Leadfield calculation --- INPUT V(n,:) : dipole position (3-D coordinate) at n-th vertex xx(n,:) : dipole current direction (unit vector) at n-th vertex pick(k, 1:3) : Sensor coil position : コイル位置, Qpick(k, 1:3) : Sensor coil direction : コイル方向 Wsensor(m,n) = n-th coil weight for m-th sensor channel basis(channel,dipole) = Wsensor * basis(coil,dipole) bem_mode: Leadfield calculation mode MEG 0 : MEG (Sphere model:Sarvas) 球モデル 1 : MEG (BEM) 境界要素法 2 : MEG (Spherical Harmonics) 球面関数展開 EEG 3 : EEG (Sphere model) 3層球モデル 4 : EEG (BEM) 境界要素法 BEM : BEM parameters Vhead : Vertex on boudary surface 境界三角面の座標 Fhead : triangle of boudary surface 境界三角面頂点インデックス Omega : Solid angle matrix 境界面係数行列 Sout : Normal of boudary surface 三角面の外向き面積法線 SPHinfo : Parameters for Spherical Harmonics --- OUTPUT basis( L * Nvertex , Npick) : Magnetic field for MEG case : Electric potential for EEG case basis( n, k ) : k-th sensor magnetic field for dipole current at n : k-th sensor electric potential for dipole current at n written by M. Sato 2005-8-8 modified by T. Sako 2006-8-23 Copyright (C) 2011, ATR All Rights Reserved. License : New BSD License(see VBMEG_LICENSE.txt)
0001 function basis = vb_calc_lead_field(V, xx, pick, Qpick, Wsensor, bem_mode, ... 0002 BEM, Vhead, Fhead, Omega, Sout, SPHinfo); 0003 % Leadfield calculation 0004 % 0005 % --- INPUT 0006 % V(n,:) : dipole position (3-D coordinate) at n-th vertex 0007 % xx(n,:) : dipole current direction (unit vector) at n-th vertex 0008 % pick(k, 1:3) : Sensor coil position : コイル位置, 0009 % Qpick(k, 1:3) : Sensor coil direction : コイル方向 0010 % 0011 % Wsensor(m,n) = n-th coil weight for m-th sensor channel 0012 % basis(channel,dipole) = Wsensor * basis(coil,dipole) 0013 % 0014 % bem_mode: Leadfield calculation mode 0015 % MEG 0016 % 0 : MEG (Sphere model:Sarvas) 球モデル 0017 % 1 : MEG (BEM) 境界要素法 0018 % 2 : MEG (Spherical Harmonics) 球面関数展開 0019 % EEG 0020 % 3 : EEG (Sphere model) 3層球モデル 0021 % 4 : EEG (BEM) 境界要素法 0022 % 0023 % BEM : BEM parameters 0024 % Vhead : Vertex on boudary surface 境界三角面の座標 0025 % Fhead : triangle of boudary surface 境界三角面頂点インデックス 0026 % Omega : Solid angle matrix 境界面係数行列 0027 % Sout : Normal of boudary surface 三角面の外向き面積法線 0028 % SPHinfo : Parameters for Spherical Harmonics 0029 % 0030 % --- OUTPUT 0031 % basis( L * Nvertex , Npick) 0032 % : Magnetic field for MEG case 0033 % : Electric potential for EEG case 0034 % basis( n, k ) 0035 % : k-th sensor magnetic field for dipole current at n 0036 % : k-th sensor electric potential for dipole current at n 0037 % 0038 % written by M. Sato 2005-8-8 0039 % modified by T. Sako 2006-8-23 0040 % 0041 % Copyright (C) 2011, ATR All Rights Reserved. 0042 % License : New BSD License(see VBMEG_LICENSE.txt) 0043 0044 % Basis mode definition 0045 global vbmeg_inst; 0046 const = vbmeg_inst.const; 0047 % const.BASIS_MEG_SPHERE; % MEG Sphere model (Sarvas) 0048 % const.BASIS_MEG_BEM; % MEG BEM 0049 % const.BASIS_MEG_HARMONICS;% MEG Spherical harmonics expansion 0050 % const.BASIS_EEG_SPHERE; % EEG 3-shell Sphere model 0051 % const.BASIS_EEG_BEM; % EEG 3-shell BEM 0052 0053 NV = size(V,1); 0054 Npick = size(pick,1); 0055 Nsensor = size(Wsensor,1); 0056 0057 % 0058 % --- Leadfield calculation (リードフィールド行列計算) 0059 % 0060 tic; 0061 0062 switch bem_mode 0063 case const.BASIS_MEG_SPHERE, 0064 % MEG Sphere model (Sarvas) 0065 fprintf('--- MEG Sphere model (Sarvas) '); 0066 0067 BB = zeros(NV , Npick); 0068 0069 for i=1:Npick, 0070 BB(:,i) = vb_sarvas_new(V, xx, pick(i,:), Qpick(i,:) ); 0071 end 0072 0073 %basis = BB(:,1:Nsensor) - BB(:,(Nsensor+1):end); 0074 basis = BB * Wsensor'; 0075 0076 case const.BASIS_MEG_BEM, 0077 % MEG BEM 0078 fprintf('--- MEG BEM-linear-Galerkin '); 0079 0080 BB = vb_bem_meg(Omega,Sout,Vhead,Fhead,V,xx,pick,Qpick); 0081 0082 %basis = BB(1:Nsensor,:) - BB((Nsensor+1):end,:); 0083 basis = Wsensor * BB; 0084 basis = basis'; 0085 0086 case const.BASIS_MEG_HARMONICS, 0087 % MEG Spherical harmonics expansion 0088 fprintf('--- MEG Spherical Harmonics expansion '); 0089 N = SPHinfo.func_order; 0090 Rmax = SPHinfo.radius; 0091 A = SPHinfo.harmo_coef; 0092 0093 BB = vb_spherical_harmo_magnetic_field(V, xx, N, Rmax, pick, Qpick, A); 0094 0095 %basis = BB(:,1:Nsensor) - BB(:,(Nsensor+1):end); 0096 basis = BB * Wsensor'; 0097 0098 case const.BASIS_EEG_SPHERE, 0099 % EEG 3-shell Sphere model 0100 fprintf('--- EEG potential (Sphere model) '); 0101 0102 Rmax = max(sqrt(sum(V.^2,2))); 0103 Rpick = mean(sqrt(sum(pick.^2,2))); 0104 R = BEM.R * Rpick; 0105 0106 % Electric potential 0107 basis = vb_eeg_legendre(V, xx, pick, R, BEM.sigma); 0108 0109 case const.BASIS_EEG_BEM, 0110 % EEG 3-shell BEM 0111 fprintf('--- EEG potential (BEM) '); 0112 0113 % Electric potential 0114 basis = vb_bem_eeg(Omega,Sout,Vhead,Fhead,V,xx,pick,BEM); 0115 0116 basis = basis'; 0117 end 0118 0119 fprintf('%f[sec]\n',toc);