EEG potential for multilayer model by BEM 多層境界面モデルのポテンシャルと磁場 P = vb_bem_eeg(D,Sout,Vhead,Fhead,V,J,Xeeg,BEM) OUTPUTS: P : EEG Lead field matrix ( Neeg x NP) P : 多層モデルのEEG センサポテンシャル INPUTS (Required): V : dipole location(in meters) NP x 3 J : dipole moment NP x 3 Xeeg : EEG sensors(in meters) on the scalp Neeg x 3 D : solid angle matrix for boundary surfaces Vhead : vertex coordinate for boundary surfaces Fhead : patch index for boundary surfaces Sout : normal vector for boundary surfaces D : 境界面係数行列 Vhead : 境界三角面の座標 Fhead : 境界三角面頂点インデックス Sout : 三角面の外向き面積法線 V : 電流双極子座標 (m) J : 電流双極子モーメント (A) Xeeg : EEG センサ座標 BEM.Nvertex = 各境界面の頂点インデックス = [ start_id(1) end_id(1) ; ... ; start_id(Nsurf) end_id(Nsurf) ] = [開始インデックス, 終了インデックス] innermost to outermost BEM.Npatch = 各境界面の三角面インデックス = [開始インデックス, 終了インデックス] BEM.sigma = 各領域の伝導率 = [sigma(1), ..., sigma(Nsurf), 0] 2004-12-26 M. Sato 2007-12-6 M. Sato EEG positions are mapped to outermost head surface first 2007-12-25 M. Sato EEG = Weighted average of potentials at nearest head positions 2008-10-8 M. Sato Multisurface correction & Auto solid angle calculation are done before saving omega Copyright (C) 2011, ATR All Rights Reserved. License : New BSD License(see VBMEG_LICENSE.txt)
0001 function P = vb_bem_eeg(D,Sout,Vhead,Fhead,V,J,Xeeg,BEM) 0002 % EEG potential for multilayer model by BEM 0003 % 多層境界面モデルのポテンシャルと磁場 0004 % P = vb_bem_eeg(D,Sout,Vhead,Fhead,V,J,Xeeg,BEM) 0005 % OUTPUTS: 0006 % P : EEG Lead field matrix ( Neeg x NP) 0007 % P : 多層モデルのEEG センサポテンシャル 0008 % 0009 % INPUTS (Required): 0010 % V : dipole location(in meters) NP x 3 0011 % J : dipole moment NP x 3 0012 % Xeeg : EEG sensors(in meters) on the scalp Neeg x 3 0013 % D : solid angle matrix for boundary surfaces 0014 % Vhead : vertex coordinate for boundary surfaces 0015 % Fhead : patch index for boundary surfaces 0016 % Sout : normal vector for boundary surfaces 0017 % 0018 % D : 境界面係数行列 0019 % Vhead : 境界三角面の座標 0020 % Fhead : 境界三角面頂点インデックス 0021 % Sout : 三角面の外向き面積法線 0022 % V : 電流双極子座標 (m) 0023 % J : 電流双極子モーメント (A) 0024 % Xeeg : EEG センサ座標 0025 % 0026 % BEM.Nvertex = 各境界面の頂点インデックス 0027 % = [ start_id(1) end_id(1) ; 0028 % ... ; 0029 % start_id(Nsurf) end_id(Nsurf) ] 0030 % = [開始インデックス, 終了インデックス] 0031 % innermost to outermost 0032 % 0033 % BEM.Npatch = 各境界面の三角面インデックス 0034 % = [開始インデックス, 終了インデックス] 0035 % BEM.sigma = 各領域の伝導率 0036 % = [sigma(1), ..., sigma(Nsurf), 0] 0037 % 0038 % 2004-12-26 M. Sato 0039 % 2007-12-6 M. Sato 0040 % EEG positions are mapped to outermost head surface first 0041 % 2007-12-25 M. Sato 0042 % EEG = Weighted average of potentials at nearest head positions 0043 % 2008-10-8 M. Sato 0044 % Multisurface correction & Auto solid angle calculation 0045 % are done before saving omega 0046 % 0047 % Copyright (C) 2011, ATR All Rights Reserved. 0048 % License : New BSD License(see VBMEG_LICENSE.txt) 0049 0050 0051 %% 多層モデル補正 0052 %fprintf('--- Multisurface correction\n'); 0053 %[D, Keeg] = vb_bem_matrix(D,BEM); 0054 % 0055 %% 自己重み係数:基底関数を三角面上で積分した寄与を加える 0056 %fprintf('--- Auto solid angle calculation\n'); 0057 %D = vb_solid_auto_grk(D,Fhead,Sout); 0058 0059 % Inter-vertex distance of head surface 0060 dnext = vb_next_vertex_distance(Fhead,Vhead); 0061 dmax = max(dnext); 0062 dmean = mean(dnext); 0063 0064 fprintf('Mean inter-vertex distance of head surface = %f\n',dmean); 0065 fprintf('MAX inter-vertex distance of head surface = %f\n',dmax); 0066 0067 % EEG positions are mapped to outermost head surface 0068 % IXeeg : Vertex index nearest to EEG sensor [Neeg x Np] 0069 % distance : Distance from EEG sensor to head vertex [Neeg x Np] 0070 Np = 3; % number of close surface points for each sensor 0071 [Xeeg, IXeeg, distance] = vb_map_eeg_to_head(Xeeg, Vhead, BEM, dmax, Np); 0072 0073 fprintf('Mean EEG-Vertex minimum distance = %f\n', mean(distance(:,1))); 0074 0075 % Extract valid vertex 0076 ix_next = find( IXeeg > 0 ); 0077 0078 fprintf('# of neighbor vertex for EEG = %d\n', length(ix_next)) 0079 fprintf('Mean EEG-Vertex distance = %f\n', mean(distance(ix_next))); 0080 0081 tic 0082 % Peeg : 双極子ポテンシャルからEEGポテンシャルへの変換行列 0083 fprintf('--- EEG transfer matrix calculation\n'); 0084 [Peeg] = vb_bem_inverse_eeg(D, IXeeg(ix_next), BEM); 0085 0086 vb_ptime(toc); 0087 tic 0088 0089 % ポテンシャル計算 0090 fprintf('--- EEG potential calculation\n'); 0091 Phead = vb_bem_eeg_linear_grk(Peeg,Vhead,Fhead,V,J); 0092 0093 vb_ptime(toc); 0094 0095 % EEG = Weighted average of potentials at close head positions 0096 % Phead : head vertex potential (Nhead x Ndipole) 0097 % P : EEG Lead field matrix ( Neeg x Ndipole) 0098 0099 fprintf('--- Interpolate head potential to EEG sensor position\n'); 0100 0101 [Neeg, Np] = size(IXeeg); % Neeg = # of EEG sensor 0102 NP = size(V,1); % # of dipoles 0103 P = zeros(Neeg,NP); 0104 0105 % Weight factor for average 0106 Ptemp = zeros(Neeg, Np); 0107 weight = vb_linear_interpolate3(distance); 0108 0109 for n = 1:NP 0110 Ptemp(ix_next) = Phead(:,n); 0111 P(:,n) = sum(Ptemp .* weight , 2); 0112 end