Scilabにおける3D plotter表示

Scilabを用いた3D plotter表示プログラムの説明です。 // preparation MEG data
stacksize(50000000);
loadmatfile('.../xxx.mat');
getf('.../CalcSSP2.sce');
getf('.../FiltfiltIIR.sce');
meg=CalcSSP2(ave_0001(:,1:306));
time=ave_0001(:,$);
smp=length(time);
meg=FiltfiltIIR(meg,time,4,'bp','ellip',[0.5,100],[0.1,0.1]);
// preparation sensor array
P=MEG2HEAD(1:3,1:3)*ChInform(1:3,1:3:304)+MEG2HEAD(1:3,4)*ones(1,102);
X=MEG2HEAD(1:3,1:3)*ChInform(4:6,1:3:304);
Y=MEG2HEAD(1:3,1:3)*ChInform(7:9,1:3:304);
Z=MEG2HEAD(1:3,1:3)*ChInform(10:12,1:3:304);
wd=0.014;
V=zeros(3,102,5);
V(:,:,1)=P+wd*X+wd*Y;
V(:,:,2)=P+wd*X-wd*Y;
V(:,:,3)=P-wd*X-wd*Y;
V(:,:,4)=P-wd*X+wd*Y;
V(:,:,5)=V(:,:,1);
// draw process
scf();hf=gcf();hf.visible='off';
for t1=1:smp;if time(t1)>=-50;break;end;end;
for t2=t1:smp;if time(t2)>=150;break;end;end;
M=-5e-5*meg(t1:t2,:)';
M(3:3:306,:)=M(3:3:306,:)/5;
smpt=size(M,2);
W1=zeros(3,102,smpt);
W2=W1;
W3=W1;
CL=P-wd*Y;CR=P+wd*Y;
dX=7/1000*X;
for sns=1:102;...
Tx=linspace(CL(1,sns),CR(1,sns),smpt);...
Ty=linspace(CL(2,sns),CR(2,sns),smpt);...
Tz=linspace(CL(3,sns),CR(3,sns),smpt);...
W1(:,sns,:)=[Tx;Ty;Tz]+X(:,sns)*M(sns*3-2,:)+dX(:,sns)*ones(1,smpt);...
W2(:,sns,:)=[Tx;Ty;Tz]+X(:,sns)*M(sns*3-1,:)-dX(:,sns)*ones(1,smpt);...
W3(:,sns,:)=[Tx;Ty;Tz]+X(:,sns)*M(sns*3,:);
end;
for sns=1:102;...
param3d(matrix(V(1,sns,:),1,5),matrix(V(2,sns,:),1,5),matrix(V(3,sns,:),1,5));...
param3d(matrix(W1(1,sns,:),smpt,1),matrix(W1(2,sns,:),smpt,1),matrix(W1(3,sns,:),smpt,1));...
param3d(matrix(W2(1,sns,:),smpt,1),matrix(W2(2,sns,:),smpt,1),matrix(W2(3,sns,:),smpt,1));...
param3d(matrix(W3(1,sns,:),smpt,1),matrix(W3(2,sns,:),smpt,1),matrix(W3(3,sns,:),smpt,1));...
end; ha=gca();ha.box='off';ha.isoview='on';
ha.margins=[0,0,0,0];
ha.axes_visible='off';
ha.children(1:4:408).foreground=13;
ha.children(2:4:408).foreground=2;
ha.children(3:4:408).foreground=2;
ha.rotation_angles=[30,180];
hf.visible='on';


Scilabの3D陰線処理はbugが残っています。そこでセンサー裏面を表示させないことにします。
x=[30,180]/180*%pi; //予め角度を決めておく
yaxis=[cos(x(1)),0,sin(x(1));0,1,0;-sin(x(1)),0,cos(x(1))];
zaxis=[cos(x(2)),-sin(x(2)),0;sin(x(2)),cos(x(2)),0;0,0,1];
angle=zaxis*yaxis*[0;0;1];
visible=-sum(Z.*(angle*ones(1,102)),1);
scf();kf=gcf();kf.visible='off';
for sns=1:102;...
if visible(sns)<=0; param3d(matrix(V(1,sns,:),1,5),matrix(V(2,sns,:),1,5),matrix(V(3,sns,:),1,5));...
param3d(matrix(W1(1,sns,:),smpt,1),matrix(W1(2,sns,:),smpt,1),matrix(W1(3,sns,:),smpt,1));...
param3d(matrix(W2(1,sns,:),smpt,1),matrix(W2(2,sns,:),smpt,1),matrix(W2(3,sns,:),smpt,1));...
param3d(matrix(W3(1,sns,:),smpt,1),matrix(W3(2,sns,:),smpt,1),matrix(W3(3,sns,:),smpt,1));...
end;end; ka=gca();ka.box='off';ka.isoview='on';
ka.margins=[0,0,0,0];ka.axes_visible='off';
ka.children(1:4:$).foreground=13;
ka.children(2:4:$).foreground=2;
ka.children(3:4:$).foreground=2;
ka.rotation_angles=x/%pi*180;
kf.visible='on';

いろいろ試しましたが、はじめから裏面に該当するセンサーを計算させないようにするしか、いい方法はないようです。 MATLABを使えばこんな苦労はありません。ScilabではopenGLやDirectXを使ってないのかもしれません。