001: /*
002: * The Unified Mapping Platform (JUMP) is an extensible, interactive GUI
003: * for visualizing and manipulating spatial features with geometry and attributes.
004: *
005: * Copyright (C) 2003 Vivid Solutions
006: *
007: * This program is free software; you can redistribute it and/or
008: * modify it under the terms of the GNU General Public License
009: * as published by the Free Software Foundation; either version 2
010: * of the License, or (at your option) any later version.
011: *
012: * This program is distributed in the hope that it will be useful,
013: * but WITHOUT ANY WARRANTY; without even the implied warranty of
014: * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
015: * GNU General Public License for more details.
016: *
017: * You should have received a copy of the GNU General Public License
018: * along with this program; if not, write to the Free Software
019: * Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
020: *
021: * For more information, contact:
022: *
023: * Vivid Solutions
024: * Suite #1A
025: * 2328 Government Street
026: * Victoria BC V8T 5G5
027: * Canada
028: *
029: * (250)385-6040
030: * www.vividsolutions.com
031: */
032:
033: package com.vividsolutions.jump.coordsys.impl;
034:
035: import com.vividsolutions.jump.coordsys.*;
036: import com.vividsolutions.jump.coordsys.Geographic;
037: import com.vividsolutions.jump.coordsys.Planar;
038:
039: /**
040: * Implements the Polyconic projection.
041: * *
042: * @author $Author: javamap $
043: * @version $Revision: 4 $
044: * <pre>
045: * $Id: Polyconic.java 4 2005-06-16 15:27:48Z javamap $
046: * $Date: 2005-06-16 08:27:48 -0700 (Thu, 16 Jun 2005) $
047: *
048: *
049: * $Log$
050: * Revision 1.1 2005/06/16 15:25:29 javamap
051: * *** empty log message ***
052: *
053: * Revision 1.2 2005/05/03 15:23:55 javamap
054: * *** empty log message ***
055: *
056: * Revision 1.2 2003/11/05 05:14:18 dkim
057: * Added global header; cleaned up Javadoc.
058: *
059: * Revision 1.1 2003/09/15 20:26:12 jaquino
060: * Reprojection
061: *
062: * Revision 1.2 2003/07/25 17:01:04 gkostadinov
063: * Moved classses reponsible for performing the basic projection to a new
064: * package -- base.
065: *
066: * Revision 1.1 2003/07/24 23:14:44 gkostadinov
067: * adding base projection classes
068: *
069: * Revision 1.1 2003/06/20 18:34:31 gkostadinov
070: * Entering the source code into the CVS.
071: *
072: * </pre>
073: *
074: */
075:
076: public class Polyconic extends Projection {
077:
078: // private:
079: double L0;// central meridian
080: double k0;// scale factor
081: double phi1;// 1st standard parallel
082: double phi2;// 2nd standard parallel
083: double phi0;// Latitude of projection
084: double X0;// false Easting
085: double Y0;// false Northing
086: int zone;// UTMzone
087: MeridianArcLength S = new MeridianArcLength();
088: Geographic q = new Geographic();
089:
090: public Polyconic() {
091: super ();
092: }
093:
094: public void setParameters(double originLatitude,
095: double originLongitude) {
096: // Polyconic projection
097: L0 = originLongitude / 180.0 * Math.PI;
098: phi0 = originLatitude / 180.0 * Math.PI;
099: }
100:
101: public Planar asPlanar(Geographic q0, Planar p) {
102: q.lat = q0.lat / 180.0 * Math.PI;
103: q.lon = q0.lon / 180.0 * Math.PI;
104: forward(q, p);
105: return p;
106: }
107:
108: public Geographic asGeographic(Planar p, Geographic q) {
109: inverse(p, q);
110: q.lat = q.lat * 180.0 / Math.PI;
111: q.lon = q.lon * 180.0 / Math.PI;
112: return q;
113: }
114:
115: public void forward(Geographic q, Planar p) {
116: double M;
117: double M0;
118: S.compute(currentSpheroid, q.lat, 0);
119: M = S.s;
120: S.compute(currentSpheroid, phi0, 0);
121: M0 = S.s;
122: double a;
123: double e;
124: double e2;
125: a = currentSpheroid.a;
126: e = currentSpheroid.e;
127: e2 = e * e;
128: double N;
129: double t;
130: t = Math.sin(q.lat);
131: N = a / Math.sqrt(1.0 - e2 * t * t);
132: double E;
133: E = (q.lon - L0) * Math.sin(q.lat);
134: t = 1.0 / Math.tan(q.lat);
135: p.x = N * t * Math.sin(E);
136: p.y = M - M0 + N * t * (1.0 - Math.cos(E));
137: }
138:
139: public void inverse(Planar p, Geographic q) {
140: double a;
141: double e;
142: double es;
143: a = currentSpheroid.getA();
144: e = currentSpheroid.getE();
145: es = e * e;
146: double A;
147: double B;
148: double M0;
149: S.compute(currentSpheroid, phi0, 0);
150: M0 = S.s;
151: A = (M0 + p.y) / a;
152: B = (p.x * p.x) / (a * a) + A * A;
153: double C;
154: double phiN;
155: double M;
156: double Mp;
157: double Ma;
158: double Ma2;
159: double s2p;
160: q.lat = A;
161: int count = 0;
162: do {
163: phiN = q.lat;
164: C = Math.sqrt(1.0 - es * Math.sin(phiN) * Math.sin(phiN))
165: * Math.tan(phiN);
166: S.compute(currentSpheroid, phiN, 0);
167: M = S.s;
168: Ma = M / a;
169: Ma2 = Ma * Ma;
170: S.compute(currentSpheroid, phiN, 1);
171: Mp = S.s;
172: s2p = Math.sin(2.0 * phiN);
173: q.lat = q.lat
174: - (A * (C * Ma + 1.0) - Ma - 0.5 * (Ma2 + B) * C)
175: / (es * s2p * (Ma2 + B - 2.0 * A * Ma) / 4.0 * C
176: + (A - Ma) * (C * Mp - 2.0 / s2p) - Mp);
177: } while (Math.abs(q.lat - phiN) > 1.0e-6 && count++ < 100);//1.0e-12);
178: q.lon = Math.asin(p.x * C / a) / Math.sin(q.lat) + L0;
179: }
180:
181: }
|