Make whole brain data (cortex sheet) by reducing BV/FS surface data. (VBMEG public function) [syntax] vb_job_brain(proj_root,brain_parm); [input] proj_root : <<string>> VBMEG project root directory. brain_parm: <<struct>> Parameters for creating cortical surface model --- fields of brain_parm brain_file : <<string>> Cortical surface model file (.brain.mat). act_file : <<string>> Cortical activity map file (.act.mat). area_file : <<string>> Cortical area file (.area.mat). Nvertex : <<int>> Total number of vertices for cortical surface model Priority: Nvertex > reduce_ratio reduce_ratio: <<float>> Reduce ratio for cortical surface model Rmax : <<float>> Maximum radius for neighbor search N_step : <<int>> Division number of z axis in matching between cortical surface model and MRI points display : <<int>> Display progress for this steps make_reduce : <bool?> Making reduced cortex analyze_file: <optional> <<string>> MRI filename (.hdr) spm_normalization_file: <optional> <<string>> SPM normalization file (_sn.mat) (parameters for BrainVoyager cortical surface model) <<string>> BV_left_file : GM/WM boundary file (.srf) for left hemisphere BV_right_file : GM/WM boundary file (.srf) for right hemisphere (parameters for FreeSurfer cortical surface model) <<string>> FS_left_file : GM/WM boundary file (.smoothwm.asc) for LH FS_right_file: GM/WM boundary file (.smoothwm.asc) for RH --- [output] The following files are created: VBMEG cortical surface model file (.brain.mat) VBMEG area file (.area.mat) VBMEG activity map file (.act.mat) The following variables are stored in .brain.mat file: V : Cortical vertex point cordinate (SPM_Right_m) F : Patch index structure xx : Normal vector to cortical surface xxA : area asigned for each vertex ([m^2]) BV_index : original vertex index in BV/FS corresponding to brain Vinfo : Vertex dimension structure xxF{i} : Next-neighbor index for the vertex-i xxD{i} : Distance between next-neighbor and the vertex-i xxT{n} : vertex index for n-th triangle nextIX{i} : Neighbor index list for the vertex-i nextDD{i} : Distance from the vertex-i The following areas are automatically created and stored in .area.mat file: areakey = 'Cortex' : whole brain areakey = 'reduced_cortex' : reduced cortex without corpus region areakey = 'cortex_region' : cortex region without corpus area areakey = 'corpus_area' : corpus area [history] Line 117 was modified by Kawawaki, Dai at 06/04/07 2006/07/21 M. Sato Make reduced cortex for soft normal constraint 2007/03/05 O. Yamashita Add lines to compute MNI and Talairach coordinate and two area files derived from anatomical divisions. These lines are only valid when spm_normalization_file is specified. 2008-12-24 Taku Yoshioka Reduce model supression mode vb_disp() is used for displaying message Show figures for checking standard brain mask 2010-05-28 Taku Yoshioka Minor change 2011-01-28 rhayashi call vb_job_inflate() when brain_parm have fields for inflate model. Copyright (C) 2011, ATR All Rights Reserved. License : New BSD License(see VBMEG_LICENSE.txt)
0001 function [V,F,xx,BV_index,Vinfo] = vb_job_brain(proj_root,brain_parm, mode) 0002 % Make whole brain data (cortex sheet) by reducing BV/FS surface data. 0003 % (VBMEG public function) 0004 % 0005 % [syntax] 0006 % vb_job_brain(proj_root,brain_parm); 0007 % 0008 % [input] 0009 % proj_root : <<string>> VBMEG project root directory. 0010 % brain_parm: <<struct>> Parameters for creating cortical surface model 0011 % --- fields of brain_parm 0012 % brain_file : <<string>> Cortical surface model file (.brain.mat). 0013 % act_file : <<string>> Cortical activity map file (.act.mat). 0014 % area_file : <<string>> Cortical area file (.area.mat). 0015 % Nvertex : <<int>> Total number of vertices for cortical surface model 0016 % Priority: Nvertex > reduce_ratio 0017 % reduce_ratio: <<float>> Reduce ratio for cortical surface model 0018 % Rmax : <<float>> Maximum radius for neighbor search 0019 % N_step : <<int>> Division number of z axis in matching between 0020 % cortical surface model and MRI points 0021 % display : <<int>> Display progress for this steps 0022 % make_reduce : <bool?> Making reduced cortex 0023 % analyze_file: <optional> <<string>> MRI filename (.hdr) 0024 % spm_normalization_file: <optional> <<string>> SPM normalization file 0025 % (_sn.mat) 0026 % 0027 % (parameters for BrainVoyager cortical surface model) <<string>> 0028 % BV_left_file : GM/WM boundary file (.srf) for left hemisphere 0029 % BV_right_file : GM/WM boundary file (.srf) for right hemisphere 0030 % 0031 % (parameters for FreeSurfer cortical surface model) <<string>> 0032 % FS_left_file : GM/WM boundary file (.smoothwm.asc) for LH 0033 % FS_right_file: GM/WM boundary file (.smoothwm.asc) for RH 0034 % 0035 % --- 0036 % 0037 % [output] 0038 % The following files are created: 0039 % VBMEG cortical surface model file (.brain.mat) 0040 % VBMEG area file (.area.mat) 0041 % VBMEG activity map file (.act.mat) 0042 % 0043 % The following variables are stored in .brain.mat file: 0044 % V : Cortical vertex point cordinate (SPM_Right_m) 0045 % F : Patch index structure 0046 % xx : Normal vector to cortical surface 0047 % xxA : area asigned for each vertex ([m^2]) 0048 % 0049 % BV_index : original vertex index in BV/FS corresponding to brain 0050 % Vinfo : Vertex dimension structure 0051 % 0052 % xxF{i} : Next-neighbor index for the vertex-i 0053 % xxD{i} : Distance between next-neighbor and the vertex-i 0054 % xxT{n} : vertex index for n-th triangle 0055 % nextIX{i} : Neighbor index list for the vertex-i 0056 % nextDD{i} : Distance from the vertex-i 0057 % 0058 % The following areas are automatically created and stored in .area.mat 0059 % file: 0060 % areakey = 'Cortex' : whole brain 0061 % areakey = 'reduced_cortex' : reduced cortex without corpus region 0062 % areakey = 'cortex_region' : cortex region without corpus area 0063 % areakey = 'corpus_area' : corpus area 0064 % 0065 % [history] 0066 % Line 117 was modified by Kawawaki, Dai at 06/04/07 0067 % 2006/07/21 M. Sato 0068 % Make reduced cortex for soft normal constraint 0069 % 2007/03/05 O. Yamashita 0070 % Add lines to compute MNI and Talairach coordinate and two area files 0071 % derived from anatomical divisions. These lines are only valid when 0072 % spm_normalization_file is specified. 0073 % 2008-12-24 Taku Yoshioka 0074 % Reduce model supression mode 0075 % vb_disp() is used for displaying message 0076 % Show figures for checking standard brain mask 0077 % 2010-05-28 Taku Yoshioka 0078 % Minor change 0079 % 2011-01-28 rhayashi 0080 % call vb_job_inflate() when brain_parm have fields for inflate model. 0081 % 0082 % Copyright (C) 2011, ATR All Rights Reserved. 0083 % License : New BSD License(see VBMEG_LICENSE.txt) 0084 0085 global vbmeg_inst 0086 const = vb_define_verbose; 0087 VERBOSE_LEVEL_NOTICE = const.VERBOSE_LEVEL_NOTICE; 0088 0089 % Set BV-SPM version (mode = 0 : SBI-compatible) 0090 if ~exist('mode','var'), mode = 1; end; 0091 0092 % Define constant 0093 define = vbmeg_inst.const; 0094 0095 proj_root = vb_rm_trailing_slash(proj_root); 0096 0097 % 0098 % Do not modify following lines 0099 % 0100 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0101 % 0102 0103 % Output file path 0104 brain_file = [proj_root filesep brain_parm.brain_file]; 0105 area_file = [proj_root filesep brain_parm.area_file ]; 0106 act_file = [proj_root filesep brain_parm.act_file ]; 0107 0108 % ID for MRI image data 0109 [udir, fname] = fileparts(brain_parm.analyze_file); 0110 MRI_ID = fname; 0111 0112 %------------------------------------------------------- 0113 % Make cortex coordinate from Brain-Voyager srf-file 0114 %------------------------------------------------------- 0115 switch mode 0116 case 0 0117 % BV-Vox-SPM version (SBI-compatible) 0118 [V,F,xx,BV_index,Vinfo] = vb_make_brain_data(brain_parm); 0119 DEBUG_MODE = 0; 0120 case 1 0121 % BV-SPM version 0122 [V,F,xx,BV_index,Vinfo] = vb_make_brain_data2(brain_parm); 0123 DEBUG_MODE = 0; 0124 case 2 0125 % BV-Vox-SPM version (SBI-compatible) 0126 [V,F,xx,BV_index,Vinfo] = vb_make_brain_data(brain_parm); 0127 DEBUG_MODE = 2; 0128 case 3 0129 % BV-SPM version 0130 [V,F,xx,BV_index,Vinfo] = vb_make_brain_data2(brain_parm); 0131 DEBUG_MODE = 2; 0132 end 0133 0134 if DEBUG_MODE==2, return; end; 0135 0136 %---------------------------------------- 0137 % Search next-point index and distance 0138 %---------------------------------------- 0139 tic; 0140 vb_disp_nonl('Search next-point index and distance ', ... 0141 VERBOSE_LEVEL_NOTICE); 0142 0143 [xxD, xxF, xxT] = vb_next_distance( F.F3, V ); % unit : m 0144 0145 vb_disp(sprintf('%f[sec]',toc), VERBOSE_LEVEL_NOTICE); 0146 0147 % total number of vertex 0148 Ndipole = Vinfo.Ndipole; 0149 Vindx = 1:Ndipole; 0150 0151 %--------------------------------------- 0152 % Smoothing normal vectors (Double mean) 0153 %--------------------------------------- 0154 tic 0155 vb_disp_nonl('Mean of patch normal ', VERBOSE_LEVEL_NOTICE); 0156 0157 % Mean of normal vector for neighboring patch 0158 xx = vb_mean_normal_vector_1(Vindx ,F,V,xx); 0159 % Mean of normal vector for neighbor vertex 0160 xx = vb_mean_normal_vector_4(Vindx ,F,V,xx,xxF); 0161 0162 vb_disp(sprintf('%f[sec]',toc), VERBOSE_LEVEL_NOTICE); 0163 0164 %---------------------------------------------- 0165 % Calculate area assigned to the vertex 0166 %---------------------------------------------- 0167 tic; 0168 [xxA] = vb_calc_patch_area(V, F.F3, xxT); 0169 vb_disp(sprintf('%f[sec]',toc), VERBOSE_LEVEL_NOTICE); 0170 0171 % 0172 % Save cortical surface model 0173 % 0174 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 0175 vb_disp(['Save brain model '], VERBOSE_LEVEL_NOTICE); 0176 vb_disp(sprintf(' filename = %s', brain_file), ... 0177 VERBOSE_LEVEL_NOTICE); 0178 vb_save([brain_file],... 0179 'F','V','xx','xxF','xxD','xxA','Vinfo','BV_index','MRI_ID'); 0180 0181 if DEBUG_MODE==1, return; end; 0182 0183 %---------------------------------------------- 0184 % Search neighbor points along cortex sheet 0185 %---------------------------------------------- 0186 Rmax = brain_parm.R_max; 0187 Display = brain_parm.display; 0188 0189 [nextIX , nextDD] = vb_find_neighbor_all(Rmax, xxF, xxD, Vindx, Display); 0190 0191 vb_disp('Save neighbor index of cortical vertex', ... 0192 VERBOSE_LEVEL_NOTICE); 0193 vb_save([brain_file],'nextIX','nextDD'); 0194 %save([brain_file],'nextIX','nextDD','-append'); 0195 0196 %-------------------------------------------------------------------------- 0197 % Calculate MNI coordinate and Talairach coordinate and 0198 % two area files derived from two anatomical divisons (added by OY 07/03/05) 0199 %-------------------------------------------------------------------------- 0200 if isfield(brain_parm, 'spm_normalization_file') ... 0201 & ~isempty(brain_parm.spm_normalization_file) 0202 0203 atlas_dir = 'MNI_atlas_templates/'; 0204 EXT_brain = '.brain.mat'; 0205 brain_id = brain_file(1:findstr(brain_file,EXT_brain)-1); 0206 0207 %% MNI coordinate 0208 tic 0209 vb_disp(sprintf('-- Computing MNI and Talairach coordinate '), ... 0210 VERBOSE_LEVEL_NOTICE); 0211 0212 mniparm.mask_file = ['betmask.img']; % fixed 0213 % mniparm.brainfile = brain_parm.brain_file; 0214 mniparm.brainfile = brain_file; 0215 mniparm.snfile = brain_parm.spm_normalization_file; 0216 0217 [Vmni, Vtal, derr, IMGinfo, pXYZmni0] = vb_calc_mni_coord(mniparm); 0218 0219 vb_save([brain_file],'Vmni','Vtal','derr'); 0220 vb_disp_nonl(sprintf('Done...%f[sec]\n\n',toc), VERBOSE_LEVEL_NOTICE); 0221 0222 % Display brain mask superimposed on MRI 0223 %[V, F] = vb_load_cortex(brainfile); 0224 %load(brainfile,'Info','Vmni'); 0225 %load(brainfile,'Vmni'); 0226 %Vmni0 = Vmni; 0227 %Vmni = Info.Vmni; 0228 %ixL = Info.ixL; 0229 %ixR = Info.ixR; 0230 0231 % Load analyze data 0232 [B, Vdim, Vsize] = vb_load_analyze_to_right(brain_parm.analyze_file); 0233 0234 % Transform vertex coordinate to analyze right-handed coordinate 0235 pXYZmni = vb_spm_right_to_analyze_right(pXYZmni0./1000,Vdim,Vsize); 0236 0237 Msize = 1; % marker size 0238 Mtype = 'y-'; % marker type 0239 dz = 5; % search radius for intersection triangle 0240 vcut = ['x','y','z','z']; % slice cut direction 0241 indx = [80 120 130 170 ]; % slice number to cut 0242 xymode = 0; 0243 0244 label = {'Sagittal cut'; 'Coronal cut'; 'Axial cut'; 'Axial cut'}; 0245 0246 Nfig = length(vcut); 0247 NX = 2; 0248 NY = 2; 0249 0250 figure('Name','Brain mask superimposed on MR images'); 0251 0252 for n= 1:Nfig 0253 subplot(NY,NX,n) 0254 vdim = vcut(n); 0255 [strX,strY] = vb_plot_3d_image(B, indx(n), vcut(n), xymode); 0256 xlabel(strX); 0257 ylabel(strY); 0258 0259 hold on 0260 vb_plot_vertex(pXYZmni,vcut(n),indx(n),dz,Msize,Mtype,xymode); 0261 title(label{n}) 0262 end 0263 0264 % Display histogram of x-values of MNI coordinate 0265 ixL = unique(F.F3L(:)); 0266 ixR = unique(F.F3R(:)); 0267 %save('tmp.mat','pXYZmni0','pXYZmni','Vmni'); 0268 %error('aho'); 0269 %Info.Vmni = pXYZmni; 0270 0271 figure('Name','Histogram of MNI coordinate (x-axis)'); 0272 x = -0.1:0.001:0.1; 0273 subplot(2,1,1); 0274 title('LH') 0275 hist(Vmni(ixL,1),x); 0276 hold on; 0277 plot([0 0],[0 800],'r'); 0278 subplot(2,1,2); 0279 title('RH'); 0280 hist(Vmni(ixR,1),x); 0281 hold on; 0282 plot([0 0],[0 800],'r'); 0283 0284 %% AAL atlas 0285 tic 0286 vb_disp_nonl(['Creating an area file consisting of AAL anatomical' ... 0287 'division: '], VERBOSE_LEVEL_NOTICE); 0288 atlas_id = 'aal'; 0289 0290 mniparm.atlas_text = [atlas_dir atlas_id '.txt']; 0291 mniparm.atlas_file = [atlas_dir atlas_id '.img']; 0292 mniparm.atlas_id = atlas_id; 0293 mniparm.atlas_dir = atlas_dir; 0294 0295 0296 mniparm.save_areafile = [brain_id '_' atlas_id '.area.mat']; 0297 mniparm.save_atlasfile = [brain_id '_atlas.act.mat']; % Atlas label 0298 0299 vb_atlas2vb_aal(mniparm) 0300 0301 vb_disp_nonl(sprintf('Done...%f[sec]\n\n',toc)); 0302 0303 0304 %% Brodmann atlas 0305 tic 0306 vb_disp_nonl(['Creating an area file consisting of Brodmann ' ... 0307 'anatomical division: ']); 0308 0309 atlas_id = 'brodmann'; 0310 0311 mniparm.atlas_text = [atlas_dir atlas_id '.txt']; 0312 mniparm.atlas_file = [atlas_dir atlas_id '.img']; 0313 mniparm.atlas_id = atlas_id; 0314 0315 mniparm.save_areafile = [brain_id '_' atlas_id '.area.mat']; 0316 0317 vb_atlas2vb(mniparm) 0318 0319 vb_disp_nonl(sprintf('Done...%f[sec]\n\n',toc), VERBOSE_LEVEL_NOTICE); 0320 else 0321 vb_disp('MNI and Talairach coordinate are not computed.', ... 0322 VERBOSE_LEVEL_NOTICE); 0323 vb_disp(['You can obtain MNI and Talairach coordinate by ' ... 0324 'running vb_atlas2vb.m later '], VERBOSE_LEVEL_NOTICE); 0325 end 0326 0327 %-------------------------------------- 0328 % Make new area and act files 0329 %-------------------------------------- 0330 act_new.key = 'Uniform'; 0331 act_new.xxP = ones(Ndipole,1); 0332 act_new.comment = 'homogeneous fMRI constraint'; 0333 0334 vb_add_act([act_file], act_new, MRI_ID, OFF); 0335 0336 % Make Area file 0337 AreaNew.key = 'Cortex'; 0338 AreaNew.Iextract = [1:Ndipole]'; %modified by Kawawaki, Dai at 06/04/07 0339 0340 vb_add_area([area_file], AreaNew, MRI_ID, OFF); 0341 0342 %----------------------------------------- 0343 % Check brain model by 3D & MRI image 0344 %----------------------------------------- 0345 vb_check_brain_model(proj_root,brain_parm); 0346 0347 %--------------------------------------------------- 0348 % Find original vertex near the cortical vertex 0349 %--------------------------------------------------- 0350 normal_stat = vb_original_normal_statics(proj_root,brain_parm); 0351 0352 vb_disp('Save neighbor statics in the original brain', ... 0353 VERBOSE_LEVEL_NOTICE); 0354 vb_save(brain_file,'normal_stat'); 0355 0356 %-------------------------- 0357 % Make reduced cortex 0358 %-------------------------- 0359 if isfield(brain_parm,'make_reduce') & brain_parm.make_reduce, 0360 vb_job_reduced_cortex(proj_root,brain_parm); 0361 end 0362 0363 %-------------------------- 0364 % Import inflate model 0365 %-------------------------- 0366 if has_inflate_model(brain_parm) 0367 vb_job_inflate(proj_root, brain_parm); 0368 end 0369 0370 %-------------------------- 0371 % Save project file 0372 %-------------------------- 0373 proj_file = get_project_filename; 0374 if isempty(proj_file) 0375 return; 0376 end 0377 0378 project_file_mgr('load', proj_file); 0379 project_file_mgr('add', 'brain_parm', brain_parm); 0380 0381 0382 function [result] = has_inflate_model(brain_parm) 0383 % Check if brain_parm contains valid fields for inflate model. 0384 % result : return bit pattern 0385 % = 01(1 as int): BV inflate data is set. 0386 % = 10(2 as int): FS inflate data is set. 0387 % = 11(3 as int): Both inflate model were set 0388 result = 0; 0389 BV_BIT = 1; FS_BIT = 2; 0390 % bit all on 0391 result = bitset(result, BV_BIT, 1); 0392 result = bitset(result, FS_BIT, 1); 0393 0394 BV_fields = {'BV_left_infl_file', 'BV_right_infl_file'}; 0395 FS_fields = {'FS_left_infl_file', 'FS_right_infl_file', ... 0396 'FS_left_curv_file', 'FS_right_curv_file'}; 0397 for k=1:length(BV_fields) 0398 if ~isfield(brain_parm, BV_fields{k}) || ... 0399 exist(brain_parm.(BV_fields{k}), 'file') ~= 2 0400 % bit off 0401 result = bitset(result, BV_BIT, 0); 0402 break; 0403 end 0404 end 0405 for k=1:length(FS_fields) 0406 if ~isfield(brain_parm, FS_fields{k}) || ... 0407 exist(brain_parm.(FS_fields{k}), 'file') ~= 2 0408 % bit off 0409 result = bitset(result, FS_BIT, 0); 0410 break; 0411 end 0412 end