
n = 20;   % plots an n by n grid of quadrilateral cells from certain slices
          % of 3d grid.  

% Plot images of 6 sides of computational cube - shows full sphere

clf

[x1,x2] = meshgrid(linspace(-1,1,n+1), linspace(-1,1,n+1));
one = ones(size(x1));

xc = x1;
yc = x2;
zc = one;
[xp,yp,zp] = mapc2p(xc,yc,zc);
surf(xp,yp,zp)

hold on
zc = -one;
[xp,yp,zp] = mapc2p(xc,yc,zc);
surf(xp,yp,zp)

xc = x1;
yc = one;
zc = x2;
[xp,yp,zp] = mapc2p(xc,yc,zc);
surf(xp,yp,zp)

yc = -one;
[xp,yp,zp] = mapc2p(xc,yc,zc);
surf(xp,yp,zp)

xc = one;
yc = x1;
zc = x2;
[xp,yp,zp] = mapc2p(xc,yc,zc);
surf(xp,yp,zp)

xc = -one;
[xp,yp,zp] = mapc2p(xc,yc,zc);
surf(xp,yp,zp)


colormap([1 1 1]);   
daspect([1 1 1])
axis([-1 1 -1 1 -1 1])
hold off

