[fia,theta]=meshgrid([linspace(0,pi,100),pi]);
x=sin(theta).*cos(fia);
y=sin(theta).*sin(fia);
z=cos(theta);
surf(x,y,z)
shading interp
axis equal