Adaptive beamformer3
電子技術総合研究所の浅野太氏のHP
(たぶん今は独立行政法人産業技術総合研究所)を参考にして、
とし、γをBBT の最大特異値の1%と最大固有値の1%で試してみました。
//preparation//
stacksize(50000000);
loadmatfile('.../xxx.mat');
loadmatfile('.../xxx-bem_pointsetXX.mat');
getf('.../CalcLFM.sce');
getf('.../CalcSSP2.sce');
Ch=ConversionCoordination(MEG2HEAD,ChInform);
Center=[0,0,40]/1000;
R=CalcTangentialVectors(pointset.points,Center);
wm=1;// gradiometer : magnetometer
LFM=CalcVectorViewLFM(R,Ch,Center,wm);
//MEG data, noise data, variance//
meg=CalcSSP2(ave_0001(:,1:306));
meg=meg';time=ave_0001(:,$);
Cmeg=meg*meg';
//adaptive beamformer (minimum variance)//
[Usvd,Ssvd,Vsvd]=svd(Cmeg);
Tikhonov=0.01*Ssvd(1,1)*eye(Usvd);// 1% of max singular value
//Tikhonov=0.01*Ssvd(1,1)^2*eye(Usvd);// 1% of max eigen value
CS=inv(Cmeg+Tikhonov);
H=zeros(LFM);
for x=1:size(LFM,2);H(:,x)=CS*LFM(:,x)/(LFM(:,x)'*CS*LFM(:,x));end;
Q=H'*meg;
//Graph//
scf();
set(gcf(),'visible','off');xset('use color',0);
subplot(211);plot2d(time,meg',rect=[-50,-600,500,600]);
subplot(212);plot2d(time,Q',rect=[-50,-5E-9,500,5E-9]);
set(gcf(),'visible','on');
なんかいい感じです。
因みにγを最大固有値の1%としたものです。
格子点の電流が30nA以上になっています。
-50〜500ミリ秒のGIF画像を作成し、動画化してみました。
color scaleは0〜4nAです。
Q=Q*1E+9;
scale=10;
scf();
ncolor=32;P=pointset.points(:,1:3)';nq=size(P,2);
set(gcf(),'color_map',jetcolormap(ncolor));
for t=1:length(time);...
set(gcf(),'visible','off');
q=(Q(1:nq,t).^2+Q(nq+1:$,t).^2).^0.5;
subplot(121);...
param3d(P(1,:),P(2,:),P(3,:));...
h=gca();h.box='off';...
h.children.mark_mode='on';h.children.mark_size=1;...
h.children.foreground=1;...
for nn=1:ncolor-1;...
PP=P;n=1:nq;n(q < nn/ncolor*scale)=[];PP=PP(:,n);...
if size(PP,2)==0;break;end;
param3d(PP(1,:),PP(2,:),PP(3,:));...
h.children(1).mark_mode='on';...
h.children(1).mark_size=round(nn/32*5+3);...
h.children(1).foreground=nn+1;...
end;
h.rotation_angles=[70,200];...
h.axes_visible='off';h.isoview='on';...
subplot(122);Matplot((ncolor:-1:1)');...
hh=gca();hh.axes_visible='off';hh.box='off';...
hh.tight_limits='on';hh.axes_bounds=[0.75,0.2,0.05,0.6];...
h.axes_bounds=[0,0,0.9,1];h.margins=[0,0,0,0];...
xstring(-6,-3,sprintf('%0.1fmsec',time(t)));...
set(gce(),'font_size',3);...
set(gcf(),'visible','on');...
xs2gif(get(gcf(),'figure_id'),sprintf('AB_%d.gif',t));...
clf();...
end;
GIF動画(1.6MB) MPEG1動画(2.5MB)
最大固有値の1%の時の動画です。color scaleは30nAです。
GIF動画(1.8MB) MPEG1動画(2.2MB)
点の大きさを電流の大きさに反映させているはずなのですが、
xs2gifでGIF画像として描き出すときに無視されるようです。
bugです。
L2ノルムでは格子点が格子点付近の電流活動を代表しますが、adaptive beamformerではお構いなしです。
すなわち格子点間にもっと大きな電流が存在するかもしれません!
格子点を5mm間隔としました。グラフは処理速度が余りにも遅いので描いていません。
γを最大特異値の1%としたものです。color scaleは0〜4nAです。
GIF動画(6.2MB) MPEG1動画(5.2MB)
γを最大固有値の1%としたものです。color scaleは0〜40nAです。
GIF動画(4.4MB) MPEG1動画(4.2MB)
なんとなくbeamformerっぽい雰囲気が出ていると思います。
しかし、adaptive beamformerって、parameterの設定でどうにでもなるような気がします。