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 = 境インデ% = [ start_id(1) end_id(1) ; 0045 % ... ; 0046 % start_id(Nsurf) end_id(Nsurf) ] 0047 % = [始ンデ 終了ンデ 0048 % innermost to outermost 0049 % 0050 % Npatch = 境角インデ 0051 % = [始ンデ 終了ンデ 0052 % Sigma = 0053 % = [sigma(1), ..., sigma(Nsurf)] 0054 % 0055 % 2008-10-9 Masa-aki Sato 0056 % structure of program is changed (algorithm is not changed) 0057 0058 % pre-processing 0059 [Basis_mode, pick, Qpick, Wsensor, V0, V, xx] = ... 0060 vb_lf_preprocess(proj_root, basis_parm); 0061 0062 % check parameters and set default values if not be set 0063 basis_parm = vb_util_check_head_shell_info(proj_root, basis_parm); 0064 0065 [V,xx] = vb_current_vector(V, xx, Basis_mode); 0066 0067 if iscell(basis_parm.head_file) 0068 headfile = basis_parm.head_file; 0069 joinfile = [vb_join_cell_file_name(headfile) '.head.mat']; 0070 joinpath = fullfile(proj_root, joinfile); 0071 0072 if ~exist(joinpath,'file') 0073 sigma=[];layer=[];mode = 3; % MAX mode for radius 0074 vb_head_join_files(headfile,proj_root,joinpath, sigma,layer,mode); 0075 end 0076 basis_parm.head_file = joinfile; 0077 end 0078 0079 new_headfile = fullfile(proj_root, basis_parm.head_file); 0080 omega_file = [vb_get_basename(new_headfile,'.head.mat') '.omega.mat']; 0081 0082 fprintf('Load %s\n',new_headfile) 0083 load(new_headfile); 0084 % Vhead : (# of vertex) x 3 0085 % Fhead : (# of patch) x 3 0086 % Nvertex : (# of surface) x 2 0087 % Npatch : (# of surface) x 2 0088 % Sigma : 1 x (# of surface) 0089 0090 fprintf('--- BEM for EEG\n') 0091 fprintf('# of vertices = %d\n',size(Vhead,1)); 0092 0093 if ~exist('Nvertex', 'var'), Nvertex = [1 length(Vhead)]; end; 0094 if ~exist('Npatch', 'var'), Npatch = [1 length(Fhead)]; end; 0095 0096 Nsurf = size(Nvertex,1); 0097 fprintf('# of surfaces = %d\n',Nsurf); 0098 0099 Sigma = basis_parm.sigma; 0100 0101 BEM.sigma = [Sigma(1:Nsurf) 0]; 0102 BEM.Nvertex = Nvertex; 0103 BEM.Npatch = Npatch; 0104 0105 if exist(omega_file,'file') 0106 fprintf('Load %s\n',omega_file) 0107 load(omega_file) 0108 end 0109 0110 if ~exist('Omega','var') || basis_parm.Recalc == ON, 0111 tic; 0112 fprintf('--- Solid angle calculation \n'); 0113 0114 % (Normal vectors (XXhead) are outwarded here) 0115 [Omega, Sout] = vb_solid_angle_grk(Vhead,Fhead,XXhead); 0116 0117 vb_ptime(toc); 0118 vb_fsave(omega_file,'Omega', 'Sout'); 0119 fprintf('Saved omega file: %s\n',omega_file); 0120 end 0121 0122 tic; 0123 0124 % 多層モデfprintf('--- Multisurface correction\n'); 0125 [Omega] = vb_bem_matrix(Omega,BEM); 0126 0127 % 数を三角寄 0128 fprintf('--- Auto solid angle calculation\n'); 0129 Omega = vb_solid_auto_grk(Omega,Fhead,Sout); 0130 0131 vb_ptime(toc); 0132 0133 fprintf('--- EEG potential (BEM) \n'); 0134 0135 % Electric potential 0136 basis = vb_bem_eeg( Omega, Sout, Vhead, Fhead, V, xx, pick, BEM ); 0137 basis = basis'; 0138 0139 %%% END OF FILE %%%