function [xp,yp,zp] = mapc2m(xc,yc);
%
% Specifies the mapping from computational coordinates to a 2d manifold.
% Maps the rectangle [-3,1] x [-1,1] to the sphere of radius Rsphere.
%  [-1,1]x[-1,1]  is mapped to the upper hemisphere.
%  [-3,-1]x[-1,1] is mapped to the lower hemisphere.
%
%  Wave Propagation Software, Computational Science, and Reproducible Research
%  by R. J. LeVeque, http://www.amath.washington.edu/~rjl/pubs/icm06

Rsphere = 1;
r1 = Rsphere;
sgnz = ones(size(xc));
ij = find(xc < -1);
xc(ij) = -2 - xc(ij);
sgnz(ij) = -1;
xc1 = abs(xc);
yc1 = abs(yc);
d = max(xc1,yc1);
d = max(d, 1e-10);
D = r1*d.*(2-d) / sqrt(2);
R = r1*ones(size(d));

center = D - sqrt(R.^2 - D.^2);
xp = D./d .* xc1;
yp = D./d .* yc1;

ij = find(yc1==d);
yp(ij) = center(ij) + sqrt(R(ij).^2 - xp(ij).^2);
ij = find(xc1==d);
xp(ij) = center(ij) + sqrt(R(ij).^2 - yp(ij).^2);

xp = sign(xc) .* xp;
yp = sign(yc) .* yp;
zp = sgnz .* sqrt(r1^2 - (xp.^2 + yp.^2));