function [coordVERTICES,varargout] = READ_stlascii(stlFILENAME)
% READ_stlascii  Read mesh data in the form of an ascii <*.stl> file
%==========================================================================
% FILENAME:          READ_stlascii.m
% AUTHOR:            Adam H. Aitkenhead
% INSTITUTION:       The Christie NHS Foundation Trust
% CONTACT:           adam.aitkenhead@physics.cr.man.ac.uk
% DATE:              29th March 2010
% PURPOSE:           Read mesh data in the form of an ascii <*.stl> file.
%
% USAGE:             [coordVERTICES] = READ_stlascii('file.stl')
%                    where coordVERTICES is an Nx3x3 array defining the
%                    vertex positions for each facet, with:
%                       - 1 row for each facet
%                       - 3 cols for the x,y,z coordinates
%                       - 3 pages for the three vertices
%
%                    [coordVERTICES,coordNORMALS] = READ_stlascii('file.stl')
%                    also outputs coordNORMALS, an Nx3 array defining the
%                    normal vector for each facet, with:
%                       - 1 row for each facet
%                       - 3 cols for the x,y,z components of the vector
%
%                    [coordVERTICES,coordNORMALS,STLname] = READ_stlascii('file.stl')
%                    also provides the string STLname defining the name of
%                    the stl mesh.
%
% STL REQUIREMENTS:  The STL data must meet the following criteria:
%                    1. All facets are triangular.
%==========================================================================

%==========================================================================
% VERSION:  USER:  CHANGES:
% -------   -----  --------
% 100329    AHA    Original version
% 100513    AHA    Totally reworked the code.  Now use textscan to read the
%                  file all at once, rather than one line at a time with
%                  fgetl.  Major speed improvment.
%==========================================================================

%==========================================================================
% STL ascii file format
%======================
% ASCII STL files have the following structure.  Technically each facet
% could be any 2D shape, but in practice only triangular facets tend to be
% used.  The present code ONLY works for meshes composed of triangular
% facets.
%
% solid object_name
% facet normal x y z
%   outer loop
%     vertex x y z
%     vertex x y z
%     vertex x y z
%   endloop
% endfacet
%
% <Repeat for all facets...>
%
% endsolid object_name
%==========================================================================

%================================
% Read the ascii STL file
%================================

fidIN = fopen(stlFILENAME);
fidCONTENTcell = textscan(fidIN,'%s','delimiter','\n');                  %Read all the file
fidCONTENT = fidCONTENTcell{:}(logical(~strcmp(fidCONTENTcell{:},'')));  %Remove all blank lines
fclose(fidIN);

%================================
% Read the STL name
%================================

if nargout == 3
  line1 = char(fidCONTENT(1));
  if (size(line1,2) >= 7)
    STLname = line1(7:end);
  else
    STLname = 'unnamed_object'; 
  end
end

%================================
% Read the vector normals
%================================

if nargout >= 2
  stringNORMALS = char(fidCONTENT(logical(strncmp(fidCONTENT,'facet normal',12))));
  coordNORMALS  = str2num(stringNORMALS(:,13:end));
end

%================================
% Read the vertex coordinates
%================================

facetTOTAL       = sum(strcmp(fidCONTENT,'endfacet'));

stringVERTICES   = char(fidCONTENT(logical(strncmp(fidCONTENT,'vertex',6))));
coordVERTICESall = str2num(stringVERTICES(:,7:end));
cotemp           = zeros(3,facetTOTAL,3);
cotemp(:)        = coordVERTICESall;
coordVERTICES    = shiftdim(cotemp,1);

%================================
% Prepare the output arguments
%================================

if nargout == 2
  varargout(1) = {coordNORMALS};
elseif nargout == 3
  varargout(1) = {coordNORMALS};
  varargout(2) = {STLname};
end
