All Downloads are FREE. Search and download functionalities are using the official Maven repository.
Many resources are needed to download a project. Please understand that we have to compensate our server costs. Thank you in advance.
Project price only 1 $
You can buy this project and download/modify it how often you want.
org.locationtech.proj4j.proj.LambertAzimuthalEqualAreaProjection Maven / Gradle / Ivy
/*******************************************************************************
* Copyright 2009, 2017 Martin Davis
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
/*
* This file was semi-automatically converted from the public-domain USGS PROJ source.
*/
package org.locationtech.proj4j.proj;
import org.locationtech.proj4j.ProjCoordinate;
import org.locationtech.proj4j.ProjectionException;
import org.locationtech.proj4j.util.ProjectionMath;
public class LambertAzimuthalEqualAreaProjection extends Projection {
private static final int N_POLE = 0;
private static final int S_POLE = 1;
private static final int EQUIT = 2;
private static final int OBLIQ = 3;
private int mode = 0;
private double phi0;
private double sinb1;
private double cosb1;
private double xmf;
private double ymf;
private double mmf;
private double qp;
private double dd;
private double rq;
private double[] apa;
private double sinph0;
private double cosph0;
public LambertAzimuthalEqualAreaProjection() {
this( false );
}
public LambertAzimuthalEqualAreaProjection( boolean south ) {
//minLatitude = Math.toRadians(0);
//maxLatitude = Math.toRadians(90);
//projectionLatitude1 = south ? -ProjectionMath.QUARTERPI : ProjectionMath.QUARTERPI;
//projectionLatitude2 = south ? -ProjectionMath.HALFPI : ProjectionMath.HALFPI;
//initialize();
}
public void initialize() {
super.initialize();
double t;
phi0 = projectionLatitude;
if (Math.abs((t = Math.abs(phi0)) - ProjectionMath.HALFPI) < ProjectionMath.EPS10) {
mode = phi0 < 0. ? S_POLE : N_POLE;
}
else if (Math.abs(t) < ProjectionMath.EPS10) {
mode = EQUIT;
}
else {
mode = OBLIQ;
}
if (!spherical) {
double sinphi;
e = Math.sqrt(es);
qp = ProjectionMath.qsfn(1., e, one_es);
mmf = .5 / (1. - es);
apa = ProjectionMath.authset(es);
switch (mode) {
case N_POLE:
case S_POLE:
dd = 1.;
break;
case EQUIT:
dd = 1. / (rq = Math.sqrt(.5 * qp));
xmf = 1.;
ymf = .5 * qp;
break;
case OBLIQ:
rq = Math.sqrt(.5 * qp);
sinphi = Math.sin(phi0);
sinb1 = ProjectionMath.qsfn(sinphi, e, one_es) / qp;
cosb1 = Math.sqrt(1. - sinb1 * sinb1);
dd = Math.cos(phi0) / (Math.sqrt(1. - es * sinphi * sinphi) *
rq * cosb1);
ymf = (xmf = rq) / dd;
xmf *= dd;
break;
}
} else {
if (mode == OBLIQ) {
sinph0 = Math.sin(phi0);
cosph0 = Math.cos(phi0);
}
}
}
public ProjCoordinate project(double lplam, double lpphi, ProjCoordinate out) {
if (spherical) {
project_s(lplam, lpphi, out);
}
else {
project_e(lplam, lpphi, out);
}
return out;
}
public void project_s(double lplam, double lpphi, ProjCoordinate out) {
double coslam, cosphi, sinphi;
sinphi = Math.sin(lpphi);
cosphi = Math.cos(lpphi);
coslam = Math.cos(lplam);
switch (mode) {
case EQUIT:
case OBLIQ:
if (mode == EQUIT)
out.y = 1. + cosphi * coslam;
else
out.y = 1. + sinph0 * sinphi + cosph0 * cosphi * coslam;
if (out.y <= EPS10) throw new ProjectionException("F");
out.x = (out.y = Math.sqrt(2. / out.y)) * cosphi * Math.sin(lplam);
out.y *= mode == EQUIT ? sinphi :
cosph0 * sinphi - sinph0 * cosphi * coslam;
break;
case N_POLE:
coslam = -coslam;
case S_POLE:
if (Math.abs(lpphi + phi0) < EPS10) throw new ProjectionException("F");;
out.y = ProjectionMath.QUARTERPI - lpphi * .5;
out.y = 2. * (mode == S_POLE ? Math.cos(out.y) : Math.sin(out.y));
out.x = out.y * Math.sin(lplam);
out.y *= coslam;
break;
}
}
public void project_e(double lplam, double lpphi, ProjCoordinate out) {
double coslam, sinlam, sinphi, q, sinb=0.0, cosb=0.0, b=0.0;
coslam = Math.cos(lplam);
sinlam = Math.sin(lplam);
sinphi = Math.sin(lpphi);
q = ProjectionMath.qsfn(sinphi, e, one_es);
if (mode == OBLIQ || mode == EQUIT) {
sinb = q / qp;
cosb = Math.sqrt(1. - sinb * sinb);
}
switch (mode) {
case OBLIQ:
b = 1. + sinb1 * sinb + cosb1 * cosb * coslam;
break;
case EQUIT:
b = 1. + cosb * coslam;
break;
case N_POLE:
b = ProjectionMath.HALFPI + lpphi;
q = qp - q;
break;
case S_POLE:
b = lpphi - ProjectionMath.HALFPI;
q = qp + q;
break;
}
if (Math.abs(b) < EPS10) throw new ProjectionException("F");
switch (mode) {
case OBLIQ:
case EQUIT:
if (mode == OBLIQ) {
out.y = ymf * ( b = Math.sqrt(2. / b) )
* (cosb1 * sinb - sinb1 * cosb * coslam);
}
else {
out.y = (b = Math.sqrt(2. / (1. + cosb * coslam))) * sinb * ymf;
}
out.x = xmf * b * cosb * sinlam;
break;
case N_POLE:
case S_POLE:
if (q >= 0.) {
out.x = (b = Math.sqrt(q)) * sinlam;
out.y = coslam * (mode == S_POLE ? b : -b);
} else
out.x = out.y = 0.;
break;
}
}
public ProjCoordinate projectInverse(double xyx, double xyy, ProjCoordinate out) {
if (spherical) {
projectInverse_s(xyx, xyy, out);
}
else {
projectInverse_e(xyx, xyy, out);
}
return out;
}
public void projectInverse_s(double xyx, double xyy, ProjCoordinate out) {
double cosz=0.0, rh, sinz=0.0;
double lpphi, lplam;
rh = Math.hypot(xyx, xyy);
if ((lpphi = rh * .5 ) > 1.) throw new ProjectionException("I_ERROR");
lpphi = 2. * Math.asin(lpphi);
if (mode == OBLIQ || mode == EQUIT) {
sinz = Math.sin(lpphi);
cosz = Math.cos(lpphi);
}
switch (mode) {
case EQUIT:
lpphi = Math.abs(rh) <= EPS10 ? 0. : Math.asin(xyy * sinz / rh);
xyx *= sinz;
xyy = cosz * rh;
break;
case OBLIQ:
lpphi = Math.abs(rh) <= EPS10 ? phi0 :
Math.asin(cosz * sinph0 + xyy * sinz * cosph0 / rh);
xyx *= sinz * cosph0;
xyy = (cosz - Math.sin(lpphi) * sinph0) * rh;
break;
case N_POLE:
xyy = -xyy;
lpphi = ProjectionMath.HALFPI - lpphi;
break;
case S_POLE:
lpphi -= ProjectionMath.HALFPI;
break;
}
lplam = (xyy == 0. && (mode == EQUIT || mode == OBLIQ)) ?
0. : Math.atan2(xyx, xyy);
out.x = lplam;
out.y = lpphi;
}
public void projectInverse_e(double xyx, double xyy, ProjCoordinate out) {
double lpphi, lplam;
double cCe, sCe, q, rho, ab=0.0;
switch (mode) {
case EQUIT:
case OBLIQ:
if ((rho = Math.hypot(xyx /= dd, xyy *= dd)) < EPS10) {
lplam = 0.;
lpphi = phi0;
out.x = lplam;
out.y = lpphi;
return;
}
cCe = Math.cos(sCe = 2. * Math.asin(.5 * rho / rq));
xyx *= (sCe = Math.sin(sCe));
if (mode == OBLIQ) {
q = qp * (ab = cCe * sinb1 + xyy * sCe * cosb1 / rho);
xyy = rho * cosb1 * cCe - xyy * sinb1 * sCe;
} else {
q = qp * (ab = xyy * sCe / rho);
xyy = rho * cCe;
}
break;
case N_POLE:
xyy = -xyy;
case S_POLE:
if (0 == (q = (xyx * xyx + xyy * xyy)) ) {
lplam = 0.;
lpphi = phi0;
out.x = lplam;
out.y = lpphi;
return;
}
/*
q = P->qp - q;
*/
ab = 1. - q / qp;
if (mode == S_POLE)
ab = - ab;
break;
}
lplam = Math.atan2(xyx, xyy);
lpphi = ProjectionMath.authlat(Math.asin(ab), apa);
out.x = lplam;
out.y = lpphi;
}
/**
* Returns true if this projection is equal area
*/
public boolean isEqualArea() {
return true;
}
public boolean hasInverse() {
return true;
}
public String toString() {
return "Lambert Azimuthal Equal Area";
}
}
