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 [V,F] = vb_load_cortex(parm.brain_file); 0053 0054 % convert cortex to analyze_right_mm coordinate 0055 [Vdim, Vsize] = analyze_hdr_read(parm.analyze_file); 0056 V = vb_spm_right_to_analyze_right_mm(V, Vdim, Vsize); 0057 0058 if isfield(parm,'gray_file') && ~isempty(parm.gray_file) 0059 % Load gray matter segment image 0060 [Bg, vdim, vsize] = vb_load_analyze_to_right(parm.gray_file); 0061 0062 if sum(abs(Vsize-vsize))>=1e-3 || sum(abs(Vdim-vdim))>=1e-3 0063 error('MRI image and segment image is not the same size') 0064 end 0065 0066 % Make mask pattern from gray matter segment file 0067 Bg = vb_image_to_mask2(Bg, Glevel, step, Vsize); 0068 Bg = vb_morphology_operation(Bg, Radius_gray, step); 0069 Dim = size(Bg); 0070 else 0071 Bg = []; 0072 Dim = ceil(Vdim.*Vsize/step); 0073 end 0074 0075 % fill inside of cortex surface 0076 % B : voxcel mask pattern for cortex 0077 fprintf('Fill Cortex \n') 0078 [B] = vb_cortex_fill(V,F,step,Radius_fill,'LR', Dim); 0079 0080 if mode < 0, return; end; 0081 0082 if ~isempty(Bg) 0083 % merge mask pattern 0084 if vb_matlab_version > 6, 0085 B = int8(B) + int8(Bg); 0086 else 0087 B = double(B) + double(Bg); 0088 end 0089 end 0090 0091 if nargout==1, 0092 % Morphology smoothing 0093 ix = find(B(:) > 0); 0094 B(ix) = 1; 0095 B = vb_morphology_operation(B, Radius_csf, step); 0096 else 0097 % Morphology smoothing and surface extraction 0098 parm.vstep = step; 0099 parm.Radius = Radius_csf; 0100 0101 [V, F, xx, B] = vb_mask_to_surf_expand(B, parm); 0102 V = vb_analyze_right_mm_to_spm_right(V, Vdim, Vsize); 0103 end 0104 0105 if mode == 0 0106 Nslice = 9; 0107 zmax = round(max(V(:,3))/step); 0108 zmin = round(min(V(:,3))/step); 0109 zstep = fix((zmax-zmin)/Nslice); 0110 zindx = (1:Nslice)*zstep + zmin; 0111 vb_plot_slice( B, V/step, zindx); 0112 end