Make 1 shell head model from brain model and gray matter segment [B, V, F, xx] = vb_make_head_surf_model(parm) %% --- Required Parameter parm.analyze_file : MRI image file parm.brain_file : Brain model file (.brain.mat) parm.gray_file : Gray matter segment file made by SPM 'Segment' SPM2 - 'Segment' output 3 image files for gray/white/CSF Only gray matter file (*_seg1.hdr/img) is necessary SPM5 might output different naming files %% --- Optional Parameter parm.Nvertex = 5000; # of head vertex - Radius of Morphology operation [mm] parm.Radius_csf = [6 6 -4 -4 ]; Radius for CSF mask smoothing %% --- Advanced Optional Parameter parm.Radius_fill = [4]; cortex fill radius parm.Radius_gray = [-2 2]; Radius to remove noise in gray matter image parm.Nloop = 200; % Iteration number for fit boundary parm.Nlast = 10 ; % Iteration number for smoothing by spring model parm.Glevel = [0.8]; % level of gray matter selection threshold Masa-aki Sato 2009-9-30 [minor] 2011-02-22 taku-y Relax MRI size check: MATLAB returns a very little (<1e-8) value for sum(abs(Vsize-vsize)) even though Vsize=[1 1 1] and vsize=[1.0 1.0 1.0], likely to be numerical error. Copyright (C) 2011, ATR All Rights Reserved. License : New BSD License(see VBMEG_LICENSE.txt)
0001 function [B, V, F, xx] = vb_make_head_surf_model(parm) 0002 % Make 1 shell head model from brain model and gray matter segment 0003 % [B, V, F, xx] = vb_make_head_surf_model(parm) 0004 %%% --- Required Parameter 0005 % parm.analyze_file : MRI image file 0006 % parm.brain_file : Brain model file (.brain.mat) 0007 % parm.gray_file : Gray matter segment file made by SPM 'Segment' 0008 % SPM2 - 'Segment' output 3 image files for gray/white/CSF 0009 % Only gray matter file (*_seg1.hdr/img) is necessary 0010 % SPM5 might output different naming files 0011 % 0012 %%% --- Optional Parameter 0013 % parm.Nvertex = 5000; # of head vertex 0014 % 0015 % - Radius of Morphology operation [mm] 0016 % parm.Radius_csf = [6 6 -4 -4 ]; Radius for CSF mask smoothing 0017 % 0018 %%% --- Advanced Optional Parameter 0019 % parm.Radius_fill = [4]; cortex fill radius 0020 % parm.Radius_gray = [-2 2]; Radius to remove noise in gray matter image 0021 % 0022 % parm.Nloop = 200; % Iteration number for fit boundary 0023 % parm.Nlast = 10 ; % Iteration number for smoothing by spring model 0024 % 0025 % parm.Glevel = [0.8]; % level of gray matter selection threshold 0026 % 0027 % Masa-aki Sato 2009-9-30 0028 % [minor] 2011-02-22 taku-y 0029 % Relax MRI size check: MATLAB returns a very little (<1e-8) value for 0030 % sum(abs(Vsize-vsize)) even though Vsize=[1 1 1] and vsize=[1.0 1.0 0031 % 1.0], likely to be numerical error. 0032 % 0033 % Copyright (C) 2011, ATR All Rights Reserved. 0034 % License : New BSD License(see VBMEG_LICENSE.txt) 0035 0036 Glevel = [0.8]; 0037 Radius_fill = [4]; % cortex fill 0038 Radius_gray = [-2 2]; % gray matter 0039 Radius_csf = [6 6 -4 -4 ]; % mask smoothing 0040 0041 mode = 1; % = 1: make head surface, = 0: only make mask file 0042 step = 2; % subsampling step for mask image [mm] 0043 0044 if isfield(parm,'Glevel'), Glevel = parm.Glevel; end; 0045 if isfield(parm,'Radius_fill'), Radius_fill = parm.Radius_fill; end; 0046 if isfield(parm,'Radius_gray'), Radius_gray = parm.Radius_gray; end; 0047 if isfield(parm,'Radius_csf'), Radius_csf = parm.Radius_csf; end; 0048 if isfield(parm,'vstep'), step = parm.vstep; end; 0049 if isfield(parm,'test_mode'), mode = parm.test_mode; end; 0050 0051 % Load cortical surface 0052 warning('off', 'MATLAB:load:variableNotFound'); 0053 load(parm.brain_file, 'brain_parm'); 0054 warning('on', 'MATLAB:load:variableNotFound'); 0055 0056 if exist('brain_parm', 'var') 0057 % V2 format 0058 [V0L,F0L,n0L,V0R,F0R,n0R] = vb_load_orig_brain_surf(brain_parm); 0059 [FL,VL] = vb_reducepatch(F0L, V0L, 10000); 0060 [FR,VR] = vb_reducepatch(F0R, V0R, 10000); 0061 V = [VL; VR]; 0062 V = V/1000; 0063 NdipoleL = size(VL,1); 0064 F.F3 = [FL; FR + NdipoleL]; 0065 F.F3R = FR; 0066 F.F3L = FL; 0067 F.NdipoleL = NdipoleL; 0068 0069 else 0070 % V1 format 0071 [V,F] = vb_load_cortex(parm.brain_file); 0072 end 0073 0074 % convert cortex to analyze_right_mm coordinate 0075 [Vdim, Vsize] = analyze_hdr_read(parm.analyze_file); 0076 V = vb_spm_right_to_analyze_right_mm(V, Vdim, Vsize); 0077 0078 if isfield(parm,'gray_file') && ~isempty(parm.gray_file) 0079 % Load gray matter segment image 0080 [Bg, vdim, vsize] = vb_load_analyze_to_right(parm.gray_file); 0081 0082 if sum(abs(Vsize-vsize))>=1e-3 || sum(abs(Vdim-vdim))>=1e-3 0083 error('MRI image and segment image is not the same size') 0084 end 0085 0086 % Make mask pattern from gray matter segment file 0087 Bg = vb_image_to_mask2(Bg, Glevel, step, Vsize); 0088 Bg = vb_morphology_operation(Bg, Radius_gray, step); 0089 Dim = size(Bg); 0090 else 0091 Bg = []; 0092 Dim = ceil(Vdim.*Vsize/step); 0093 end 0094 0095 % fill inside of cortex surface 0096 % B : voxcel mask pattern for cortex 0097 fprintf('Fill Cortex \n') 0098 [B] = vb_cortex_fill(V,F,step,Radius_fill,'LR', Dim); 0099 0100 if mode < 0, return; end; 0101 0102 if ~isempty(Bg) 0103 % merge mask pattern 0104 if vb_matlab_version > 6, 0105 B = int8(B) + int8(Bg); 0106 else 0107 B = double(B) + double(Bg); 0108 end 0109 end 0110 0111 if nargout==1, 0112 % Morphology smoothing 0113 ix = find(B(:) > 0); 0114 B(ix) = 1; 0115 B = vb_morphology_operation(B, Radius_csf, step); 0116 else 0117 % Morphology smoothing and surface extraction 0118 parm.vstep = step; 0119 parm.Radius = Radius_csf; 0120 0121 [V, F, xx, B] = vb_mask_to_surf_expand(B, parm); 0122 V = vb_analyze_right_mm_to_spm_right(V, Vdim, Vsize); 0123 end 0124 0125 if mode == 0 0126 Nslice = 9; 0127 zmax = round(max(V(:,3))/step); 0128 zmin = round(min(V(:,3))/step); 0129 zstep = fix((zmax-zmin)/Nslice); 0130 zindx = (1:Nslice)*zstep + zmin; 0131 vb_plot_slice( B, V/step, zindx); 0132 end