EEG potencial for multilayer model by BEM [basis, basis_parm] = vb_job_leadfield_eeg_bem(proj_root, basis_parm) basis : EEG Lead field matrix ( NP x NS) basis_parm.brain_file : brain file (.brain.mat) basis_parm.area_file : area file (.area.mat) basis_parm.meg_file : MEG data file (.meg.mat) basis_parm.Basis_mode : number of independent curent direction = 1 : current vextor = xx(n,:) = 2 : current vextor = orthogonal vectors to V(n,:) = 3 : current vextor = xx(n,:) & orthogonal vectors to xx(n,:) basis_parm.bem_mode : Leadfield calculation method MEG 0 : MEG spherical model (Sarvas) 1 : MEG BEM 1 shell model 2 : MEG spherical harmonics expansion EEG 3 : EEG spherical 3 shell model 4 : EEG BEM 3 shell model basis_parm.sigma : Conductivity of each region basis_parm.head_file : Head surface file (.head.mat) --- Optional Input basis_parm.normal_mode : = 0 : normal vector at the vertex (Defualt) = 1 : average normal vector in the neighbor of BV original brain basis_parm.area_key : Area key to specify possible current region basis_parm.Recalc : Recalculation force flag Copyright (C) 2011, ATR All Rights Reserved. License : New BSD License(see VBMEG_LICENSE.txt)
0001 function [basis, basis_parm] = vb_job_leadfield_eeg_bem(proj_root, basis_parm) 0002 % EEG potencial for multilayer model by BEM 0003 % [basis, basis_parm] = vb_job_leadfield_eeg_bem(proj_root, basis_parm) 0004 % basis : EEG Lead field matrix ( NP x NS) 0005 % 0006 % basis_parm.brain_file : brain file (.brain.mat) 0007 % basis_parm.area_file : area file (.area.mat) 0008 % basis_parm.meg_file : MEG data file (.meg.mat) 0009 % basis_parm.Basis_mode : number of independent curent direction 0010 % = 1 : current vextor = xx(n,:) 0011 % = 2 : current vextor = orthogonal vectors to V(n,:) 0012 % = 3 : current vextor = xx(n,:) & orthogonal vectors to xx(n,:) 0013 % basis_parm.bem_mode : Leadfield calculation method 0014 % MEG 0015 % 0 : MEG spherical model (Sarvas) 0016 % 1 : MEG BEM 1 shell model 0017 % 2 : MEG spherical harmonics expansion 0018 % EEG 0019 % 3 : EEG spherical 3 shell model 0020 % 4 : EEG BEM 3 shell model 0021 % 0022 % basis_parm.sigma : Conductivity of each region 0023 % basis_parm.head_file : Head surface file (.head.mat) 0024 % --- Optional Input 0025 % basis_parm.normal_mode : 0026 % = 0 : normal vector at the vertex (Defualt) 0027 % = 1 : average normal vector in the neighbor of BV original brain 0028 % 0029 % basis_parm.area_key : Area key to specify possible current region 0030 % basis_parm.Recalc : Recalculation force flag 0031 % 0032 % Copyright (C) 2011, ATR All Rights Reserved. 0033 % License : New BSD License(see VBMEG_LICENSE.txt) 0034 0035 % --- 0036 % V : dipole location(in meters) NP x 3 0037 % xx : dipole moment direction NP x 3 0038 % pick : EEG sensors(in meters) on the scalp NS x 3 0039 % Omega : solid angle matrix for boundary surfaces 0040 % Vhead : vertex coordinate for boundary surfaces 0041 % Fhead : patch index for boundary surfaces 0042 % Sout : normal vector for boundary surfaces 0043 % 0044 % Nvertex = 各境界面の頂点インデックス 0045 % = [ start_id(1) end_id(1) ; 0046 % ... ; 0047 % start_id(Nsurf) end_id(Nsurf) ] 0048 % = [開始インデックス, 終了インデックス] 0049 % innermost to outermost 0050 % 0051 % Npatch = 各境界面の三角面インデックス 0052 % = [開始インデックス, 終了インデックス] 0053 % Sigma = 各領域の伝導率 0054 % = [sigma(1), ..., sigma(Nsurf)] 0055 % 0056 % 2008-10-9 Masa-aki Sato 0057 % structure of program is changed (algorithm is not changed) 0058 0059 % pre-processing 0060 [Basis_mode, pick, Qpick, Wsensor, V0, V, xx] = ... 0061 vb_lf_preprocess(proj_root, basis_parm); 0062 0063 % check parameters and set default values if not be set 0064 basis_parm = vb_util_check_head_shell_info(proj_root, basis_parm); 0065 0066 [V,xx] = vb_current_vector(V, xx, Basis_mode); 0067 0068 if iscell(basis_parm.head_file) 0069 headfile = basis_parm.head_file; 0070 joinfile = [vb_join_cell_file_name(headfile) '.head.mat']; 0071 joinpath = [proj_root '/' joinfile]; 0072 0073 if ~exist(joinpath,'file') 0074 sigma=[];layer=[];mode = 3; % MAX mode for radius 0075 vb_head_join_files(headfile,proj_root,joinpath, sigma,layer,mode); 0076 end 0077 basis_parm.head_file = joinfile; 0078 end 0079 0080 new_headfile = [proj_root filesep basis_parm.head_file]; 0081 omega_file = [vb_get_basename(new_headfile,'.head.mat') '.omega.mat']; 0082 0083 fprintf('Load %s\n',new_headfile) 0084 load(new_headfile); 0085 % Vhead : (# of vertex) x 3 0086 % Fhead : (# of patch) x 3 0087 % Nvertex : (# of surface) x 2 0088 % Npatch : (# of surface) x 2 0089 % Sigma : 1 x (# of surface) 0090 0091 fprintf('--- BEM for EEG\n') 0092 fprintf('# of vertices = %d\n',size(Vhead,1)); 0093 0094 if ~exist('Nvertex', 'var'), Nvertex = [1 length(Vhead)]; end; 0095 if ~exist('Npatch', 'var'), Npatch = [1 length(Fhead)]; end; 0096 0097 Nsurf = size(Nvertex,1); 0098 fprintf('# of surfaces = %d\n',Nsurf); 0099 0100 Sigma = basis_parm.sigma; 0101 0102 BEM.sigma = [Sigma(1:Nsurf) 0]; 0103 BEM.Nvertex = Nvertex; 0104 BEM.Npatch = Npatch; 0105 0106 if exist(omega_file,'file') 0107 fprintf('Load %s\n',omega_file) 0108 load(omega_file) 0109 end 0110 0111 if ~exist('Omega','var') || basis_parm.Recalc == ON, 0112 tic; 0113 fprintf('--- Solid angle calculation \n'); 0114 0115 % (Normal vectors (XXhead) are outwarded here) 0116 [Omega, Sout] = vb_solid_angle_grk(Vhead,Fhead,XXhead); 0117 0118 vb_ptime(toc); 0119 vb_fsave(omega_file,'Omega', 'Sout'); 0120 fprintf('Saved omega file: %s\n',omega_file); 0121 end 0122 0123 tic; 0124 0125 % 多層モデル補正 0126 fprintf('--- Multisurface correction\n'); 0127 [Omega] = vb_bem_matrix(Omega,BEM); 0128 0129 % 自己重み係数:基底関数を三角面上で積分した寄与を加える 0130 fprintf('--- Auto solid angle calculation\n'); 0131 Omega = vb_solid_auto_grk(Omega,Fhead,Sout); 0132 0133 vb_ptime(toc); 0134 0135 fprintf('--- EEG potential (BEM) \n'); 0136 0137 % Electric potential 0138 basis = vb_bem_eeg( Omega, Sout, Vhead, Fhead, V, xx, pick, BEM ); 0139 basis = basis'; 0140 0141 %%% END OF FILE %%%