function sSpline2Ply(sM, filename)
%  sSpline2Ply - writes B-spline spherical surface into a ascii PLY file
%
%      sSpline2Ply(sM, filename)
%
%      sM - matrix with B-spline knots, i.e. returned by sFit function
%      filename - a string with a filename
%
%  see also: sFit, sFitCart, sSuperQuad, sSurfM
%

sM = sExtend(sValueM(sM),[0 1],[0 0]);
[m,n] = size(sM);
fi = (((1:m)-1)/(m-1)*2*pi-pi)' * ones(1,n);
teta= ones(m,1) * (((1:n)-1)/(n-1)*pi-pi/2);

x = sM.*cos(teta).*cos(fi); 
y = sM.*cos(teta).*sin(fi); 
z = sM.*sin(teta);

points = [ x(:) y(:) z(:)];

t = 1:(m-1)*n;
upper = t' * ones(1,3) + ones((m-1)*n,1)*[1 0 m];
lower = t' * ones(1,3) + ones((m-1)*n,1)*[1 m (m+1)];
triangles = ones(2* (m-1)*n,5);
triangles( 2*t-1, 2:4 ) = upper;
triangles( 2*t  , 2:4 ) = lower;

clear upper lower
triangles( :,1) = 3*ones(size(triangles,1),1);

FID = fopen(filename,'wt');
fprintf( FID, 'ply\n');
fprintf( FID, 'format ascii 1.0\n');
fprintf( FID, 'element vertex %d\n', m*n);
fprintf( FID, 'property float x\n');
fprintf( FID, 'property float y\n');
fprintf( FID, 'property float z\n');
fprintf( FID, 'element face %d\n', 2*(m-1)*(n-1));
fprintf( FID, 'property list uchar int vertex_indices\n');
fprintf( FID, 'property int surf_type\n');
fprintf( FID, 'end_header\n');

fprintf( FID, '%10.5f %10.5f %10.5f \n', points');
fprintf( FID, '%10d   %10d   %10d   %10d   %10d \n', triangles');
fclose(FID);

return;