public ProjCoordinate project(double lplam, double lpphi, ProjCoordinate out) {
double al, al2, g, g2, p2;
p2 = Math.abs(lpphi / ProjectionMath.HALFPI);
if ((p2 - TOL) > 1.) throw new ProjectionException("F");
if (p2 > 1.)
p2 = 1.;
if (Math.abs(lpphi) <= TOL) {
out.x = lplam;
out.y = 0.;
} else if (Math.abs(lplam) <= TOL || Math.abs(p2 - 1.) < TOL) {
out.x = 0.;
out.y = Math.PI * Math.tan(.5 * Math.asin(p2));
if (lpphi < 0.) out.y = -out.y;
} else {
al = .5 * Math.abs(Math.PI / lplam - lplam / Math.PI);
al2 = al * al;
g = Math.sqrt(1. - p2 * p2);
g = g / (p2 + g - 1.);
g2 = g * g;
p2 = g * (2. / p2 - 1.);
p2 = p2 * p2;
out.x = g - p2; g = p2 + al2;
out.x = Math.PI * (al * out.x + Math.sqrt(al2 * out.x * out.x - g * (g2 - p2))) / g;
if (lplam < 0.) out.x = -out.x;
out.y = Math.abs(out.x / Math.PI);
out.y = 1. - out.y * (out.y + 2. * al);
if (out.y < -TOL) throw new ProjectionException("F");
if (out.y < 0.)
out.y = 0.;
else
out.y = Math.sqrt(out.y) * (lpphi < 0. ? -Math.PI : Math.PI);
}