mne_read_surfaces

 [surfs] = mne_read_surfaces(surfname,read_curv,read_left,read_right,subject,subjects_dir,add_info)

 Reads FreeSurfer surface files for both hemispheres
 as well as curvatures if requested.

 Adds the derived geometry information to the surfaces

 surfname     - The name of the surface to read, e.g., 'pial'
 read_curv    - read the curvatures as well
 read_left    - read the left hemisphere (default true)
 read_right   - read the right hemisphere (default true)
 subject      - The name of the subject (defaults to SUBJECT environment
                variable)
 subjects_dir - The name of the MRI data directory (defaults to
                SUBJECTS_DIR environment variable)
 add_info     - Add auxilliary information to the surfaces
                (vertex normals, triangle centroids, triangle normals, triangle
                areas) (default true)

 surfs          - Output structure containing the surfaces

0001 function [surfs] = mne_read_surfaces(surfname,read_curv,read_left,read_right,subject,subjects_dir,add_info)
0002 %
0003 % [surfs] = mne_read_surfaces(surfname,read_curv,read_left,read_right,subject,subjects_dir,add_info)
0004 %
0005 % Reads FreeSurfer surface files for both hemispheres
0006 % as well as curvatures if requested.
0007 %
0008 % Adds the derived geometry information to the surfaces
0009 %
0010 % surfname     - The name of the surface to read, e.g., 'pial'
0011 % read_curv    - read the curvatures as well
0012 % read_left    - read the left hemisphere (default true)
0013 % read_right   - read the right hemisphere (default true)
0014 % subject      - The name of the subject (defaults to SUBJECT environment
0015 %                variable)
0016 % subjects_dir - The name of the MRI data directory (defaults to
0017 %                SUBJECTS_DIR environment variable)
0018 % add_info     - Add auxilliary information to the surfaces
0019 %                (vertex normals, triangle centroids, triangle normals, triangle
0020 %                areas) (default true)
0021 %
0022 % surfs          - Output structure containing the surfaces
0023 %
0024 
0025 %
0026 %
0027 %   Author : Matti Hamalainen, MGH Martinos Center
0028 %   License : BSD 3-clause
0029 %
0030 %
0031 %   Revision 1.8  2006/09/27 13:07:25  msh
0032 %   Added vertex normals to BEM surfaces.
0033 %   Fixed some data type casts in mne_read_surfaces
0034 %
0035 %   Revision 1.7  2006/09/25 19:48:16  msh
0036 %   Added projection item kinds to fiff_define_constants
0037 %   Changed some fields to int32 in surface structures
0038 %
0039 %   Revision 1.6  2006/09/15 16:03:02  msh
0040 %   Added an extra argument to mne_read_surfaces to avoid generation of auxilliary data
0041 %
0042 %   Revision 1.5  2006/05/30 18:29:21  msh
0043 %   Return an empty variable if neither read_left nor read_right is set
0044 %
0045 %   Revision 1.4  2006/05/30 18:24:56  msh
0046 %   Added options to read only one hemisphere.
0047 %
0048 %   Revision 1.3  2006/05/22 11:01:47  msh
0049 %   Deleted superfluous text from the comments.
0050 %
0051 %   Revision 1.2  2006/05/22 10:55:02  msh
0052 %   Fixed help text.
0053 %
0054 %   Revision 1.1  2006/05/22 10:44:44  msh
0055 %   Added surface and curvature reading routines.
0056 %   Fixed error in mne_read_source_spaces: triangle normals were not normalized
0057 %
0058 
0059 me='MNE:mne_read_surfaces';
0060 
0061 global FIFF;
0062 if isempty(FIFF)
0063     FIFF = fiff_define_constants();
0064 end
0065 
0066 if nargin < 7
0067     add_info = true;
0068     if nargin < 6
0069         subjects_dir=getenv('SUBJECTS_DIR');
0070         if isempty(subjects_dir)
0071             error(me,'SUBJECTS_DIR not set');
0072         end
0073         if nargin < 5
0074             subject=getenv('SUBJECT');
0075             if isempty(subject)
0076                 error(me,'SUBJECT not set');
0077             end
0078         end
0079         if nargin < 4
0080             read_right = true;
0081         end
0082         if nargin < 3
0083             read_left = true;
0084         end
0085     end
0086 elseif nargin ~= 7
0087     error(me,'Incorrect number of arguments');
0088 end
0089 %
0090 %   Read both hemispheres
0091 %
0092 if read_left
0093     [ lhvert, lhtri ] = mne_read_surface(sprintf('%s/%s/surf/lh.%s',subjects_dir,subject,surfname));
0094     fprintf(1,'\n');
0095 end
0096 if read_right
0097     [ rhvert, rhtri ] = mne_read_surface(sprintf('%s/%s/surf/rh.%s',subjects_dir,subject,surfname));
0098     fprintf(1,'\n');
0099 end
0100 %
0101 if read_curv
0102     if read_left
0103         lhcurv = mne_read_curvature(sprintf('%s/%s/surf/lh.curv', ...
0104             subjects_dir,subject));
0105     end
0106     if read_right
0107         rhcurv = mne_read_curvature(sprintf('%s/%s/surf/rh.curv', ...
0108             subjects_dir,subject));
0109     end
0110     fprintf(1,'\n');
0111 else
0112     lhcurv = [];
0113     rhcurv = [];
0114 end
0115 %
0116 %   Compose the source space structure
0117 %
0118 nsurf = 0;
0119 if read_left
0120     this.id           = int32(FIFF.FIFFV_MNE_SURF_LEFT_HEMI);
0121     this.np           = int32(size(lhvert,1));
0122     this.ntri         = int32(size(lhtri,1));
0123     this.coord_frame  = int32(FIFF.FIFFV_COORD_MRI);
0124     this.rr           = lhvert;
0125     this.nn           = [];
0126     this.tris         = cast(lhtri,'int32');
0127     this.nuse         = this.np;
0128     this.inuse        = ones(1,this.np);
0129     this.vertno       = [ 1 : this.np ];
0130     this.curv         = lhcurv;
0131     %
0132     %   Add the derived geometry data
0133     %
0134     hemi='left';
0135     if add_info == true
0136         complete_surface_info;
0137     end
0138     nsurf = nsurf + 1;
0139     surfs(nsurf) = this;
0140     clear('this');
0141 end
0142 if read_right
0143     this.id           = int32(FIFF.FIFFV_MNE_SURF_RIGHT_HEMI);
0144     this.np           = int32(size(rhvert,1));
0145     this.ntri         = int32(size(rhtri,1));
0146     this.coord_frame  = int32(FIFF.FIFFV_COORD_MRI);
0147     this.rr           = rhvert;
0148     this.nn           = [];
0149     this.tris         = cast(rhtri,'int32');
0150     this.nuse         = int32(this.np);
0151     this.inuse        = int32(ones(1,this.np));
0152     this.vertno       = int32([ 1 : this.np ]);
0153     this.curv         = rhcurv;
0154     hemi='right';
0155     if add_info == true
0156         complete_surface_info;
0157     end
0158     nsurf = nsurf + 1;
0159     surfs(nsurf) = this;
0160     clear('this');
0161 end
0162 
0163 if nsurf == 0
0164     surfs = [];
0165 end
0166 
0167 return;
0168 
0169     function complete_surface_info
0170         %
0171         %   Main triangulation
0172         %
0173         fprintf(1,'\tCompleting triangulation info for the %s hemisphere...',hemi);
0174         this.tri_area = zeros(1,this.ntri);
0175         fprintf(1,'triangle normals...');
0176         r1 = this.rr(this.tris(:,1),:);
0177         r2 = this.rr(this.tris(:,2),:);
0178         r3 = this.rr(this.tris(:,3),:);
0179         this.tri_cent = (r1+r2+r3)/3.0;
0180         this.tri_nn   = cross((r2-r1),(r3-r1));
0181         %
0182         %   Triangle normals and areas
0183         %
0184         for p = 1:this.ntri
0185             size = sqrt(this.tri_nn(p,:)*this.tri_nn(p,:)');
0186             this.tri_area(p) = size/2.0;
0187             if size > 0.0
0188                 this.tri_nn(p,:) = this.tri_nn(p,:)/size;
0189             end
0190         end
0191         %
0192         %   Accumulate the vertex normals
0193         %
0194         fprintf(1,'vertex normals...');
0195         this.nn = zeros(this.np,3);
0196         for p = 1:this.ntri
0197             this.nn(this.tris(p,:),:) = this.nn(this.tris(p,:),:) + kron(ones(3,1),this.tri_nn(p,:));
0198         end
0199         %
0200         %   Compute the lengths of the vertex normals and scale
0201         %
0202         fprintf(1,'normalize...');
0203         for p = 1:this.np
0204             size = sqrt(this.nn(p,:)*this.nn(p,:)');
0205             if size > 0
0206                 this.nn(p,:) = this.nn(p,:)/size;
0207             end
0208         end
0209         fprintf(1,'[done]\n');
0210     end
0211 end