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 Albers projection.
041: *
042: *
043: * @author $Author: javamap $
044: * @version $Revision: 4 $
045: *
046: * <pre>
047: * $Id: Albers.java 4 2005-06-16 15:27:48Z javamap $
048: * $Date: 2005-06-16 08:27:48 -0700 (Thu, 16 Jun 2005) $
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:10:31 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: * </pre>
072: *
073: */
074:
075: public class Albers extends Projection {
076:
077: double L0;// central meridian
078: double k0;// scale factor
079: double phi1;// 1st standard parallel
080: double phi2;// 2nd standard parallel
081: double phi0;// Latitude of projection
082: double X0;// false Easting
083: double Y0;// false Northing
084: double A_n, A_C, A_p0;// variables for Albers, see constructor
085: Geographic q = new Geographic();
086:
087: public Albers() {
088: super ();
089: }
090:
091: public void setParameters(double centralMeridian,
092: double firstStandardParallel,
093: double secondStandardParallel, double latitudeOfProjection,
094: double falseEasting, double falseNorthing) {
095: L0 = centralMeridian * Math.PI / 180.0;
096: phi1 = firstStandardParallel * Math.PI / 180.0;
097: phi2 = secondStandardParallel * Math.PI / 180.0;
098: phi0 = latitudeOfProjection * Math.PI / 180.0;
099: X0 = falseEasting;
100: Y0 = falseNorthing;
101: double m1;
102: double m2;
103: double q1;
104: double q2;
105: double q0;
106: m1 = albersM(phi1);
107: m2 = albersM(phi2);
108: q1 = albersQ(phi1);
109: q2 = albersQ(phi2);
110: q0 = albersQ(phi0);
111: A_n = (m1 * m1 - m2 * m2) / (q2 - q1);
112: A_C = m1 * m1 + A_n * q1;
113: double a;
114: a = currentSpheroid.getA();
115: A_p0 = (a * Math.sqrt(A_C - A_n * q0)) / A_n;
116: }// END - constructor for Albers projection plane
117:
118: public Planar asPlanar(Geographic q0, Planar p) {
119: q.lat = q0.lat / 180.0 * Math.PI;
120: q.lon = q0.lon / 180.0 * Math.PI;
121: forward(q, p);
122: return p;
123: }
124:
125: public Geographic asGeographic(Planar p, Geographic q) {
126: inverse(p, q);
127: q.lat = q.lat * 180.0 / Math.PI;
128: q.lon = q.lon * 180.0 / Math.PI;
129: return q;
130: }
131:
132: void forward(Geographic q, Planar p) {
133: double que;
134: double theta;
135: double pee;
136: double a;
137: a = currentSpheroid.getA();
138: que = albersQ(q.lat);
139: theta = A_n * (q.lon - L0);
140: pee = (a * Math.sqrt(A_C - A_n * que)) / A_n;
141: p.x = pee * Math.sin(theta) + X0;
142: p.y = A_p0 - pee * Math.cos(theta) + Y0;
143: }
144:
145: void inverse(Planar p, Geographic q) {
146: double es;
147: double x;
148: double y;
149: double theta;
150: double pee;
151: double que;
152: double a;
153: double e;
154: a = currentSpheroid.getA();
155: e = currentSpheroid.getE();
156: es = e * e;
157: x = p.x - X0;
158: y = p.y - Y0;
159: theta = Math.atan2(x, A_p0 - y);
160: pee = Math.sqrt(x * x + Math.pow((A_p0 - y), 2.0));
161: que = (A_C - (pee * pee * A_n * A_n) / (a * a)) / A_n;
162: q.lon = L0 + theta / A_n;
163: double li;
164: double delta;
165: double j1;
166: double k1;
167: double k2;
168: double k3;
169: double lip1;
170: //li = Math.sin(que / 2.0); -- transcription error
171: li = Math.asin(que / 2.0);
172: delta = 10e010;
173: do {
174: j1 = Math
175: .pow((1.0 - es * Math.pow(Math.sin(li), 2.0)), 2.0)
176: / (2.0 * Math.cos(li));
177: k1 = que / (1.0 - es);
178: k2 = Math.sin(li)
179: / (1.0 - es * Math.pow(Math.sin(li), 2.0));
180: k3 = (1.0 / (2.0 * e))
181: * Math.log((1.0 - e * Math.sin(li))
182: / (1.0 + e * Math.sin(li)));
183: lip1 = li + j1 * (k1 - k2 + k3);
184: delta = Math.abs(lip1 - li);
185: li = lip1;
186: } while (delta > 1.0e-012);
187: q.lat = li;
188: }
189:
190: double albersQ(double lat) {
191: double q;
192: double e;
193: e = currentSpheroid.getE();
194: q = (1.0 - e * e)
195: * (Math.sin(lat)
196: / (1.0 - e * e * Math.pow(Math.sin(lat), 2.0)) - (1.0 / (2.0 * e))
197: * Math.log((1.0 - e * Math.sin(lat))
198: / (1.0 + e * Math.sin(lat))));
199: return q;
200: }
201:
202: double albersM(double lat) {
203: double m;
204: double e;
205: e = currentSpheroid.getE();
206: m = Math.cos(lat)
207: / Math.sqrt(1.0 - e * e * Math.pow(Math.sin(lat), 2.0));
208: return m;
209: }
210:
211: }
|