// Two examples using contour3d
// Tony Roberts, July 2006
getf("contour3d.sci");

// first draw spheres centred on (1,0,0)
nx=3;ny=4;nz=5;
x=linspace(0,1,nx); 
y=linspace(0,1,ny); 
z=linspace(0,1,nz); 
// make row/column/?? of hypermat
xx=hypermat([nx,1,1],x);
yy=hypermat([1,ny,1],y);
zz=hypermat([1,1,nz],z);
xx=xx(:,ones(1,ny),ones(1,nz));
yy=yy(ones(1,nx),:,ones(1,nz));
zz=zz(ones(1,nx),ones(1,ny),:);
// compute function and draw
f=sqrt((xx-1).^2+yy.^2+zz.^2)-0.5;
xbasc()
xset("colormap",hotcolormap(6))
fs=contour3d(x,y,z,f,6)

// second, use cylindrical coordinates
nr=4;nt=5;nz=3;
r=linspace(0.25,1,nr); 
t=linspace(0,%pi,nt); 
z=linspace(0,1,nz); 
// make row/column/?? of hypermat
r=hypermat([nr,1,1],r);
t=hypermat([1,nt,1],t);
z=hypermat([1,1,nz],z);
r=r(:,ones(1,nt),ones(1,nz));
t=t(ones(1,nr),:,ones(1,nz));
z=z(ones(1,nr),ones(1,nt),:);
// compute function
f=r-z+t/2;
// get xy-coords of polars and draw
x=r.*cos(t);
y=r.*sin(t);
scf();
xset("colormap",jetcolormap(6))
fs=contour3d(x,y,z,f,6,34,47)