Add utm point representation

This commit is contained in:
Victor Shcherb 2014-12-06 17:14:01 +01:00
parent 9f1a9c0f3b
commit 02233788fd
9 changed files with 3373 additions and 0 deletions

View file

@ -82,6 +82,7 @@
- junidecode-0.1.jar - BSD-4-Clause-UC (http://sourceforge.net/projects/junidecode/) - junidecode-0.1.jar - BSD-4-Clause-UC (http://sourceforge.net/projects/junidecode/)
- kxml2-2.3.0.jar - BSD license (http://www.kxml.org/) - kxml2-2.3.0.jar - BSD license (http://www.kxml.org/)
- tuprolog.jar - LGPL (http://apice.unibo.it/xwiki/bin/view/Tuprolog/) - tuprolog.jar - LGPL (http://apice.unibo.it/xwiki/bin/view/Tuprolog/)
- OpenMap framework - Apache License (https://code.google.com/p/android-openmap-framework/)
* Pull-requests and translations: * Pull-requests and translations:
- All pull requests are accepted under MIT License (most honorable contributors are mentioned in AUTHORS list) - All pull requests are accepted under MIT License (most honorable contributors are mentioned in AUTHORS list)

View file

@ -0,0 +1,167 @@
// **********************************************************************
//
// <copyright>
//
// BBN Technologies
// 10 Moulton Street
// Cambridge, MA 02138
// (617) 873-8000
//
// Copyright (C) BBNT Solutions LLC. All rights reserved.
//
// </copyright>
// **********************************************************************
package com.jwetherell.openmap.common;
public class Ellipsoid {
/** "Airy" */
public final static Ellipsoid AIRY = new Ellipsoid("Airy", 6377563.0d, 0.00667054d);
/** "Australian National" */
public final static Ellipsoid AUSTRALIAN_NATIONAL = new Ellipsoid("Australian National", 6378160.0d, 0.006694542d);
/** "Bessel 1841" */
public final static Ellipsoid BESSEL_1841 = new Ellipsoid("Bessel 1841", 6377397.0d, 0.006674372d);
/** "Bessel 1841 (Nambia) " */
public final static Ellipsoid BESSEL_1841_NAMIBIA = new Ellipsoid("Bessel 1841 Namibia", 6377484.0d, 0.006674372d);
/** "Clarke 1866" */
public final static Ellipsoid CLARKE_1866 = new Ellipsoid("Clarke 1866", 6378206.0d, 0.006768658d);
/** "Clarke 1880" */
public final static Ellipsoid CLARKE_1880 = new Ellipsoid("Clarke 1880", 6378249.0d, 0.006803511d);
/** "Everest" */
public final static Ellipsoid EVEREST = new Ellipsoid("Everest", 6377276.0d, 0.006637847d);
/** "Fischer 1960 (Mercury) " */
public final static Ellipsoid FISHER_1960_MERCURY = new Ellipsoid("Fisher 1960 Mercury", 6378166.0d, 0.006693422d);
/** "Fischer 1968" */
public final static Ellipsoid FISHER_1968 = new Ellipsoid("Fisher 1968", 6378150.0d, 0.006693422d);
/** "GRS 1967" */
public final static Ellipsoid GRS_1967 = new Ellipsoid("GRS 1967", 6378160.0d, 0.006694605d);
/** "GRS 1980" */
public final static Ellipsoid GRS_1980 = new Ellipsoid("GRS 1980", 6378137.0d, 0.081819191d, 0.00669438d, 6356752.3141d);
/** "Helmert 1906" */
public final static Ellipsoid HELMERT_1906 = new Ellipsoid("Helmert 1906", 6378200.0d, 0.006693422d);
/** "Hough" */
public final static Ellipsoid HOUGH = new Ellipsoid("Hough", 6378270.0d, 0.00672267d);
/** "International" */
public final static Ellipsoid INTERNATIONAL = new Ellipsoid("International", 6378388.0d, 0.08199189, 0.00672267d, 6356911.946d);
/** "Krassovsky" */
public final static Ellipsoid KRASSOVSKY = new Ellipsoid("Krassovsky", 6378245.0d, 0.006693422d);
/** "Modified Airy" */
public final static Ellipsoid MODIFIED_AIRY = new Ellipsoid("Modified Airy", 6377340.0d, 0.00667054d);
/** "Modified Everest" */
public final static Ellipsoid MODIFIED_EVEREST = new Ellipsoid("Modified Everest", 6377304.0d, 0.006637847d);
/** "Modified Fischer 1960" */
public final static Ellipsoid MODIFIED_FISCHER_1960 = new Ellipsoid("Modified Fischer", 6378155.0d, 0.006693422d);
/** "South American 1969" */
public final static Ellipsoid SOUTH_AMERICAN_1969 = new Ellipsoid("South American 1969", 6378160.0d, 0.006694542d);
/** "WGS 60" */
public final static Ellipsoid WGS_60 = new Ellipsoid("WSG 60", 6378165.0d, 0.006693422d);
/** "WGS 66" */
public final static Ellipsoid WGS_66 = new Ellipsoid("WGS 66", 6378145.0d, 0.006694542d);
/** "WGS-72" */
public final static Ellipsoid WGS_72 = new Ellipsoid("WGS 72", 6378135.0d, 0.006694318d);
/** "WGS-84" */
public final static Ellipsoid WGS_84 = new Ellipsoid("WGS 84", 6378137.0d, 0.081819191d, 0.00669438d, 6356752.3142d);
/**
* The display name for this ellipsoid.
*/
public final String name;
/**
* The equitorial radius for this ellipsoid.
*/
public final double radius;
/**
* The polar radius for this ellipsoid.
*/
public final double polarRadius;
/**
* The ellipsoid's eccentricity.
*/
public final double ecc;
/**
* The square of this ellipsoid's eccentricity.
*/
public final double eccsq;
/**
* Constructs a new Ellipsoid instance.
*
* @param radius
* The earth radius for this ellipsoid.
* @param eccsq
* The square of the eccentricity for this ellipsoid.
*/
public Ellipsoid(String name, double radius, double eccsq) {
this(name, radius, eccsq, Double.NaN, Double.NaN);
}
/**
* Constructs a new Ellipsoid instance.
*
* @param name
* The name of the ellipsoid.
* @param equitorialRadius
* The earth equitorial radius for this ellipsoid.
* @param ecc
* The eccentricity for this ellipsoid.
* @param eccsq
* The square of the eccentricity for this ellipsoid.
* @param polarRadius
* The earth polar radius for this ellipsoid.
*/
public Ellipsoid(String name, double equitorialRadius, double ecc, double eccsq, double polarRadius) {
this.name = name;
this.radius = equitorialRadius;
this.ecc = ecc;
this.eccsq = eccsq;
this.polarRadius = polarRadius;
}
/**
* Returns an array of all available ellipsoids in alphabetical order by
* name.
*
* @return An Ellipsoid[] array containing all the available ellipsoids
*/
public static Ellipsoid[] getAllEllipsoids() {
Ellipsoid[] all = { AIRY, AUSTRALIAN_NATIONAL, BESSEL_1841, BESSEL_1841_NAMIBIA, CLARKE_1866, CLARKE_1880, EVEREST, FISHER_1960_MERCURY, FISHER_1968,
GRS_1967, GRS_1980, HELMERT_1906, HOUGH, INTERNATIONAL, KRASSOVSKY, MODIFIED_AIRY, MODIFIED_EVEREST, MODIFIED_FISCHER_1960,
SOUTH_AMERICAN_1969, WGS_60, WGS_66, WGS_72, WGS_84 };
return all;
}
/**
* Given the name of an Ellipsoid, find the object for it out of the
* possible selections. Returns null if the Ellipsoid isn't found.
*
* @param name
* @return Ellipsoid
*/
public static Ellipsoid getByName(String name) {
Ellipsoid[] all = getAllEllipsoids();
if (name != null && name.length() > 0) {
name = name.replace('_', ' ');
for (int i = 0; i < all.length; i++) {
if (name.equalsIgnoreCase(all[i].name)) {
return all[i];
}
}
}
return null;
}
/**
* {@inheritDoc}
*/
@Override
public String toString() {
return "Ellipsoid[name=" + name + ", radius=" + radius + ", eccsq=" + eccsq + "]";
}
}

View file

@ -0,0 +1,640 @@
// **********************************************************************
//
// <copyright>
//
// BBN Technologies
// 10 Moulton Street
// Cambridge, MA 02138
// (617) 873-8000
//
// Copyright (C) BBNT Solutions LLC. All rights reserved.
//
// </copyright>
// **********************************************************************
package com.jwetherell.openmap.common;
public abstract class GreatCircle {
/**
* Calculate spherical arc distance between two points.
* <p>
* Computes arc distance `c' on the sphere. equation (5-3a). (0 &lt;= c
* &lt;= PI)
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param phi
* latitude in radians of end point
* @param lambda
* longitude in radians of end point
* @return float arc distance `c'
*
*/
public static final float sphericalDistance(float phi1, float lambda0, float phi, float lambda) {
return (float) sphericalDistance((double) phi1, (double) lambda0, (double) phi, (double) lambda);
}
/**
* Calculate spherical arc distance between two points with double
* precision.
* <p>
* Computes arc distance `c' on the sphere. equation (5-3a). (0 &lt;= c
* &lt;= PI)
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param phi
* latitude in radians of end point
* @param lambda
* longitude in radians of end point
* @return float arc distance `c'
*/
public static final double sphericalDistance(double phi1, double lambda0, double phi, double lambda) {
double pdiff = Math.sin(((phi - phi1) / 2.0));
double ldiff = Math.sin((lambda - lambda0) / 2.0);
double rval = Math.sqrt((pdiff * pdiff) + Math.cos(phi1) * Math.cos(phi) * (ldiff * ldiff));
return 2.0 * Math.asin(rval);
}
/**
* Calculate spherical azimuth between two points.
* <p>
* Computes the azimuth `Az' east of north from phi1, lambda0 bearing toward
* phi and lambda. (5-4b). (-PI &lt;= Az &lt;= PI).
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param phi
* latitude in radians of end point
* @param lambda
* longitude in radians of end point
* @return float azimuth east of north `Az'
*
*/
public static final float sphericalAzimuth(float phi1, float lambda0, float phi, float lambda) {
return (float) sphericalAzimuth((double) phi1, (double) lambda0, (double) phi, (double) lambda);
}
/**
* Calculate spherical azimuth between two points with double precision.
* <p>
* Computes the azimuth `Az' east of north from phi1, lambda0 bearing toward
* phi and lambda. (5-4b). (-PI &lt;= Az &lt;= PI).
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param phi
* latitude in radians of end point
* @param lambda
* longitude in radians of end point
* @return float azimuth east of north `Az'
*
*/
public static final double sphericalAzimuth(double phi1, double lambda0, double phi, double lambda) {
double ldiff = lambda - lambda0;
double cosphi = Math.cos(phi);
return Math.atan2(cosphi * Math.sin(ldiff), (Math.cos(phi1) * Math.sin(phi) - Math.sin(phi1) * cosphi * Math.cos(ldiff)));
}
/**
* Calculate point at azimuth and distance from another point, with double
* precision.
* <p>
* Returns a LatLonPoint.Double at arc distance `c' in direction `Az' from
* start point.
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param c
* arc radius in radians (0 &lt; c &lt;= PI)
* @param Az
* azimuth (direction) east of north (-PI &lt;= Az &lt; PI)
* @return LatLonPoint
*
*/
public static final LatLonPoint sphericalBetween(double phi1, double lambda0, double c, double Az) {
double cosphi1 = Math.cos(phi1);
double sinphi1 = Math.sin(phi1);
double cosAz = Math.cos(Az);
double sinAz = Math.sin(Az);
double sinc = Math.sin(c);
double cosc = Math.cos(c);
return new LatLonPoint(ProjMath.radToDeg(Math.asin(sinphi1 * cosc + cosphi1 * sinc * cosAz)), ProjMath.radToDeg(Math.atan2(sinc * sinAz, cosphi1 * cosc
- sinphi1 * sinc * cosAz)
+ lambda0));
}
/**
* Calculate point between two points.
* <p>
* Same as spherical_between() above except it calculates n equal segments
* along the length of c.
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param c
* arc radius in radians (0 &lt; c &lt;= PI)
* @param Az
* azimuth (direction) east of north (-PI &lt;= Az &lt; PI)
* @param n
* number of points along great circle edge to calculate
* @return float[n+1] radian lat,lon pairs
*
*/
public static final float[] sphericalBetween(float phi1, float lambda0, float c, float Az, int n) {
// full constants for the computation
double cosphi1 = Math.cos(phi1);
double sinphi1 = Math.sin(phi1);
double cosAz = Math.cos(Az);
double sinAz = Math.sin(Az);
int end = n << 1;
// new radian points
float[] points = new float[end + 2];
points[0] = phi1;
points[1] = lambda0;
float inc = c / n;
c = inc;
for (int i = 2; i <= end; i += 2, c += inc) {
// partial constants
double sinc = Math.sin(c);
double cosc = Math.cos(c);
// generate new point
points[i] = (float) Math.asin(sinphi1 * cosc + cosphi1 * sinc * cosAz);
points[i + 1] = (float) Math.atan2(sinc * sinAz, cosphi1 * cosc - sinphi1 * sinc * cosAz) + lambda0;
}
return points;
}
/**
* Calculate point between two points with double precision.
* <p>
* Same as spherical_between() above except it calculates n equal segments
* along the length of c.
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param c
* arc radius in radians (0 &lt; c &lt;= PI)
* @param Az
* azimuth (direction) east of north (-PI &lt;= Az &lt; PI)
* @param n
* number of points along great circle edge to calculate
* @return double[n+1] radian lat,lon pairs
*
*/
public static final double[] sphericalBetween(double phi1, double lambda0, double c, double Az, int n) {
// full constants for the computation
double cosphi1 = Math.cos(phi1);
double sinphi1 = Math.sin(phi1);
double cosAz = Math.cos(Az);
double sinAz = Math.sin(Az);
int end = n << 1;
// new radian points
double[] points = new double[end + 2];
points[0] = phi1;
points[1] = lambda0;
double inc = c / n;
c = inc;
for (int i = 2; i <= end; i += 2, c += inc) {
// partial constants
double sinc = Math.sin(c);
double cosc = Math.cos(c);
// generate new point
points[i] = Math.asin(sinphi1 * cosc + cosphi1 * sinc * cosAz);
points[i + 1] = Math.atan2(sinc * sinAz, cosphi1 * cosc - sinphi1 * sinc * cosAz) + lambda0;
}
return points;
}
/**
* Calculate great circle between two points on the sphere.
* <p>
* Folds all computation (distance, azimuth, points between) into one
* function for optimization. returns n or n+1 pairs of lat,lon on great
* circle between lat-lon pairs.
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param phi
* latitude in radians of end point
* @param lambda
* longitude in radians of end point
* @param n
* number of segments
* @param include_last
* return n or n+1 segments
* @return float[n] or float[n+1] radian lat,lon pairs
*
*/
public static final float[] greatCircle(float phi1, float lambda0, float phi, float lambda, int n, boolean include_last) {
// number of points to generate
int end = include_last ? n + 1 : n;
end <<= 1;// *2 for pairs
// calculate a bunch of stuff for later use
double cosphi = Math.cos(phi);
double cosphi1 = Math.cos(phi1);
double sinphi1 = Math.sin(phi1);
double ldiff = lambda - lambda0;
double p2diff = Math.sin(((phi - phi1) / 2));
double l2diff = Math.sin((ldiff) / 2);
// calculate spherical distance
double c = 2.0f * Math.asin(Math.sqrt(p2diff * p2diff + cosphi1 * cosphi * l2diff * l2diff));
// calculate spherical azimuth
double Az = Math.atan2(cosphi * Math.sin(ldiff), (cosphi1 * Math.sin(phi) - sinphi1 * cosphi * Math.cos(ldiff)));
double cosAz = Math.cos(Az);
double sinAz = Math.sin(Az);
// generate the great circle line
float[] points = new float[end];
points[0] = phi1;
points[1] = lambda0;
double inc = c / n;
c = inc;
for (int i = 2; i < end; i += 2, c += inc) {
// partial constants
double sinc = Math.sin(c);
double cosc = Math.cos(c);
// generate new point
points[i] = (float) Math.asin(sinphi1 * cosc + cosphi1 * sinc * cosAz);
points[i + 1] = (float) Math.atan2(sinc * sinAz, cosphi1 * cosc - sinphi1 * sinc * cosAz) + lambda0;
}
return points;
}
/**
* Calculate great circle between two points on the sphere with double
* precision.
* <p>
* Folds all computation (distance, azimuth, points between) into one
* function for optimization. returns n or n+1 pairs of lat,lon on great
* circle between lat-lon pairs.
* <p>
*
* @param phi1
* latitude in radians of start point
* @param lambda0
* longitude in radians of start point
* @param phi
* latitude in radians of end point
* @param lambda
* longitude in radians of end point
* @param n
* number of segments
* @param include_last
* return n or n+1 segments
* @return double[n] or double[n+1] radian lat,lon pairs
*
*/
public static final double[] greatCircle(double phi1, double lambda0, double phi, double lambda, int n, boolean include_last) {
// number of points to generate
int end = include_last ? n + 1 : n;
end <<= 1;// *2 for pairs
// calculate a bunch of stuff for later use
double cosphi = Math.cos(phi);
double cosphi1 = Math.cos(phi1);
double sinphi1 = Math.sin(phi1);
double ldiff = lambda - lambda0;
double p2diff = Math.sin(((phi - phi1) / 2));
double l2diff = Math.sin((ldiff) / 2);
// calculate spherical distance
double c = 2.0f * Math.asin(Math.sqrt(p2diff * p2diff + cosphi1 * cosphi * l2diff * l2diff));
// calculate spherical azimuth
double Az = Math.atan2(cosphi * Math.sin(ldiff), (cosphi1 * Math.sin(phi) - sinphi1 * cosphi * Math.cos(ldiff)));
double cosAz = Math.cos(Az);
double sinAz = Math.sin(Az);
// generate the great circle line
double[] points = new double[end];
points[0] = phi1;
points[1] = lambda0;
double inc = c / n;
c = inc;
for (int i = 2; i < end; i += 2, c += inc) {
// partial constants
double sinc = Math.sin(c);
double cosc = Math.cos(c);
// generate new point
points[i] = Math.asin(sinphi1 * cosc + cosphi1 * sinc * cosAz);
points[i + 1] = Math.atan2(sinc * sinAz, cosphi1 * cosc - sinphi1 * sinc * cosAz) + lambda0;
}
return points;
}
/**
* Calculate partial earth circle on the sphere.
* <p>
* Returns n float lat,lon pairs at arc distance c from point at
* phi1,lambda0.
* <p>
*
* @param phi1
* latitude in radians of center point
* @param lambda0
* longitude in radians of center point
* @param c
* arc radius in radians (0 &lt; c &lt; PI)
* @param s
* starting angle in radians. North up is zero.
* @param e
* angular extent in radians, clockwise right from starting
* angle.
* @param n
* number of points along circle edge to calculate
* @return float[n] radian lat,lon pairs along earth circle
*
*/
public static final float[] earthCircle(float phi1, float lambda0, float c, float s, float e, int n) {
return earthCircle(phi1, lambda0, c, s, e, n, new float[n << 1]);
}
/**
* Calculate earth circle on the sphere.
* <p>
* Returns n float lat,lon pairs at arc distance c from point at
* phi1,lambda0.
* <p>
*
* @param phi1
* latitude in radians of center point
* @param lambda0
* longitude in radians of center point
* @param c
* arc radius in radians (0 &lt; c &lt; PI)
* @param n
* number of points along circle edge to calculate
* @return float[n] radian lat,lon pairs along earth circle
*
*/
public static final float[] earthCircle(float phi1, float lambda0, float c, int n) {
return earthCircle(phi1, lambda0, c, 0.0f, MoreMath.TWO_PI, n, new float[n << 1]);
}
/**
* Calculate earth circle in the sphere.
* <p>
* Returns n float lat,lon pairs at arc distance c from point at
* phi1,lambda0.
* <p>
*
* @param phi1
* latitude in radians of center point
* @param lambda0
* longitude in radians of center point
* @param c
* arc radius in radians (0 &lt; c &lt; PI)
* @param n
* number of points along circle edge to calculate
* @param ret_val
* float[] ret_val array of n*2 number of points along circle
* edge to calculate
* @return float[n] radian lat,lon pairs along earth circle
*
*/
public static final float[] earthCircle(float phi1, float lambda0, float c, int n, float[] ret_val) {
return earthCircle(phi1, lambda0, c, 0.0f, MoreMath.TWO_PI, n, ret_val);
}
/**
* Calculate earth circle in the sphere.
* <p>
* Returns n float lat,lon pairs at arc distance c from point at
* phi1,lambda0.
* <p>
*
* @param phi1
* latitude in radians of center point.
* @param lambda0
* longitude in radians of center point.
* @param c
* arc radius in radians (0 &lt; c &lt; PI).
* @param s
* starting angle in radians. North up is zero.
* @param e
* angular extent in radians, clockwise right from starting
* angle.
* @param n
* number of points along circle edge to calculate.
* @param ret_val
* float[] ret_val array of n*2 number of points along circle
* edge to calculate.
* @return float[n] radian lat,lon pairs along earth circle.
*
*/
public static final float[] earthCircle(float phi1, float lambda0, float c, float s, float e, int n, float[] ret_val) {
double Az, cosAz, sinAz;
double cosphi1 = Math.cos(phi1);
double sinphi1 = Math.sin(phi1);
double sinc = Math.sin(c);
double cosc = Math.cos(c);
if (n < 2) n = 2; // Safety to avoid / by zero later.
int end = n << 1;// *2
// Only want to create a new return float array if there was a
// null one passed in, or if the number of desired coordinates
// is bigger than what ret_val is currently allocated for.
if (ret_val == null || end > ret_val.length) {
ret_val = new float[end];
}
double inc = e / (n - 1);
Az = s;
// generate the points in clockwise order (conforming to
// internal standard!)
for (int i = 0; i < end; i += 2, Az += inc) {
cosAz = Math.cos(Az);
sinAz = Math.sin(Az);
ret_val[i] = (float) Math.asin(sinphi1 * cosc + cosphi1 * sinc * cosAz);
ret_val[i + 1] = (float) Math.atan2(sinc * sinAz, cosphi1 * cosc - sinphi1 * sinc * cosAz) + lambda0;
}
return ret_val;
}
/**
* Calculate partial earth circle on the sphere with double precision.
* <p>
* Returns n double lat,lon pairs at arc distance c from point at
* phi1,lambda0.
* <p>
*
* @param phi1
* latitude in radians of center point
* @param lambda0
* longitude in radians of center point
* @param c
* arc radius in radians (0 &lt; c &lt; PI)
* @param s
* starting angle in radians. North up is zero.
* @param e
* angular extent in radians, clockwise right from starting
* angle.
* @param n
* number of points along circle edge to calculate
* @return double[n] radian lat,lon pairs along earth circle
*
*/
public static final double[] earthCircle(double phi1, double lambda0, double c, double s, double e, int n) {
return earthCircle(phi1, lambda0, c, s, e, n, new double[n << 1]);
}
/**
* Calculate earth circle on the sphere with double precision.
* <p>
* Returns n double lat,lon pairs at arc distance c from point at
* phi1,lambda0.
* <p>
*
* @param phi1
* latitude in radians of center point
* @param lambda0
* longitude in radians of center point
* @param c
* arc radius in radians (0 &lt; c &lt; PI)
* @param n
* number of points along circle edge to calculate
* @return double[n] radian lat,lon pairs along earth circle
*
*/
public static final double[] earthCircle(double phi1, double lambda0, double c, int n) {
return earthCircle(phi1, lambda0, c, 0.0f, MoreMath.TWO_PI_D, n, new double[n << 1]);
}
/**
* Calculate earth circle in the sphere with double precision.
* <p>
* Returns n float lat,lon pairs at arc distance c from point at
* phi1,lambda0.
* <p>
*
* @param phi1
* latitude in radians of center point
* @param lambda0
* longitude in radians of center point
* @param c
* arc radius in radians (0 &lt; c &lt; PI)
* @param n
* number of points along circle edge to calculate
* @param ret_val
* double[] ret_val array of n*2 number of points along circle
* edge to calculate
* @return double[n] radian lat,lon pairs along earth circle
*
*/
public static final double[] earthCircle(double phi1, double lambda0, double c, int n, double[] ret_val) {
return earthCircle(phi1, lambda0, c, 0.0f, MoreMath.TWO_PI_D, n, ret_val);
}
/**
* Calculate earth circle in the sphere in double precision.
* <p>
* Returns n double lat,lon pairs at arc distance c from point at
* phi1,lambda0.
* <p>
*
* @param phi1
* latitude in radians of center point.
* @param lambda0
* longitude in radians of center point.
* @param c
* arc radius in radians (0 &lt; c &lt; PI).
* @param s
* starting angle in radians. North up is zero.
* @param e
* angular extent in radians, clockwise right from starting
* angle.
* @param n
* number of points along circle edge to calculate.
* @param ret_val
* double[] ret_val array of n*2 number of points along circle
* edge to calculate.
* @return double[n] radian lat,lon pairs along earth circle.
*
*/
public static final double[] earthCircle(double phi1, double lambda0, double c, double s, double e, int n, double[] ret_val) {
double Az, cosAz, sinAz;
double cosphi1 = Math.cos(phi1);
double sinphi1 = Math.sin(phi1);
double sinc = Math.sin(c);
double cosc = Math.cos(c);
if (n < 2) n = 2; // Safety to avoid / by zero later.
int end = n << 1;// *2
// Only want to create a new return float array if there was a
// null one passed in, or if the number of desired coordinates
// is bigger than what ret_val is currently allocated for.
if (ret_val == null || end > ret_val.length) {
ret_val = new double[end];
}
double inc = e / (n - 1);
Az = s;
// generate the points in clockwise order (conforming to
// internal standard!)
for (int i = 0; i < end; i += 2, Az += inc) {
cosAz = Math.cos(Az);
sinAz = Math.sin(Az);
ret_val[i] = Math.asin(sinphi1 * cosc + cosphi1 * sinc * cosAz);
ret_val[i + 1] = Math.atan2(sinc * sinAz, cosphi1 * cosc - sinphi1 * sinc * cosAz) + lambda0;
}
return ret_val;
}
}

View file

@ -0,0 +1,339 @@
// **********************************************************************
//
// <copyright>
//
// BBN Technologies
// 10 Moulton Street
// Cambridge, MA 02138
// (617) 873-8000
//
// Copyright (C) BBNT Solutions LLC. All rights reserved.
//
// </copyright>
// **********************************************************************
package com.jwetherell.openmap.common;
public class LatLonPoint {
public final static double NORTH_POLE = 90.0;
public final static double SOUTH_POLE = -NORTH_POLE;
public final static double DATELINE = 180.0;
public final static double LON_RANGE = 360.0;
public final static double EQUIVALENT_TOLERANCE = 0.00001;
protected double lat;
protected double lon;
protected transient double radLat;
protected transient double radLon;
/**
* Default constructor, values set to 0, 0.
*/
public LatLonPoint() {
}
/**
* Set the latitude, longitude for this point in decimal degrees.
*
* @param lat
* latitude
* @param lon
* longitude.
*/
public LatLonPoint(double lat, double lon) {
setLatLon(lat, lon, false);
}
/**
* Set the latitude, longitude for this point, with the option of noting
* whether the values are in degrees or radians.
*
* @param lat
* latitude
* @param lon
* longitude.
* @param isRadians
* true of values are radians.
*/
public LatLonPoint(double lat, double lon, boolean isRadian) {
setLatLon(lat, lon, isRadian);
}
/**
* Create Double version from another LatLonPoint.
*
* @param llp
*/
public LatLonPoint(LatLonPoint llp) {
setLatLon(llp.getY(), llp.getX(), false);
}
/**
* Point2D method, inheriting signature!!
*
* @param x
* longitude value in decimal degrees.
* @param y
* latitude value in decimal degrees.
*/
public void setLocation(double x, double y) {
setLatLon(y, x, false);
}
/**
* Set latitude and longitude.
*
* @param lat
* latitude in decimal degrees.
* @param lon
* longitude in decimal degrees.
*/
public void setLatLon(double lat, double lon) {
setLatLon(lat, lon, false);
}
/**
* Set latitude and longitude.
*
* @param lat
* latitude.
* @param lon
* longitude.
* @param isRadians
* true if lat/lon values are radians.
*/
public void setLatLon(double lat, double lon, boolean isRadians) {
if (isRadians) {
radLat = lat;
radLon = lon;
this.lat = ProjMath.radToDeg(lat);
this.lon = ProjMath.radToDeg(lon);
} else {
this.lat = normalizeLatitude(lat);
this.lon = wrapLongitude(lon);
radLat = ProjMath.degToRad(lat);
radLon = ProjMath.degToRad(lon);
}
}
/**
* @return longitude in decimal degrees.
*/
public double getX() {
return lon;
}
/**
* @return latitude in decimal degrees.
*/
public double getY() {
return lat;
}
/**
* @return float latitude in decimal degrees.
*/
public float getLatitude() {
return (float) lat;
}
/**
* @return float longitude in decimal degrees.
*/
public float getLongitude() {
return (float) lon;
}
/**
* @return radian longitude.
*/
public double getRadLon() {
return radLon;
}
/**
* @return radian latitude.
*/
public double getRadLat() {
return radLat;
}
/**
* Set latitude.
*
* @param lat
* latitude in decimal degrees
*/
public void setLatitude(double lat) {
this.lat = normalizeLatitude(lat);
radLat = ProjMath.degToRad(lat);
}
/**
* Set longitude.
*
* @param lon
* longitude in decimal degrees
*/
public void setLongitude(double lon) {
this.lon = wrapLongitude(lon);
radLon = ProjMath.degToRad(lon);
}
/**
* Set location values from another lat/lon point.
*
* @param llp
*/
public void setLatLon(LatLonPoint llp) {
setLatLon(llp.getY(), llp.getX(), false);
}
/**
* Ensure latitude is between the poles.
*
* @param lat
* @return
*/
public final static float normalizeLatitude(float lat) {
return (float) normalizeLatitude((double) lat);
}
/**
* Sets latitude to something sane.
*
* @param lat
* latitude in decimal degrees
* @return float normalized latitude in decimal degrees (&minus;90&deg; &le;
* &phi; &le; 90&deg;)
*/
public final static double normalizeLatitude(double lat) {
if (lat > NORTH_POLE) {
lat = NORTH_POLE;
}
if (lat < SOUTH_POLE) {
lat = SOUTH_POLE;
}
return lat;
}
/**
* Ensure the longitude is between the date line.
*
* @param lon
* @return
*/
public final static float wrapLongitude(float lon) {
return (float) wrapLongitude((double) lon);
}
/**
* Sets longitude to something sane.
*
* @param lon
* longitude in decimal degrees
* @return float wrapped longitude in decimal degrees (&minus;180&deg; &le;
* &lambda; &le; 180&deg;)
*/
public final static double wrapLongitude(double lon) {
if ((lon < -DATELINE) || (lon > DATELINE)) {
lon += DATELINE;
lon = lon % LON_RANGE;
lon = (lon < 0) ? DATELINE + lon : -DATELINE + lon;
}
return lon;
}
/**
* Check if latitude is bogus. Latitude is invalid if lat &gt; 90&deg; or if
* lat &lt; &minus;90&deg;.
*
* @param lat
* latitude in decimal degrees
* @return boolean true if latitude is invalid
*/
public static boolean isInvalidLatitude(double lat) {
return ((lat > NORTH_POLE) || (lat < SOUTH_POLE));
}
/**
* Check if longitude is bogus. Longitude is invalid if lon &gt; 180&deg; or
* if lon &lt; &minus;180&deg;.
*
* @param lon
* longitude in decimal degrees
* @return boolean true if longitude is invalid
*/
public static boolean isInvalidLongitude(double lon) {
return ((lon < -DATELINE) || (lon > DATELINE));
}
/**
* Determines whether two LatLonPoints are equal.
*
* @param obj
* Object
* @return Whether the two points are equal up to a tolerance of 10 <sup>-5
* </sup> degrees in latitude and longitude.
*/
public boolean equals(Object obj) {
if (obj == null) {
return false;
}
if (getClass() != obj.getClass()) {
return false;
}
final LatLonPoint pt = (LatLonPoint) obj;
return (MoreMath.approximately_equal(getY(), pt.getY(), EQUIVALENT_TOLERANCE) && MoreMath.approximately_equal(getX(), pt.getX(), EQUIVALENT_TOLERANCE));
}
/**
* Find the distance to another LatLonPoint, based on a earth spherical
* model.
*
* @param toPoint
* LatLonPoint
* @return distance, in radians. You can use an com.bbn.openmap.proj.Length
* to convert the radians to other units.
*/
public double distance(LatLonPoint toPoint) {
return GreatCircle.sphericalDistance(getRadLat(), getRadLon(), toPoint.getRadLat(), toPoint.getRadLon());
}
/**
* Find the azimuth to another point, based on the spherical earth model.
*
* @param toPoint
* LatLonPoint
* @return the azimuth `Az' east of north from this point bearing toward the
* one provided as an argument.(-PI &lt;= Az &lt;= PI).
*
*/
public double azimuth(LatLonPoint toPoint) {
return GreatCircle.sphericalAzimuth(getRadLat(), getRadLon(), toPoint.getRadLat(), toPoint.getRadLon());
}
/**
* Get a new LatLonPoint a distance and azimuth from another point, based on
* the spherical earth model.
*
* @param distance
* radians
* @param azimuth
* radians
* @return
*/
public LatLonPoint getPoint(double distance, double azimuth) {
return GreatCircle.sphericalBetween(getRadLat(), getRadLon(), distance, azimuth);
}
/**
* {@inheritDoc}
*/
@Override
public String toString() {
return "Lat=" + lat + ", Lon=" + lon;
}
}

View file

@ -0,0 +1,643 @@
// **********************************************************************
//
// <copyright>
//
// BBN Technologies
// 10 Moulton Street
// Cambridge, MA 02138
// (617) 873-8000
//
// Copyright (C) BBNT Solutions LLC. All rights reserved.
//
// </copyright>
// **********************************************************************
package com.jwetherell.openmap.common;
public class MGRSPoint extends ZonedUTMPoint {
/**
* UTM zones are grouped, and assigned to one of a group of 6 sets.
*/
private final static int NUM_100K_SETS = 6;
/**
* The column letters (for easting) of the lower left value, per set.
*/
private final static int[] SET_ORIGIN_COLUMN_LETTERS = { 'A', 'J', 'S', 'A', 'J', 'S' };
/**
* The row letters (for northing) of the lower left value, per set.
*/
private final static int[] SET_ORIGIN_ROW_LETTERS = { 'A', 'F', 'A', 'F', 'A', 'F' };
public final static int ACCURACY_1_METER = 5;
public final static int ACCURACY_10_METER = 4;
public final static int ACCURACY_100_METER = 3;
public final static int ACCURACY_1000_METER = 2;
public final static int ACCURACY_10000_METER = 1;
/** The set origin column letters to use. */
private int[] originColumnLetters = SET_ORIGIN_COLUMN_LETTERS;
/** The set origin row letters to use. */
private int[] originRowLetters = SET_ORIGIN_ROW_LETTERS;
private final static int A = 'A';
private final static int I = 'I';
private final static int O = 'O';
private final static int V = 'V';
private final static int Z = 'Z';
/** The String holding the MGRS coordinate value. */
private String mgrs = null;
/**
* Controls the number of digits that the MGRS coordinate will have, which
* directly affects the accuracy of the coordinate. Default is
* ACCURACY_1_METER, which indicates that MGRS coordinates will have 10
* digits (5 easting, 5 northing) after the 100k two letter code, indicating
* 1 meter resolution.
*/
private int accuracy = ACCURACY_1_METER;
public MGRSPoint() {
}
public MGRSPoint(LatLonPoint llpoint) {
this(llpoint, Ellipsoid.WGS_84);
}
public MGRSPoint(LatLonPoint llpoint, Ellipsoid ellip) {
this();
LLtoMGRS(llpoint, ellip, this);
}
public MGRSPoint(String mgrsString) throws NumberFormatException {
this();
setMGRS(mgrsString);
}
public MGRSPoint(double northing, double easting, int zoneNumber, char zoneLetter) {
super(northing, easting, zoneNumber, zoneLetter);
}
/**
* Set the MGRS value for this Point. Will be decoded, and the UTM values
* figured out. You can call toLatLonPoint() to translate it to lat/lon
* decimal degrees.
*/
public void setMGRS(String mgrsString) throws NumberFormatException {
try {
mgrs = mgrsString.toUpperCase(); // Just to make sure.
decode(mgrs);
} catch (StringIndexOutOfBoundsException sioobe) {
throw new NumberFormatException("MGRSPoint has bad string: " + mgrsString);
} catch (NullPointerException npe) {
// Blow off
}
}
/**
* Set the UTM parameters from a MGRS string.
*
* @param mgrsString
* an UPPERCASE coordinate string is expected.
*/
protected void decode(String mgrsString) throws NumberFormatException {
if (mgrsString == null || mgrsString.length() == 0) {
throw new NumberFormatException("MGRSPoint coverting from nothing");
}
// Ensure an upper-case string
mgrsString = mgrsString.toUpperCase();
int length = mgrsString.length();
String hunK = null;
StringBuffer sb = new StringBuffer();
char testChar;
int i = 0;
// get Zone number
while (!Character.isLetter(testChar = mgrsString.charAt(i))) {
if (i >= 2) {
throw new NumberFormatException("MGRSPoint bad conversion from: " + mgrsString + ", first two characters need to be a number between 1-60.");
}
sb.append(testChar);
i++;
}
zone_number = Integer.parseInt(sb.toString());
if (zone_number < 1 || zone_number > 60) {
throw new NumberFormatException("MGRSPoint bad conversion from: " + mgrsString + ", first two characters need to be a number between 1-60.");
}
if (i == 0 || i + 3 > length) {
// A good MGRS string has to be 4-5 digits long,
// ##AAA/#AAA at least.
throw new NumberFormatException("MGRSPoint bad conversion from: " + mgrsString + ", MGRS string must be at least 4-5 digits long");
}
zone_letter = mgrsString.charAt(i++);
// Should we check the zone letter here? Why not.
if (zone_letter <= 'A' || zone_letter == 'B' || zone_letter == 'Y' || zone_letter >= 'Z' || zone_letter == 'I' || zone_letter == 'O') {
throw new NumberFormatException("MGRSPoint zone letter " + zone_letter + " not handled: " + mgrsString);
}
hunK = mgrsString.substring(i, i += 2);
// Validate, check the zone, make sure each letter is between A-Z, not I
// or O
char char1 = hunK.charAt(0);
char char2 = hunK.charAt(1);
if (char1 < 'A' || char2 < 'A' || char1 > 'Z' || char2 > 'Z' || char1 == 'I' || char2 == 'I' || char1 == 'O' || char2 == 'O') {
throw new NumberFormatException("MGRSPoint bad conversion from " + mgrsString + ", invalid 100k designator");
}
int set = get100kSetForZone(zone_number);
float east100k = getEastingFromChar(char1, set);
float north100k = getNorthingFromChar(char2, set);
// We have a bug where the northing may be 2000000 too low.
// How do we know when to roll over?
while (north100k < getMinNorthing(zone_letter)) {
north100k += 2000000;
}
// calculate the char index for easting/northing separator
int remainder = length - i;
if (remainder % 2 != 0) {
throw new NumberFormatException(
"MGRSPoint has to have an even number \nof digits after the zone letter and two 100km letters - front \nhalf for easting meters, second half for \nnorthing meters"
+ mgrsString);
}
int sep = remainder / 2;
float sepEasting = 0f;
float sepNorthing = 0f;
if (sep > 0) {
float accuracyBonus = 100000f / (float) Math.pow(10, sep);
String sepEastingString = mgrsString.substring(i, i + sep);
sepEasting = Float.parseFloat(sepEastingString) * accuracyBonus;
String sepNorthingString = mgrsString.substring(i + sep);
sepNorthing = Float.parseFloat(sepNorthingString) * accuracyBonus;
}
easting = sepEasting + east100k;
northing = sepNorthing + north100k;
}
/**
* Given the first letter from a two-letter MGRS 100k zone, and given the
* MGRS table set for the zone number, figure out the easting value that
* should be added to the other, secondary easting value.
*/
protected float getEastingFromChar(char e, int set) {
int baseCol[] = getOriginColumnLetters();
// colOrigin is the letter at the origin of the set for the
// column
int curCol = baseCol[set - 1];
float eastingValue = 100000f;
boolean rewindMarker = false;
while (curCol != e) {
curCol++;
if (curCol == I) curCol++;
if (curCol == O) curCol++;
if (curCol > Z) {
if (rewindMarker) {
throw new NumberFormatException("Bad character: " + e);
}
curCol = A;
rewindMarker = true;
}
eastingValue += 100000f;
}
return eastingValue;
}
/**
* Given the second letter from a two-letter MGRS 100k zone, and given the
* MGRS table set for the zone number, figure out the northing value that
* should be added to the other, secondary northing value. You have to
* remember that Northings are determined from the equator, and the vertical
* cycle of letters mean a 2000000 additional northing meters. This happens
* approx. every 18 degrees of latitude. This method does *NOT* count any
* additional northings. You have to figure out how many 2000000 meters need
* to be added for the zone letter of the MGRS coordinate.
*
* @param n
* second letter of the MGRS 100k zone
* @param set
* the MGRS table set number, which is dependent on the UTM zone
* number.
*/
protected float getNorthingFromChar(char n, int set) {
if (n > 'V') {
throw new NumberFormatException("MGRSPoint given invalid Northing " + n);
}
int baseRow[] = getOriginRowLetters();
// rowOrigin is the letter at the origin of the set for the
// column
int curRow = baseRow[set - 1];
float northingValue = 0f;
boolean rewindMarker = false;
while (curRow != n) {
curRow++;
if (curRow == I) curRow++;
if (curRow == O) curRow++;
// fixing a bug making whole application hang in this loop
// when 'n' is a wrong character
if (curRow > V) {
if (rewindMarker) { // making sure that this loop ends
throw new NumberFormatException("Bad character: " + n);
}
curRow = A;
rewindMarker = true;
}
northingValue += 100000f;
}
return northingValue;
}
/**
* The function getMinNorthing returns the minimum northing value of a MGRS
* zone.
*
* portted from Geotrans' c Latitude_Band_Value structure table. zoneLetter
* : MGRS zone (input)
*/
protected float getMinNorthing(char zoneLetter) throws NumberFormatException {
float northing;
switch (zoneLetter) {
case 'C':
northing = 1100000.0f;
break;
case 'D':
northing = 2000000.0f;
break;
case 'E':
northing = 2800000.0f;
break;
case 'F':
northing = 3700000.0f;
break;
case 'G':
northing = 4600000.0f;
break;
case 'H':
northing = 5500000.0f;
break;
case 'J':
northing = 6400000.0f;
break;
case 'K':
northing = 7300000.0f;
break;
case 'L':
northing = 8200000.0f;
break;
case 'M':
northing = 9100000.0f;
break;
case 'N':
northing = 0.0f;
break;
case 'P':
northing = 800000.0f;
break;
case 'Q':
northing = 1700000.0f;
break;
case 'R':
northing = 2600000.0f;
break;
case 'S':
northing = 3500000.0f;
break;
case 'T':
northing = 4400000.0f;
break;
case 'U':
northing = 5300000.0f;
break;
case 'V':
northing = 6200000.0f;
break;
case 'W':
northing = 7000000.0f;
break;
case 'X':
northing = 7900000.0f;
break;
default:
northing = -1.0f;
}
if (northing >= 0.0)
return northing;
// else
throw new NumberFormatException("Invalid zone letter: " + zone_letter);
}
/**
* Convert this MGRSPoint to a LatLonPoint, and assume a WGS_84 ellipsoid.
*/
public LatLonPoint toLatLonPoint() {
return toLatLonPoint(Ellipsoid.WGS_84, new LatLonPoint());
}
/**
* Convert this MGRSPoint to a LatLonPoint, and use the given ellipsoid.
*/
public LatLonPoint toLatLonPoint(Ellipsoid ellip) {
return toLatLonPoint(ellip, new LatLonPoint());
}
/**
* Fill in the given LatLonPoint with the converted values of this
* MGRSPoint, and use the given ellipsoid.
*/
public LatLonPoint toLatLonPoint(Ellipsoid ellip, LatLonPoint llpoint) {
return MGRStoLL(this, ellip, llpoint);
}
/**
* Create a LatLonPoint from a MGRSPoint.
*
* @param mgrsp
* to convert.
* @param ellip
* Ellipsoid for earth model.
* @param llp
* a LatLonPoint to fill in values for. If null, a new
* LatLonPoint will be returned. If not null, the new values will
* be set in this object, and it will be returned.
* @return LatLonPoint with values converted from MGRS coordinate.
*/
public static LatLonPoint MGRStoLL(MGRSPoint mgrsp, Ellipsoid ellip, LatLonPoint llp) {
return UTMtoLL(ellip, mgrsp.northing, mgrsp.easting, mgrsp.zone_number, MGRSPoint.MGRSZoneToUTMZone(mgrsp.zone_letter), llp);
}
/**
* Converts a LatLonPoint to a MGRS Point, assuming the WGS_84 ellipsoid.
*
* @return MGRSPoint, or null if something bad happened.
*/
public static MGRSPoint LLtoMGRS(LatLonPoint llpoint) {
return LLtoMGRS(llpoint, Ellipsoid.WGS_84, new MGRSPoint());
}
/**
* Create a MGRSPoint from a LatLonPoint.
*
* @param llp
* LatLonPoint to convert.
* @param ellip
* Ellipsoid for earth model.
* @param mgrsp
* a MGRSPoint to fill in values for. If null, a new MGRSPoint
* will be returned. If not null, the new values will be set in
* this object, and it will be returned.
* @return MGRSPoint with values converted from lat/lon.
*/
public static MGRSPoint LLtoMGRS(LatLonPoint llp, Ellipsoid ellip, MGRSPoint mgrsp) {
if (mgrsp == null) {
mgrsp = new MGRSPoint();
}
// Calling LLtoUTM here results in N/S zone letters! wrong!
mgrsp = (MGRSPoint) LLtoUTM(llp, ellip, mgrsp);
// Need to add this to set the right letter for the latitude.
mgrsp.zone_letter = mgrsp.getLetterDesignator(llp.getLatitude());
mgrsp.resolve();
return mgrsp;
}
/**
* Convert MGRS zone letter to UTM zone letter, N or S.
*
* @param mgrsZone
* @return N of given zone is equal or larger than N, S otherwise.
*/
public static char MGRSZoneToUTMZone(char mgrsZone) {
if (Character.toUpperCase(mgrsZone) >= 'N')
return 'N';
// else
return 'S';
}
/**
* Method that provides a check for MGRS zone letters. Returns an uppercase
* version of any valid letter passed in.
*/
public char checkZone(char zone) {
zone = Character.toUpperCase(zone);
if (zone <= 'A' || zone == 'B' || zone == 'Y' || zone >= 'Z' || zone == 'I' || zone == 'O') {
throw new NumberFormatException("Invalid MGRSPoint zone letter: " + zone);
}
return zone;
}
/**
* Set the number of digits to use for easting and northing numbers in the
* mgrs string, which reflects the accuracy of the coordinate. From 5 (1
* meter) to 1 (10,000 meter).
*/
public void setAccuracy(int value) {
accuracy = value;
mgrs = null;
}
public int getAccuracy() {
return accuracy;
}
/**
* Create the mgrs string based on the internal UTM settings, should be
* called if the accuracy changes.
*
* @param digitAccuracy
* The number of digits to use for the northing and easting
* numbers. 5 digits reflect a 1 meter accuracy, 4 - 10 meter, 3
* - 100 meter, 2 - 1000 meter, 1 - 10,000 meter.
*/
public void resolve(int digitAccuracy) {
setAccuracy(digitAccuracy);
resolve();
}
/**
* Create the mgrs string based on the internal UTM settings, using the
* accuracy set in the MGRSPoint.
*/
public void resolve() {
if (zone_letter == 'Z') {
mgrs = "Latitude limit exceeded";
} else {
StringBuffer sb = new StringBuffer(Integer.toString(zone_number)).append(zone_letter).append(get100kID(easting, northing, zone_number));
StringBuffer seasting = new StringBuffer(Integer.toString((int) easting));
StringBuffer snorthing = new StringBuffer(Integer.toString((int) northing));
while (accuracy + 1 > seasting.length()) {
seasting.insert(0, '0');
}
// We have to be careful here, the 100k values shouldn't
// be
// used for calculating stuff here.
while (accuracy + 1 > snorthing.length()) {
snorthing.insert(0, '0');
}
while (snorthing.length() > 6) {
snorthing.deleteCharAt(0);
}
try {
sb.append(seasting.substring(1, accuracy + 1)).append(snorthing.substring(1, accuracy + 1));
mgrs = sb.toString();
} catch (IndexOutOfBoundsException ioobe) {
mgrs = null;
}
}
}
/**
* Given a UTM zone number, figure out the MGRS 100K set it is in.
*/
private int get100kSetForZone(int i) {
int set = i % NUM_100K_SETS;
if (set == 0) set = NUM_100K_SETS;
return set;
}
/**
* Provided so that extensions to this class can provide different origin
* letters, in case of different ellipsoids. The int[] represents all of the
* first letters in the bottom left corner of each set box, as shown in an
* MGRS 100K box layout.
*/
private int[] getOriginColumnLetters() {
return originColumnLetters;
}
/**
* Provided so that extensions to this class can provide different origin
* letters, in case of different ellipsoids. The int[] represents all of the
* second letters in the bottom left corner of each set box, as shown in an
* MGRS 100K box layout.
*/
private int[] getOriginRowLetters() {
return originRowLetters;
}
/**
* Get the two letter 100k designator for a given UTM easting, northing and
* zone number value.
*/
private String get100kID(double easting, double northing, int zone_number) {
int set = get100kSetForZone(zone_number);
int setColumn = ((int) easting / 100000);
int setRow = ((int) northing / 100000) % 20;
return get100kID(setColumn, setRow, set);
}
/**
* Get the two-letter MGRS 100k designator given information translated from
* the UTM northing, easting and zone number.
*
* @param setColumn
* the column index as it relates to the MGRS 100k set
* spreadsheet, created from the UTM easting. Values are 1-8.
* @param setRow
* the row index as it relates to the MGRS 100k set spreadsheet,
* created from the UTM northing value. Values are from 0-19.
* @param set
* the set block, as it relates to the MGRS 100k set spreadsheet,
* created from the UTM zone. Values are from 1-60.
* @return two letter MGRS 100k code.
*/
private String get100kID(int setColumn, int setRow, int set) {
int baseCol[] = getOriginColumnLetters();
int baseRow[] = getOriginRowLetters();
// colOrigin and rowOrigin are the letters at the origin of
// the set
int colOrigin = baseCol[set - 1];
int rowOrigin = baseRow[set - 1];
// colInt and rowInt are the letters to build to return
int colInt = colOrigin + setColumn - 1;
int rowInt = rowOrigin + setRow;
boolean rollover = false;
if (colInt > Z) {
colInt = colInt - Z + A - 1;
rollover = true;
}
if (colInt == I || (colOrigin < I && colInt > I) || ((colInt > I || colOrigin < I) && rollover)) {
colInt++;
}
if (colInt == O || (colOrigin < O && colInt > O) || ((colInt > O || colOrigin < O) && rollover)) {
colInt++;
if (colInt == I) {
colInt++;
}
}
if (colInt > Z) {
colInt = colInt - Z + A - 1;
}
if (rowInt > V) {
rowInt = rowInt - V + A - 1;
rollover = true;
} else {
rollover = false;
}
if (rowInt == I || (rowOrigin < I && rowInt > I) || ((rowInt > I || rowOrigin < I) && rollover)) {
rowInt++;
}
if (rowInt == O || (rowOrigin < O && rowInt > O) || ((rowInt > O || rowOrigin < O) && rollover)) {
rowInt++;
if (rowInt == I) {
rowInt++;
}
}
if (rowInt > V) {
rowInt = rowInt - V + A - 1;
}
String twoLetter = (char) colInt + "" + (char) rowInt;
return twoLetter;
}
/**
* {@inheritDoc}
*/
@Override
public String toString() {
return mgrs;
}
}

View file

@ -0,0 +1,546 @@
// **********************************************************************
//
// <copyright>
//
// BBN Technologies
// 10 Moulton Street
// Cambridge, MA 02138
// (617) 873-8000
//
// Copyright (C) BBNT Solutions LLC. All rights reserved.
//
// </copyright>
// **********************************************************************
package com.jwetherell.openmap.common;
public abstract class MoreMath {
/**
* 2*Math.PI
*/
public static final transient float TWO_PI = (float) Math.PI * 2.0f;
/**
* 2*Math.PI
*/
public static final transient double TWO_PI_D = Math.PI * 2.0d;
/**
* Math.PI/2
*/
public static final transient float HALF_PI = (float) Math.PI / 2.0f;
/**
* Math.PI/2
*/
public static final transient double HALF_PI_D = Math.PI / 2.0d;
/**
* Checks if a ~= b. Use this to test equality of floating point numbers.
* <p>
*
* @param a
* double
* @param b
* double
* @param epsilon
* the allowable error
* @return boolean
*/
public static final boolean approximately_equal(double a, double b, double epsilon) {
return (Math.abs(a - b) <= epsilon);
}
/**
* Checks if a ~= b. Use this to test equality of floating point numbers.
* <p>
*
* @param a
* float
* @param b
* float
* @param epsilon
* the allowable error
* @return boolean
*/
public static final boolean approximately_equal(float a, float b, float epsilon) {
return (Math.abs(a - b) <= epsilon);
}
/**
* Hyperbolic arcsin.
* <p>
* Hyperbolic arc sine: log (x+sqrt(1+x^2))
*
* @param x
* float
* @return float asinh(x)
*/
public static final float asinh(float x) {
return (float) Math.log(x + Math.sqrt(x * x + 1));
}
/**
* Hyperbolic arcsin.
* <p>
* Hyperbolic arc sine: log (x+sqrt(1+x^2))
*
* @param x
* double
* @return double asinh(x)
*/
public static final double asinh(double x) {
return Math.log(x + Math.sqrt(x * x + 1));
}
/**
* Hyperbolic sin.
* <p>
* Hyperbolic sine: (e^x-e^-x)/2
*
* @param x
* float
* @return float sinh(x)
*/
public static final float sinh(float x) {
return (float) (Math.pow(Math.E, x) - Math.pow(Math.E, -x)) / 2.0f;
}
/**
* Hyperbolic sin.
* <p>
* Hyperbolic sine: (e^x-e^-x)/2
*
* @param x
* double
* @return double sinh(x)
*/
public static final double sinh(double x) {
return (Math.pow(Math.E, x) - Math.pow(Math.E, -x)) / 2.0d;
}
// HACK - are there functions that already exist?
/**
* Return sign of number.
*
* @param x
* short
* @return int sign -1, 1
*/
public static final int sign(short x) {
return (x < 0) ? -1 : 1;
}
/**
* Return sign of number.
*
* @param x
* int
* @return int sign -1, 1
*/
public static final int sign(int x) {
return (x < 0) ? -1 : 1;
}
/**
* Return sign of number.
*
* @param x
* long
* @return int sign -1, 1
*/
public static final int sign(long x) {
return (x < 0) ? -1 : 1;
}
/**
* Return sign of number.
*
* @param x
* float
* @return int sign -1, 1
*/
public static final int sign(float x) {
return (x < 0f) ? -1 : 1;
}
/**
* Return sign of number.
*
* @param x
* double
* @return int sign -1, 1
*/
public static final int sign(double x) {
return (x < 0d) ? -1 : 1;
}
/**
* Check if number is odd.
*
* @param x
* short
* @return boolean
*/
public static final boolean odd(short x) {
return !even(x);
}
/**
* Check if number is odd.
*
* @param x
* int
* @return boolean
*/
public static final boolean odd(int x) {
return !even(x);
}
/**
* Check if number is odd.
*
* @param x
* long
* @return boolean
*/
public static final boolean odd(long x) {
return !even(x);
}
/**
* Check if number is even.
*
* @param x
* short
* @return boolean
*/
public static final boolean even(short x) {
return ((x & 0x1) == 0);
}
/**
* Check if number is even.
*
* @param x
* int
* @return boolean
*/
public static final boolean even(int x) {
return ((x & 0x1) == 0);
}
/**
* Check if number is even.
*
* @param x
* long
* @return boolean
*/
public static final boolean even(long x) {
return ((x & 0x1) == 0);
}
/**
* Converts a byte in the range of -128 to 127 to an int in the range 0 -
* 255.
*
* @param b
* (-128 &lt;= b &lt;= 127)
* @return int (0 &lt;= b &lt;= 255)
*/
public static final int signedToInt(byte b) {
return (b & 0xff);
}
/**
* Converts a short in the range of -32768 to 32767 to an int in the range 0
* - 65535.
*
* @param w
* (-32768 &lt;= b &lt;= 32767)
* @return int (0 &lt;= b &lt;= 65535)
*/
public static final int signedToInt(short w) {
return (w & 0xffff);
}
/**
* Convert an int in the range of -2147483648 to 2147483647 to a long in the
* range 0 to 4294967295.
*
* @param x
* (-2147483648 &lt;= x &lt;= 2147483647)
* @return long (0 &lt;= x &lt;= 4294967295)
*/
public static final long signedToLong(int x) {
return (x & 0xFFFFFFFFL);
}
/**
* Converts an int in the range of 0 - 65535 to an int in the range of 0 -
* 255.
*
* @param w
* int (0 &lt;= w &lt;= 65535)
* @return int (0 &lt;= w &lt;= 255)
*/
public static final int wordToByte(int w) {
return w >> 8;
}
/**
* Build short out of bytes (in big endian order).
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @return short
*/
public static final short BuildShortBE(byte bytevec[], int offset) {
return (short) (((bytevec[0 + offset]) << 8) | (signedToInt(bytevec[1 + offset])));
}
/**
* Build short out of bytes (in little endian order).
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @return short
*/
public static final short BuildShortLE(byte bytevec[], int offset) {
return (short) (((bytevec[1 + offset]) << 8) | (signedToInt(bytevec[0 + offset])));
}
/**
* Build short out of bytes.
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @param MSBFirst
* BE or LE?
* @return short
*/
public static final short BuildShort(byte bytevec[], int offset, boolean MSBFirst) {
if (MSBFirst)
return (BuildShortBE(bytevec, offset));
// else
return (BuildShortLE(bytevec, offset));
}
/**
* Build short out of bytes (in big endian order).
*
* @param bytevec
* bytes
* @param MSBFirst
* BE or LE?
* @return short
*/
public static final short BuildShortBE(byte bytevec[], boolean MSBFirst) {
return BuildShortBE(bytevec, 0);
}
/**
* Build short out of bytes (in little endian order).
*
* @param bytevec
* bytes
* @param MSBFirst
* BE or LE?
* @return short
*/
public static final short BuildShortLE(byte bytevec[], boolean MSBFirst) {
return BuildShortLE(bytevec, 0);
}
/**
* Build short out of bytes.
*
* @param bytevec
* bytes
* @param MSBFirst
* BE or LE?
* @return short
*/
public static final short BuildShort(byte bytevec[], boolean MSBFirst) {
return BuildShort(bytevec, 0, MSBFirst);
}
/**
* Build int out of bytes (in big endian order).
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @return int
*/
public static final int BuildIntegerBE(byte bytevec[], int offset) {
return (((bytevec[0 + offset]) << 24) | (signedToInt(bytevec[1 + offset]) << 16) | (signedToInt(bytevec[2 + offset]) << 8) | (signedToInt(bytevec[3 + offset])));
}
/**
* Build int out of bytes (in little endian order).
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @return int
*/
public static final int BuildIntegerLE(byte bytevec[], int offset) {
return (((bytevec[3 + offset]) << 24) | (signedToInt(bytevec[2 + offset]) << 16) | (signedToInt(bytevec[1 + offset]) << 8) | (signedToInt(bytevec[0 + offset])));
}
/**
* Build int out of bytes.
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @param MSBFirst
* BE or LE?
* @return int
*/
public static final int BuildInteger(byte bytevec[], int offset, boolean MSBFirst) {
if (MSBFirst)
return BuildIntegerBE(bytevec, offset);
// else
return BuildIntegerLE(bytevec, offset);
}
/**
* Build int out of bytes (in big endian order).
*
* @param bytevec
* bytes
* @return int
*/
public static final int BuildIntegerBE(byte bytevec[]) {
return BuildIntegerBE(bytevec, 0);
}
/**
* Build int out of bytes (in little endian order).
*
* @param bytevec
* bytes
* @return int
*/
public static final int BuildIntegerLE(byte bytevec[]) {
return BuildIntegerLE(bytevec, 0);
}
/**
* Build int out of bytes.
*
* @param bytevec
* bytes
* @param MSBFirst
* BE or LE?
* @return int
*/
public static final int BuildInteger(byte bytevec[], boolean MSBFirst) {
if (MSBFirst)
return BuildIntegerBE(bytevec, 0);
//else
return BuildIntegerLE(bytevec, 0);
}
/**
* Build long out of bytes (in big endian order).
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @return long
*/
public static final long BuildLongBE(byte bytevec[], int offset) {
return (((long) signedToInt(bytevec[0 + offset]) << 56) | ((long) signedToInt(bytevec[1 + offset]) << 48)
| ((long) signedToInt(bytevec[2 + offset]) << 40) | ((long) signedToInt(bytevec[3 + offset]) << 32)
| ((long) signedToInt(bytevec[4 + offset]) << 24) | ((long) signedToInt(bytevec[5 + offset]) << 16)
| ((long) signedToInt(bytevec[6 + offset]) << 8) | (signedToInt(bytevec[7 + offset])));
}
/**
* Build long out of bytes (in little endian order).
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @return long
*/
public static final long BuildLongLE(byte bytevec[], int offset) {
return (((long) signedToInt(bytevec[7 + offset]) << 56) | ((long) signedToInt(bytevec[6 + offset]) << 48)
| ((long) signedToInt(bytevec[5 + offset]) << 40) | ((long) signedToInt(bytevec[4 + offset]) << 32)
| ((long) signedToInt(bytevec[3 + offset]) << 24) | ((long) signedToInt(bytevec[2 + offset]) << 16)
| ((long) signedToInt(bytevec[1 + offset]) << 8) | (signedToInt(bytevec[0 + offset])));
}
/**
* Build long out of bytes.
*
* @param bytevec
* bytes
* @param offset
* byte offset
* @param MSBFirst
* BE or LE?
* @return long
*/
public static final long BuildLong(byte bytevec[], int offset, boolean MSBFirst) {
if (MSBFirst)
return BuildLongBE(bytevec, offset);
// else
return BuildLongLE(bytevec, offset);
}
/**
* Build long out of bytes (in big endian order).
*
* @param bytevec
* bytes
* @return long
*/
public static final long BuildLongBE(byte bytevec[]) {
return BuildLongBE(bytevec, 0);
}
/**
* Build long out of bytes (in little endian order).
*
* @param bytevec
* bytes
* @return long
*/
public static final long BuildLongLE(byte bytevec[]) {
return BuildLongLE(bytevec, 0);
}
/**
* Build long out of bytes.
*
* @param bytevec
* bytes
* @param MSBFirst
* BE or LE?
* @return long
*/
public static final long BuildLong(byte bytevec[], boolean MSBFirst) {
if (MSBFirst)
return BuildLongBE(bytevec, 0);
// else
return BuildLongLE(bytevec, 0);
}
}

View file

@ -0,0 +1,465 @@
// **********************************************************************
//
// <copyright>
//
// BBN Technologies
// 10 Moulton Street
// Cambridge, MA 02138
// (617) 873-8000
//
// Copyright (C) BBNT Solutions LLC. All rights reserved.
//
// </copyright>
// **********************************************************************
package com.jwetherell.openmap.common;
public abstract class ProjMath {
/**
* North pole latitude in radians.
*/
public static final transient float NORTH_POLE_F = MoreMath.HALF_PI;
/**
* South pole latitude in radians.
*/
public static final transient float SOUTH_POLE_F = -NORTH_POLE_F;
/**
* North pole latitude in radians.
*/
public static final transient double NORTH_POLE_D = MoreMath.HALF_PI_D;
/**
* North pole latitude in degrees.
*/
public static final transient double NORTH_POLE_DEG_D = 90d;
/**
* South pole latitude in radians.
*/
public static final transient double SOUTH_POLE_D = -NORTH_POLE_D;
/**
* South pole latitude in degrees.
*/
public static final transient double SOUTH_POLE_DEG_D = -NORTH_POLE_DEG_D;
/**
* Dateline longitude in radians.
*/
public static final transient float DATELINE_F = (float) Math.PI;
/**
* Dateline longitude in radians.
*/
public static final transient double DATELINE_D = Math.PI;
/**
* Dateline longitude in degrees.
*/
public static final transient double DATELINE_DEG_D = 180d;
/**
* Longitude range in radians.
*/
public static final transient float LON_RANGE_F = MoreMath.TWO_PI;
/**
* Longitude range in radians.
*/
public static final transient double LON_RANGE_D = MoreMath.TWO_PI_D;
/**
* Longitude range in degrees.
*/
public static final transient double LON_RANGE_DEG_D = 360d;
public static final double DEGREES_TO_MILS = 17.77777777777777777778;
/**
* rounds the quantity away from 0.
*
* @param x
* in value
* @return double
* @see #qint(double)
*/
public static final double roundAdjust(double x) {
return qint_old(x);
}
/**
* Rounds the quantity away from 0.
*
* @param x
* value
* @return double
*/
public static final double qint(double x) {
return qint_new(x);
}
private static final double qint_old(double x) {
return (((int) x) < 0) ? (x - 0.5) : (x + 0.5);
}
private static final double qint_new(double x) {
// -1 or +1 away from zero
return (x <= 0.0) ? (x - 1.0) : (x + 1.0);
}
/**
* Calculate the shortest arc distance between two lons.
*
* @param lon1
* radians
* @param lon2
* radians
* @return float distance
*/
public static final float lonDistance(float lon1, float lon2) {
return (float) Math.min(Math.abs(lon1 - lon2), ((lon1 < 0) ? lon1 + Math.PI : Math.PI - lon1) + ((lon2 < 0) ? lon2 + Math.PI : Math.PI - lon2));
}
/**
* Convert between decimal degrees and scoords.
*
* @param deg
* degrees
* @return long scoords
*
*/
public static final long DEG_TO_SC(double deg) {
return (long) (deg * 3600000);
}
/**
* Convert between decimal degrees and scoords.
*
* @param sc
* scoords
* @return double decimal degrees
*/
public static final double SC_TO_DEG(int sc) {
return ((sc) / (60.0 * 60.0 * 1000.0));
}
/**
* Convert radians to mils.
*
* @param rad
* radians
* @return double mils
*/
public static final double radToMils(double rad) {
double degrees = Math.toDegrees(rad);
degrees = (degrees < 0) ? 360 + degrees : degrees;
return degrees * DEGREES_TO_MILS;
}
/**
* Convert radians to degrees.
*
* @param rad
* radians
* @return double decimal degrees
*/
public static final double radToDeg(double rad) {
return Math.toDegrees(rad);
}
/**
* Convert radians to degrees.
*
* @param rad
* radians
* @return float decimal degrees
*/
public static final float radToDeg(float rad) {
return (float) Math.toDegrees(rad);
}
/**
* Convert degrees to radians.
*
* @param deg
* degrees
* @return double radians
*/
public static final double degToRad(double deg) {
return Math.toRadians(deg);
}
/**
* Convert degrees to radians.
*
* @param deg
* degrees
* @return float radians
*/
public static final float degToRad(float deg) {
return (float) Math.toRadians(deg);
}
/**
* Generate a hashCode value for a lat/lon pair.
*
* @param lat
* latitude
* @param lon
* longitude
* @return int hashcode
*
*/
public static final int hashLatLon(float lat, float lon) {
if (lat == -0f) lat = 0f;// handle negative zero (anything else?)
if (lon == -0f) lon = 0f;
int tmp = Float.floatToIntBits(lat);
int hash = (tmp << 5) | (tmp >> 27);// rotate the lat bits
return hash ^ Float.floatToIntBits(lon);// XOR with lon
}
/**
* Converts an array of decimal degrees float lat/lons to float radians in
* place.
*
* @param degs
* float[] lat/lons in decimal degrees
* @return float[] lat/lons in radians
*/
public static final float[] arrayDegToRad(float[] degs) {
for (int i = 0; i < degs.length; i++) {
degs[i] = degToRad(degs[i]);
}
return degs;
}
/**
* Converts an array of radian float lat/lons to decimal degrees in place.
*
* @param rads
* float[] lat/lons in radians
* @return float[] lat/lons in decimal degrees
*/
public static final float[] arrayRadToDeg(float[] rads) {
for (int i = 0; i < rads.length; i++) {
rads[i] = radToDeg(rads[i]);
}
return rads;
}
/**
* Converts an array of decimal degrees double lat/lons to double radians in
* place.
*
* @param degs
* double[] lat/lons in decimal degrees
* @return double[] lat/lons in radians
*/
public static final double[] arrayDegToRad(double[] degs) {
for (int i = 0; i < degs.length; i++) {
degs[i] = degToRad(degs[i]);
}
return degs;
}
/**
* Converts an array of radian double lat/lons to decimal degrees in place.
*
* @param rads
* double[] lat/lons in radians
* @return double[] lat/lons in decimal degrees
*/
public static final double[] arrayRadToDeg(double[] rads) {
for (int i = 0; i < rads.length; i++) {
rads[i] = radToDeg(rads[i]);
}
return rads;
}
/**
* Normalizes radian latitude. Normalizes latitude if at or exceeds epsilon
* distance from a pole.
*
* @param lat
* float latitude in radians
* @param epsilon
* epsilon (&gt;= 0) radians distance from pole
* @return float latitude (-PI/2 &lt;= phi &lt;= PI/2)
* @see Proj#normalizeLatitude(float)
* @see com.bbn.openmap.LatLonPoint#normalizeLatitude(float)
*/
public static final float normalizeLatitude(float lat, float epsilon) {
if (lat > NORTH_POLE_F - epsilon) {
return NORTH_POLE_F - epsilon;
} else if (lat < SOUTH_POLE_F + epsilon) {
return SOUTH_POLE_F + epsilon;
}
return lat;
}
/**
* Normalizes radian latitude. Normalizes latitude if at or exceeds epsilon
* distance from a pole.
*
* @param lat
* double latitude in radians
* @param epsilon
* epsilon (&gt;= 0) radians distance from pole
* @return double latitude (-PI/2 &lt;= phi &lt;= PI/2)
* @see Proj#normalizeLatitude(float)
*/
public static final double normalizeLatitude(double lat, double epsilon) {
if (lat > NORTH_POLE_D - epsilon) {
return NORTH_POLE_D - epsilon;
} else if (lat < SOUTH_POLE_D + epsilon) {
return SOUTH_POLE_D + epsilon;
}
return lat;
}
/**
* Sets radian longitude to something sane.
*
* @param lon
* float longitude in radians
* @return float longitude (-PI &lt;= lambda &lt; PI)
*/
public static final float wrapLongitude(float lon) {
if ((lon < -DATELINE_F) || (lon > DATELINE_F)) {
lon += DATELINE_F;
lon %= LON_RANGE_F;
lon += (lon < 0) ? DATELINE_F : -DATELINE_F;
}
return lon;
}
/**
* Sets radian longitude to something sane.
*
* @param lon
* double longitude in radians
* @return double longitude (-PI &lt;= lambda &lt; PI)
* @see #wrapLongitude(float)
*/
public static final double wrapLongitude(double lon) {
if ((lon < -DATELINE_D) || (lon > DATELINE_D)) {
lon += DATELINE_D;
lon %= LON_RANGE_D;
lon += (lon < 0) ? DATELINE_D : -DATELINE_D;
}
return lon;
}
/**
* Sets degree longitude to something sane.
*
* @param lon
* double longitude in degrees
* @return double longitude (-180 &lt;= lambda &lt; 180)
*/
public static final double wrapLongitudeDeg(double lon) {
if ((lon < -DATELINE_DEG_D) || (lon > DATELINE_DEG_D)) {
lon += DATELINE_DEG_D;
lon %= LON_RANGE_DEG_D;
lon += (lon < 0) ? DATELINE_DEG_D : -DATELINE_DEG_D;
}
return lon;
}
/**
* Converts units (km, nm, miles, etc) to decimal degrees for a spherical
* planet. This does not check for arc distances &gt; 1/2 planet
* circumference, which are better represented as (2pi - calculated arc).
*
* @param u
* units float value
* @param uCircumference
* units circumference of planet
* @return float decimal degrees
*/
public static final float sphericalUnitsToDeg(float u, float uCircumference) {
return 360f * (u / uCircumference);
}
/**
* Converts units (km, nm, miles, etc) to arc radians for a spherical
* planet. This does not check for arc distances &gt; 1/2 planet
* circumference, which are better represented as (2pi - calculated arc).
*
* @param u
* units float value
* @param uCircumference
* units circumference of planet
* @return float arc radians
*/
public static final float sphericalUnitsToRad(float u, float uCircumference) {
return MoreMath.TWO_PI * (u / uCircumference);
}
/**
* Calculate the geocentric latitude given a geographic latitude. According
* to John Synder: <br>
* "The geographic or geodetic latitude is the angle which a line
* perpendicular to the surface of the ellipsoid at the given point makes
* with the plane of the equator. ...The geocentric latitude is the angle
* made by a line to the center of the ellipsoid with the equatorial plane".
* ( <i>Map Projections --A Working Manual </i>, p 13)
* <p>
* Translated from Ken Anderson's lisp code <i>Freeing the Essence of
* Computation </i>
*
* @param lat
* float geographic latitude in radians
* @param flat
* float flatening factor
* @return float geocentric latitude in radians
* @see #geographic_latitude
*/
public static final float geocentricLatitude(float lat, float flat) {
float f = 1.0f - flat;
return (float) Math.atan((f * f) * (float) Math.tan(lat));
}
/**
* Calculate the geographic latitude given a geocentric latitude. Translated
* from Ken Anderson's lisp code <i>Freeing the Essence of Computation </i>
*
* @param lat
* float geocentric latitude in radians
* @param flat
* float flatening factor
* @return float geographic latitude in radians
* @see #geocentric_latitude
*/
public static final float geographicLatitude(float lat, float flat) {
float f = 1.0f - flat;
return (float) Math.atan((float) Math.tan(lat) / (f * f));
}
/**
* Generic test for seeing if an left longitude value and a right longitude
* value seem to constitute crossing the dateline.
*
* @param leftLon
* the leftmost longitude, in decimal degrees. Expected to
* represent the location of the left side of a map window.
* @param rightLon
* the rightmost longitude, in decimal degrees. Expected to
* represent the location of the right side of a map window.
* @param projScale
* the projection scale, considered if the two values are very
* close to each other and leftLon less than rightLon.
* @return true if it seems like these two longitude values represent a
* dateline crossing.
*/
public static boolean isCrossingDateline(double leftLon, double rightLon, float projScale) {
// if the left longitude is greater than the right, we're obviously
// crossing the dateline. If they are approximately equal, we could be
// showing the whole earth, but only if the scale is significantly
// large. If the scale is small, we could be really zoomed in.
return ((leftLon > rightLon) || (MoreMath.approximately_equal(leftLon, rightLon, .001f) && projScale > 1000000f));
}
}

View file

@ -0,0 +1,461 @@
// **********************************************************************
//
// <copyright>
//
// BBN Technologies
// 10 Moulton Street
// Cambridge, MA 02138
// (617) 873-8000
//
// Copyright (C) BBNT Solutions LLC. All rights reserved.
//
// </copyright>
// **********************************************************************
package com.jwetherell.openmap.common;
public class UTMPoint {
public double northing;
public double easting;
public int zone_number;
public char zone_letter;
public UTMPoint() {
}
public UTMPoint(LatLonPoint llpoint) {
this(llpoint, Ellipsoid.WGS_84);
}
public UTMPoint(LatLonPoint llpoint, Ellipsoid ellip) {
this();
LLtoUTM(llpoint, ellip, this);
}
public UTMPoint(double northing, double easting, int zone_number, char zone_letter) {
this.northing = northing;
this.easting = easting;
this.zone_number = zone_number;
this.zone_letter = checkZone(zone_letter);
}
/**
* Method that provides a check for UTM zone letters. Returns an uppercase
* version of any valid letter passed in, 'N' or 'S'.
*
* @throws NumberFormatException
* if zone letter is invalid.
*/
protected char checkZone(char zone) {
zone = Character.toUpperCase(zone);
if (zone != 'N' && zone != 'S') {
throw new NumberFormatException("Invalid UTMPoint zone letter: " + zone);
}
return zone;
}
public boolean equals(Object obj) {
if (obj instanceof UTMPoint) {
UTMPoint pnt = (UTMPoint) obj;
return northing == pnt.northing && easting == pnt.easting && zone_number == pnt.zone_number && zone_letter == pnt.zone_letter;
}
return false;
}
/**
* Convert this UTMPoint to a LatLonPoint, and assume a WGS_84 ellipsoid.
*/
public LatLonPoint toLatLonPoint() {
return UTMtoLL(this, Ellipsoid.WGS_84, new LatLonPoint());
}
/**
* Convert this UTMPoint to a LatLonPoint, and use the given ellipsoid.
*/
public LatLonPoint toLatLonPoint(Ellipsoid ellip) {
return UTMtoLL(this, ellip, new LatLonPoint());
}
/**
* Fill in the given LatLonPoint with the converted values of this UTMPoint,
* and use the given ellipsoid.
*/
public LatLonPoint toLatLonPoint(Ellipsoid ellip, LatLonPoint llpoint) {
return UTMtoLL(this, ellip, llpoint);
}
/**
* Converts a LatLonPoint to a UTM Point, assuming the WGS_84 ellipsoid.
*
* @return UTMPoint, or null if something bad happened.
*/
public static UTMPoint LLtoUTM(LatLonPoint llpoint) {
return LLtoUTM(llpoint, Ellipsoid.WGS_84, new UTMPoint());
}
/**
* Converts a LatLonPoint to a UTM Point.
*
* @param llpoint
* the LatLonPoint to convert.
* @param utmpoint
* a UTMPoint to put the results in. If it's null, a UTMPoint
* will be allocated.
* @return UTMPoint, or null if something bad happened. If a UTMPoint was
* passed in, it will also be returned on a successful conversion.
*/
public static UTMPoint LLtoUTM(LatLonPoint llpoint, UTMPoint utmpoint) {
return LLtoUTM(llpoint, Ellipsoid.WGS_84, utmpoint);
}
/**
* Converts a set of Longitude and Latitude co-ordinates to UTM given an
* ellipsoid
*
* @param ellip
* an ellipsoid definition.
* @param llpoint
* the coordinate to be converted
* @param utmpoint
* A UTMPoint instance to put the results in. If null, a new
* UTMPoint will be allocated.
* @return A UTM class instance containing the value of <code>null</code> if
* conversion failed. If you pass in a UTMPoint, it will be returned
* as well if successful.
*/
public static UTMPoint LLtoUTM(LatLonPoint llpoint, Ellipsoid ellip, UTMPoint utmpoint) {
// find the native zone for the given latitude/longitude point
int zoneNumber = getZoneNumber(llpoint.getY(), llpoint.getX());
boolean isnorthern = (llpoint.getLatitude() >= 0f);
return LLtoUTM(llpoint, ellip, utmpoint, zoneNumber, isnorthern);
}
/**
* Converts a set of Longitude and Latitude co-ordinates to UTM given an
* ellipsoid and the UTM zone to use.
*
* @param ellip
* an ellipsoid definition.
* @param llpoint
* the coordinate to be converted
* @param utmPoint
* A UTMPoint instance to put the results in. If null, a new
* UTMPoint will be allocated.
* @param zoneNumber
* the number of the zone
* @param isNorthern
* true if zone is in northern hemispehere
* @return A UTM class instance containing the value of <code>null</code> if
* conversion failed. If you pass in a UTMPoint, it will be returned
* as well if successful.
*/
public static UTMPoint LLtoUTM(LatLonPoint llpoint, Ellipsoid ellip, UTMPoint utmPoint, int zoneNumber, boolean isNorthern) {
double a = ellip.radius;
double k0 = 0.9996;
double eccSquared = ellip.eccsq;
double eccPrimeSquared = (eccSquared) / (1 - eccSquared);
double eccSquared2 = eccSquared * eccSquared;
double eccSquared3 = eccSquared2 * eccSquared;
double N, T, C, A, M;
double LatRad = llpoint.getRadLat();
double LongRad = llpoint.getRadLon();
// in middle of zone
double LongOrigin = (zoneNumber - 1) * 6 - 180 + 3; // +3 puts origin
double LongOriginRad = Math.toRadians(LongOrigin);
double tanLatRad = Math.tan(LatRad);
double sinLatRad = Math.sin(LatRad);
double cosLatRad = Math.cos(LatRad);
N = a / Math.sqrt(1 - eccSquared * sinLatRad * sinLatRad);
T = tanLatRad * tanLatRad;
C = eccPrimeSquared * cosLatRad * cosLatRad;
A = cosLatRad * (LongRad - LongOriginRad);
M = a
* ((1 - eccSquared / 4 - 3 * eccSquared2 / 64 - 5 * eccSquared3 / 256) * LatRad
- (3 * eccSquared / 8 + 3 * eccSquared2 / 32 + 45 * eccSquared3 / 1024) * Math.sin(2 * LatRad)
+ (15 * eccSquared2 / 256 + 45 * eccSquared3 / 1024) * Math.sin(4 * LatRad) - (35 * eccSquared3 / 3072) * Math.sin(6 * LatRad));
double UTMEasting = (k0 * N * (A + (1 - T + C) * A * A * A / 6.0d + (5 - 18 * T + T * T + 72 * C - 58 * eccPrimeSquared) * A * A * A * A * A / 120.0d) + 500000.0d);
double UTMNorthing = (k0 * (M + N
* Math.tan(LatRad)
* (A * A / 2 + (5 - T + 9 * C + 4 * C * C) * A * A * A * A / 24.0d + (61 - 58 * T + T * T + 600 * C - 330 * eccPrimeSquared) * A * A * A * A
* A * A / 720.0d)));
if (!isNorthern) {
UTMNorthing += 10000000.0f; // 10000000 meter offset for
// southern hemisphere
}
if (utmPoint == null) {
utmPoint = new UTMPoint();
}
utmPoint.northing = UTMNorthing;
utmPoint.easting = UTMEasting;
utmPoint.zone_number = zoneNumber;
utmPoint.zone_letter = isNorthern ? 'N' : 'S';
return utmPoint;
}
/**
* Returns 'N' if the latitude is equal to or above the equator, 'S' if it's
* below.
*
* @param lat
* The float value of the latitude.
*
* @return A char value
*/
protected char getLetterDesignator(double lat) {
char letterDesignator = 'N';
if (lat < 0) {
letterDesignator = 'S';
}
return letterDesignator;
}
/**
* Converts UTM coords to lat/long given an ellipsoid given an instance of
* UTMPoint.
*
* @param utm_point
* A UTMPoint instance.
* @param ellip
* a ellipsoid definition.
* @param llpoint
* a LatLonPoint, if you want it to be filled in with the
* results. If null, a new LatLonPoint will be allocated.
* @return A LatLonPoint class instance containing the lat/long value, or
* <code>null</code> if conversion failed. If you pass in a
* LatLonPoint, it will be returned as well, if successful.
*/
public static LatLonPoint UTMtoLL(UTMPoint utm_point, Ellipsoid ellip, LatLonPoint llpoint) {
return UTMtoLL(ellip, utm_point.northing, utm_point.easting, utm_point.zone_number, utm_point.zone_letter, llpoint);
}
/**
* Converts UTM coords to lat/long given an ellipsoid. This is a convenience
* class where the Zone can be specified as a single string eg."61N" which
* is then broken down into the ZoneNumber and ZoneLetter.
*
* @param ellip
* an ellipsoid definition.
* @param UTMNorthing
* A float value for the northing to be converted.
* @param UTMEasting
* A float value for the easting to be converted.
* @param UTMZone
* A String value for the UTM zone eg."61N".
* @param llpoint
* a LatLonPoint, if you want it to be filled in with the
* results. If null, a new LatLonPoint will be allocated.
* @return A LatLonPoint class instance containing the lat/long value, or
* <code>null</code> if conversion failed. If you pass in a
* LatLonPoint, it will be returned as well, if successful.
*/
public static LatLonPoint UTMtoLL(Ellipsoid ellip, double UTMNorthing, double UTMEasting, String UTMZone, LatLonPoint llpoint) {
// without the zone we can't calculate the Lat and Long
if (UTMZone == null || UTMZone.length() == 0) {
return null;
}
int ZoneNumber = 1;
char ZoneLetter = 'N'; // northern hemisphere by default if no
// character is found
// Break out the Zone number and zone letter from the UTMZone
// string We assume the string is a valid zone with a number
// followed by a zone letter If there is no Letter we assume
// that it's the Northern hemisphere
int ln = UTMZone.length() - 1;
if (ln > 0) {
// If it's Zero then there is only one character and it
// must be the Zone number
ZoneLetter = UTMZone.charAt(ln);
if (!Character.isLetter(ZoneLetter)) {
// No letter so assume it's missing & default to 'N'
ZoneLetter = 'N';
ln++;
}
}
// convert the number but catch the exception if it's not
// valid
try {
ZoneNumber = Integer.parseInt(UTMZone.substring(0, ln));
} catch (NumberFormatException nfe) {
return null;
}
return UTMtoLL(ellip, UTMNorthing, UTMEasting, ZoneNumber, ZoneLetter, llpoint);
}
/**
* Converts UTM coords to lat/long given an ellipsoid. This is a convenience
* class where the exact Zone letter is not known. Instead only the
* hemisphere needs to be indicated.
*
* @param ellip
* an ellipsoid definition.
* @param UTMNorthing
* A float value for the northing to be converted.
* @param UTMEasting
* A float value for the easting to be converted.
* @param ZoneNumber
* An int value indicating the float number.
* @param isNorthern
* A boolean which is true for the northern hemisphere otherwise
* false for the southern.
* @param llpoint
* a LatLonPoint, if you want it to be filled in with the
* results. If null, a new LatLonPoint will be allocated.
* @return A LatLonPoint class instance containing the lat/long value, or
* <code>null</code> if conversion failed. If you pass in a
* LatLonPoint, it will be returned as well, if successful.
*/
public static LatLonPoint UTMtoLL(Ellipsoid ellip, double UTMNorthing, double UTMEasting, int ZoneNumber, boolean isNorthern, LatLonPoint llpoint) {
return UTMtoLL(ellip, UTMNorthing, UTMEasting, ZoneNumber, (isNorthern) ? 'N' : 'S', llpoint);
}
/**
* Converts UTM coords to lat/long given an ellipsoid.
* <p>
* Equations from USGS Bulletin 1532 <br>
* East Longitudes are positive, West longitudes are negative. <br>
* North latitudes are positive, South latitudes are negative. <br>
*
* @param ellip
* an ellipsoid definition.
* @param UTMNorthing
* A float value for the northing to be converted.
* @param UTMEasting
* A float value for the easting to be converted.
* @param zoneNumber
* An int value specifiying the UTM zone number.
* @param zoneLetter
* A char value specifying the ZoneLetter within the ZoneNumber.
* @param llpoint
* a LatLonPoint, if you want it to be filled in with the
* results. If null, a new LatLonPoint will be allocated.
* @return A LatLonPoint class instance containing the lat/long value, or
* <code>null</code> if conversion failed. If you pass in a
* LatLonPoint, it will be returned as well, if successful.
*/
public static LatLonPoint UTMtoLL(Ellipsoid ellip, double UTMNorthing, double UTMEasting, int zoneNumber, char zoneLetter, LatLonPoint llpoint) {
// check the ZoneNummber is valid
if (zoneNumber < 0 || zoneNumber > 60) {
return null;
}
double k0 = 0.9996;
double a = ellip.radius;
double eccSquared = ellip.eccsq;
double eccPrimeSquared;
double e1 = (1 - Math.sqrt(1 - eccSquared)) / (1 + Math.sqrt(1 - eccSquared));
double N1, T1, C1, R1, D, M;
double LongOrigin;
double mu, phi1Rad;
// remove 500,000 meter offset for longitude
double x = UTMEasting - 500000.0d;
double y = UTMNorthing;
// We must know somehow if we are in the Northern or Southern
// hemisphere, this is the only time we use the letter So even
// if the Zone letter isn't exactly correct it should indicate
// the hemisphere correctly
if (zoneLetter == 'S') {
y -= 10000000.0d;// remove 10,000,000 meter offset used
// for southern hemisphere
}
// There are 60 zones with zone 1 being at West -180 to -174
LongOrigin = (zoneNumber - 1) * 6 - 180 + 3; // +3 puts origin
// in middle of
// zone
eccPrimeSquared = (eccSquared) / (1 - eccSquared);
M = y / k0;
mu = M / (a * (1 - eccSquared / 4 - 3 * eccSquared * eccSquared / 64 - 5 * eccSquared * eccSquared * eccSquared / 256));
phi1Rad = mu + (3 * e1 / 2 - 27 * e1 * e1 * e1 / 32) * Math.sin(2 * mu) + (21 * e1 * e1 / 16 - 55 * e1 * e1 * e1 * e1 / 32) * Math.sin(4 * mu)
+ (151 * e1 * e1 * e1 / 96) * Math.sin(6 * mu);
// double phi1 = ProjMath.radToDeg(phi1Rad);
N1 = a / Math.sqrt(1 - eccSquared * Math.sin(phi1Rad) * Math.sin(phi1Rad));
T1 = Math.tan(phi1Rad) * Math.tan(phi1Rad);
C1 = eccPrimeSquared * Math.cos(phi1Rad) * Math.cos(phi1Rad);
R1 = a * (1 - eccSquared) / Math.pow(1 - eccSquared * Math.sin(phi1Rad) * Math.sin(phi1Rad), 1.5);
D = x / (N1 * k0);
double lat = phi1Rad
- (N1 * Math.tan(phi1Rad) / R1)
* (D * D / 2 - (5 + 3 * T1 + 10 * C1 - 4 * C1 * C1 - 9 * eccPrimeSquared) * D * D * D * D / 24 + (61 + 90 * T1 + 298 * C1 + 45 * T1 * T1 - 252
* eccPrimeSquared - 3 * C1 * C1)
* D * D * D * D * D * D / 720);
lat = ProjMath.radToDeg(lat);
double lon = (D - (1 + 2 * T1 + C1) * D * D * D / 6 + (5 - 2 * C1 + 28 * T1 - 3 * C1 * C1 + 8 * eccPrimeSquared + 24 * T1 * T1) * D * D * D * D * D
/ 120)
/ Math.cos(phi1Rad);
lon = LongOrigin + ProjMath.radToDeg(lon);
if (llpoint != null) {
llpoint.setLatLon(lat, lon);
return llpoint;
}
return new LatLonPoint(lat, lon);
}
/**
* Find zone number based on the given latitude and longitude in *degrees*.
*
* @param lat
* @param lon
* @return
*/
private static int getZoneNumber(double lat, double lon) {
int zoneNumber = (int) ((lon + 180) / 6) + 1;
// Make sure the longitude 180.00 is in Zone 60
if (lon == 180) {
zoneNumber = 60;
}
// Special zone for Norway
if (lat >= 56.0f && lat < 64.0f && lon >= 3.0f && lon < 12.0f) {
zoneNumber = 32;
}
// Special zones for Svalbard
if (lat >= 72.0f && lat < 84.0f) {
if (lon >= 0.0f && lon < 9.0f) zoneNumber = 31;
else if (lon >= 9.0f && lon < 21.0f) zoneNumber = 33;
else if (lon >= 21.0f && lon < 33.0f) zoneNumber = 35;
else if (lon >= 33.0f && lon < 42.0f) zoneNumber = 37;
}
return zoneNumber;
}
/**
* {@inheritDoc}
*/
@Override
public String toString() {
return "Zone_number=" + zone_number + ", Hemisphere=" + zone_letter + ", Northing=" + northing + ", Easting=" + easting;
}
}

View file

@ -0,0 +1,111 @@
// **********************************************************************
//
// <copyright>
//
// BBN Technologies
// 10 Moulton Street
// Cambridge, MA 02138
// (617) 873-8000
//
// Copyright (C) BBNT Solutions LLC. All rights reserved.
//
// </copyright>
// **********************************************************************
package com.jwetherell.openmap.common;
public class ZonedUTMPoint extends UTMPoint {
public ZonedUTMPoint() {
super();
}
public ZonedUTMPoint(LatLonPoint llpoint) {
super(llpoint);
zone_letter = getLetterDesignator(llpoint.getLatitude());
}
public ZonedUTMPoint(LatLonPoint llpoint, Ellipsoid ellip) {
super(llpoint, ellip);
zone_letter = getLetterDesignator(llpoint.getLatitude());
}
public ZonedUTMPoint(double northing, double easting, int zone_number, char zone_letter) {
super(northing, easting, zone_number, MGRSPoint.MGRSZoneToUTMZone(zone_letter));
this.zone_letter = zone_letter;
}
/**
* Converts UTM coords to lat/long given an ellipsoid.
* <p>
* Equations from USGS Bulletin 1532 <br>
* East Longitudes are positive, West longitudes are negative. <br>
* North latitudes are positive, South latitudes are negative. <br>
*
* @param ellip
* an ellipsoid definition.
* @param UTMNorthing
* A float value for the northing to be converted.
* @param UTMEasting
* A float value for the easting to be converted.
* @param ZoneNumber
* An int value specifiying the UTM zone number.
* @param ZoneLetter
* A char value specifying the ZoneLetter within the ZoneNumber,
* letter being MGRS zone.
* @param llpoint
* a LatLonPoint, if you want it to be filled in with the
* results. If null, a new LatLonPoint will be allocated.
* @return A LatLonPoint class instance containing the lat/long value, or
* <code>null</code> if conversion failed. If you pass in a
* LatLonPoint, it will be returned as well, if successful.
*/
public static LatLonPoint ZonedUTMtoLL(Ellipsoid ellip, double UTMNorthing, double UTMEasting, int ZoneNumber, char ZoneLetter, LatLonPoint llpoint) {
return UTMPoint.UTMtoLL(ellip, UTMNorthing, UTMEasting, ZoneNumber, MGRSPoint.MGRSZoneToUTMZone(ZoneLetter), llpoint);
}
/**
* Determines the correct MGRS letter designator for the given latitude
* returns 'Z' if latitude is outside the MGRS limits of 84N to 80S.
*
* @param lat
* The float value of the latitude.
*
* @return A char value which is the MGRS zone letter.
*/
protected char getLetterDesignator(double lat) {
// This is here as an error flag to show that the Latitude is
// outside MGRS limits
char LetterDesignator = 'Z';
if ((84 >= lat) && (lat >= 72)) LetterDesignator = 'X';
else if ((72 > lat) && (lat >= 64)) LetterDesignator = 'W';
else if ((64 > lat) && (lat >= 56)) LetterDesignator = 'V';
else if ((56 > lat) && (lat >= 48)) LetterDesignator = 'U';
else if ((48 > lat) && (lat >= 40)) LetterDesignator = 'T';
else if ((40 > lat) && (lat >= 32)) LetterDesignator = 'S';
else if ((32 > lat) && (lat >= 24)) LetterDesignator = 'R';
else if ((24 > lat) && (lat >= 16)) LetterDesignator = 'Q';
else if ((16 > lat) && (lat >= 8)) LetterDesignator = 'P';
else if ((8 > lat) && (lat >= 0)) LetterDesignator = 'N';
else if ((0 > lat) && (lat >= -8)) LetterDesignator = 'M';
else if ((-8 > lat) && (lat >= -16)) LetterDesignator = 'L';
else if ((-16 > lat) && (lat >= -24)) LetterDesignator = 'K';
else if ((-24 > lat) && (lat >= -32)) LetterDesignator = 'J';
else if ((-32 > lat) && (lat >= -40)) LetterDesignator = 'H';
else if ((-40 > lat) && (lat >= -48)) LetterDesignator = 'G';
else if ((-48 > lat) && (lat >= -56)) LetterDesignator = 'F';
else if ((-56 > lat) && (lat >= -64)) LetterDesignator = 'E';
else if ((-64 > lat) && (lat >= -72)) LetterDesignator = 'D';
else if ((-72 > lat) && (lat >= -80)) LetterDesignator = 'C';
return LetterDesignator;
}
/**
* {@inheritDoc}
*/
@Override
public String toString() {
return "Zone_number=" + zone_number + ", Hemisphere=" + zone_letter + ", Northing=" + northing + ", Easting=" + easting;
}
}