svn-gvsig-desktop / branches / F2 / libraries / libJCRS / src / org / geotools / referencing / operation / projection / Mercator.java @ 12260
History | View | Annotate | Download (17.9 KB)
1 |
/*
|
---|---|
2 |
* Geotools 2 - OpenSource mapping toolkit
|
3 |
* (C) 2003, Geotools Project Managment Committee (PMC)
|
4 |
* (C) 2001, Institut de Recherche pour le Dveloppement
|
5 |
* (C) 1999, Fisheries and Oceans Canada
|
6 |
*
|
7 |
* This library is free software; you can redistribute it and/or
|
8 |
* modify it under the terms of the GNU Lesser General Public
|
9 |
* License as published by the Free Software Foundation; either
|
10 |
* version 2.1 of the License, or (at your option) any later version.
|
11 |
*
|
12 |
* This library is distributed in the hope that it will be useful,
|
13 |
* but WITHOUT ANY WARRANTY; without even the implied warranty of
|
14 |
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
|
15 |
* Lesser General Public License for more details.
|
16 |
*
|
17 |
* You should have received a copy of the GNU Lesser General Public
|
18 |
* License along with this library; if not, write to the Free Software
|
19 |
* Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
|
20 |
*
|
21 |
*
|
22 |
* This package contains formulas from the PROJ package of USGS.
|
23 |
* USGS's work is fully acknowledged here.
|
24 |
*/
|
25 |
package org.geotools.referencing.operation.projection; |
26 |
|
27 |
// J2SE dependencies and extensions
|
28 |
import java.awt.geom.Point2D; |
29 |
import java.util.Collection; |
30 |
import javax.units.NonSI; |
31 |
|
32 |
// OpenGIS dependencies
|
33 |
import org.opengis.parameter.ParameterDescriptor; |
34 |
import org.opengis.parameter.ParameterDescriptorGroup; |
35 |
import org.opengis.parameter.ParameterNotFoundException; |
36 |
import org.opengis.parameter.ParameterValueGroup; |
37 |
import org.opengis.referencing.operation.CylindricalProjection; |
38 |
import org.opengis.referencing.operation.MathTransform; |
39 |
|
40 |
// Geotools dependencies
|
41 |
import org.geotools.measure.Latitude; |
42 |
import org.geotools.metadata.iso.citation.CitationImpl; |
43 |
import org.geotools.referencing.NamedIdentifier; |
44 |
import org.geotools.referencing.operation.MathTransformProvider; |
45 |
import org.geotools.resources.cts.ResourceKeys; |
46 |
import org.geotools.resources.cts.Resources; |
47 |
|
48 |
|
49 |
/**
|
50 |
* Mercator Cylindrical Projection. The parallels and the meridians are straight lines and
|
51 |
* cross at right angles; this projection thus produces rectangular charts. The scale is true
|
52 |
* along the equator (by default) or along two parallels equidistant of the equator (if a scale
|
53 |
* factor other than 1 is used). This projection is used to represent areas close to the equator.
|
54 |
* It is also often used for maritime navigation because all the straight lines on the chart are
|
55 |
* <em>loxodrome</em> lines, i.e. a ship following this line would keep a constant azimuth on its
|
56 |
* compass.
|
57 |
* <br><br>
|
58 |
*
|
59 |
* This implementation handles both the 1 and 2 stardard parallel cases.
|
60 |
* For <code>Mercator_1SP</code> (EPSG code 9804), the line of contact is the equator.
|
61 |
* For <code>Mercator_2SP</code> (EPSG code 9805) lines of contact are symmetrical
|
62 |
* about the equator.
|
63 |
* <br><br>
|
64 |
*
|
65 |
* <strong>References:</strong><ul>
|
66 |
* <li>John P. Snyder (Map Projections - A Working Manual,<br>
|
67 |
* U.S. Geological Survey Professional Paper 1395, 1987)</li>
|
68 |
* <li>"Coordinate Conversions and Transformations including Formulas",<br>
|
69 |
* EPSG Guidence Note Number 7, Version 19.</li>
|
70 |
* </ul>
|
71 |
*
|
72 |
* @see <A HREF="http://mathworld.wolfram.com/MercatorProjection.html">Mercator projection on MathWorld</A>
|
73 |
* @see <A HREF="http://www.remotesensing.org/geotiff/proj_list/mercator_1sp.html">"mercator_1sp" on Remote Sensing</A>
|
74 |
* @see <A HREF="http://www.remotesensing.org/geotiff/proj_list/mercator_2sp.html">"mercator_2sp" on Remote Sensing</A>
|
75 |
*
|
76 |
* @version $Id: Mercator.java 12260 2007-06-21 09:48:20Z jlgomez $
|
77 |
* @author Andr Gosselin
|
78 |
* @author Martin Desruisseaux
|
79 |
* @author Rueben Schulz
|
80 |
*/
|
81 |
public class Mercator extends MapProjection { |
82 |
/**
|
83 |
* Standard Parallel used for the <code>Mercator_2SP</code> case.
|
84 |
* Set to {@link Double#NaN} for the <code>Mercator_1SP</code> case.
|
85 |
*/
|
86 |
protected final double standardParallel; |
87 |
//protected final double latitudeOfOrigin;
|
88 |
|
89 |
/**
|
90 |
* The {@link MathTransformProvider} for a {@link Mercator} 1SP projection.
|
91 |
*
|
92 |
* @see <A HREF="http://www.remotesensing.org/geotiff/proj_list/mercator_1sp.html">"mercator_1sp" on Remote Sensing</A>
|
93 |
* @see org.geotools.referencing.operation.DefaultMathTransformFactory
|
94 |
*
|
95 |
* @version $Id: Mercator.java 12260 2007-06-21 09:48:20Z jlgomez $
|
96 |
* @author Martin Desruisseaux
|
97 |
* @author Rueben Schulz
|
98 |
*/
|
99 |
public static final class Provider1SP extends AbstractProvider { |
100 |
/**
|
101 |
* The parameters group.
|
102 |
*/
|
103 |
public static final ParameterDescriptor LATITUDE_OF_ORIGIN = createDescriptor( |
104 |
new NamedIdentifier[] { |
105 |
new NamedIdentifier(CitationImpl.OGC, "latitude_of_origin"), |
106 |
new NamedIdentifier(CitationImpl.EPSG, "CenterLat"), |
107 |
new NamedIdentifier(CitationImpl.EPSG, "Latitude of projection centre"), |
108 |
new NamedIdentifier(CitationImpl.GEOTIFF, "NatOriginLat"), |
109 |
new NamedIdentifier(CitationImpl.EPSG, "FalseOriginLat"), |
110 |
new NamedIdentifier(CitationImpl.EPSG, "Latitude of false origin"), |
111 |
new NamedIdentifier(CitationImpl.EPSG, "Latitude of natural origin"), |
112 |
new NamedIdentifier(CitationImpl.EPSG, "Latitude of projection centre"), |
113 |
new NamedIdentifier(CitationImpl.EPSG, "ProjCenterLat") |
114 |
}, 0.0, -90.0, 90.0, NonSI.DEGREE_ANGLE); |
115 |
|
116 |
static final ParameterDescriptorGroup PARAMETERS = createDescriptorGroup(new NamedIdentifier[] { |
117 |
new NamedIdentifier(CitationImpl.OGC, "Mercator_1SP"), |
118 |
new NamedIdentifier(CitationImpl.EPSG, "Mercator (1SP)"), |
119 |
new NamedIdentifier(CitationImpl.EPSG, "9804"), |
120 |
new NamedIdentifier(CitationImpl.GEOTIFF, "CT_Mercator"), |
121 |
new NamedIdentifier(CitationImpl.GEOTOOLS, Resources.formatInternational(
|
122 |
ResourceKeys.CYLINDRICAL_MERCATOR_PROJECTION)) |
123 |
}, new ParameterDescriptor[] { |
124 |
SEMI_MAJOR, SEMI_MINOR, |
125 |
CENTRAL_MERIDIAN, SCALE_FACTOR, |
126 |
LATITUDE_OF_ORIGIN, |
127 |
FALSE_EASTING, FALSE_NORTHING |
128 |
}); |
129 |
|
130 |
/**
|
131 |
* Constructs a new provider.
|
132 |
*/
|
133 |
public Provider1SP() {
|
134 |
super(PARAMETERS);
|
135 |
} |
136 |
|
137 |
/**
|
138 |
* Returns the operation type for this map projection.
|
139 |
*/
|
140 |
protected Class getOperationType() { |
141 |
return CylindricalProjection.class;
|
142 |
} |
143 |
|
144 |
/**
|
145 |
* Creates a transform from the specified group of parameter values.
|
146 |
*
|
147 |
* @param parameters The group of parameter values.
|
148 |
* @return The created math transform.
|
149 |
* @throws ParameterNotFoundException if a required parameter was not found.
|
150 |
*/
|
151 |
public MathTransform createMathTransform(final ParameterValueGroup parameters) |
152 |
throws ParameterNotFoundException
|
153 |
{ |
154 |
final Collection descriptors = PARAMETERS.descriptors(); |
155 |
if (isSpherical(parameters)) {
|
156 |
return new Spherical(parameters, descriptors); |
157 |
} else {
|
158 |
return new Mercator (parameters, descriptors); |
159 |
} |
160 |
} |
161 |
} |
162 |
|
163 |
/**
|
164 |
* The {@link MathTransformProvider} for a {@link Mercator} 2SP projection.
|
165 |
*
|
166 |
* @see <A HREF="http://www.remotesensing.org/geotiff/proj_list/mercator_2sp.html">"mercator_2sp" on Remote Sensing</A>
|
167 |
* @see org.geotools.referencing.operation.DefaultMathTransformFactory
|
168 |
*
|
169 |
* @version $Id: Mercator.java 12260 2007-06-21 09:48:20Z jlgomez $
|
170 |
* @author Martin Desruisseaux
|
171 |
* @author Rueben Schulz
|
172 |
*/
|
173 |
public static final class Provider2SP extends AbstractProvider { |
174 |
/**
|
175 |
* The operation parameter descriptor for the {@link #standardParallel standard parallel}
|
176 |
* parameter value. Valid values range is from -90 to 90. Default value is 0.
|
177 |
*/
|
178 |
public static final ParameterDescriptor STANDARD_PARALLEL = createDescriptor( |
179 |
new NamedIdentifier[] { |
180 |
new NamedIdentifier(CitationImpl.OGC, "standard_parallel_1"), |
181 |
new NamedIdentifier(CitationImpl.EPSG, "Latitude of 1st standard parallel"), |
182 |
new NamedIdentifier(CitationImpl.GEOTIFF, "StdParallel1") |
183 |
}, |
184 |
0, -90, 90, NonSI.DEGREE_ANGLE); |
185 |
|
186 |
/**
|
187 |
* The parameters group.
|
188 |
*/
|
189 |
static final ParameterDescriptorGroup PARAMETERS = createDescriptorGroup(new NamedIdentifier[] { |
190 |
new NamedIdentifier(CitationImpl.OGC, "Mercator_2SP"), |
191 |
new NamedIdentifier(CitationImpl.EPSG, "Mercator (2SP)"), |
192 |
new NamedIdentifier(CitationImpl.EPSG, "9805"), |
193 |
new NamedIdentifier(CitationImpl.GEOTIFF, "CT_Mercator"), |
194 |
new NamedIdentifier(CitationImpl.ESRI, "Mercator"), |
195 |
new NamedIdentifier(CitationImpl.GEOTOOLS, Resources.formatInternational(
|
196 |
ResourceKeys.CYLINDRICAL_MERCATOR_PROJECTION)) |
197 |
}, new ParameterDescriptor[] { |
198 |
SEMI_MAJOR, SEMI_MINOR, |
199 |
CENTRAL_MERIDIAN, STANDARD_PARALLEL, |
200 |
FALSE_EASTING, FALSE_NORTHING |
201 |
}); |
202 |
|
203 |
/**
|
204 |
* Constructs a new provider.
|
205 |
*/
|
206 |
public Provider2SP() {
|
207 |
super(PARAMETERS);
|
208 |
} |
209 |
|
210 |
/**
|
211 |
* Returns the operation type for this map projection.
|
212 |
*/
|
213 |
protected Class getOperationType() { |
214 |
return CylindricalProjection.class;
|
215 |
} |
216 |
|
217 |
/**
|
218 |
* Creates a transform from the specified group of parameter values.
|
219 |
*
|
220 |
* @param parameters The group of parameter values.
|
221 |
* @return The created math transform.
|
222 |
* @throws ParameterNotFoundException if a required parameter was not found.
|
223 |
*/
|
224 |
public MathTransform createMathTransform(final ParameterValueGroup parameters) |
225 |
throws ParameterNotFoundException
|
226 |
{ |
227 |
final Collection descriptors = PARAMETERS.descriptors(); |
228 |
if (isSpherical(parameters)) {
|
229 |
return new Spherical(parameters, descriptors); |
230 |
} else {
|
231 |
return new Mercator (parameters, descriptors); |
232 |
} |
233 |
} |
234 |
} |
235 |
|
236 |
|
237 |
/**
|
238 |
* Constructs a new map projection from the supplied parameters.
|
239 |
*
|
240 |
* @param parameters The parameter values in standard units.
|
241 |
* @throws ParameterNotFoundException if a mandatory parameter is missing.
|
242 |
*/
|
243 |
protected Mercator(final ParameterValueGroup parameters) |
244 |
throws ParameterNotFoundException
|
245 |
{ |
246 |
this(parameters, getDescriptor(parameters).descriptors());
|
247 |
} |
248 |
|
249 |
/**
|
250 |
* Work around for RFE #4093999 in Sun's bug database
|
251 |
* ("Relax constraint on placement of this()/super() call in constructors").
|
252 |
*/
|
253 |
private static ParameterDescriptorGroup getDescriptor(final ParameterValueGroup parameters) { |
254 |
try {
|
255 |
parameters.parameter(Provider2SP.STANDARD_PARALLEL.getName().getCode()); |
256 |
return Provider2SP.PARAMETERS;
|
257 |
} catch (ParameterNotFoundException ignore) {
|
258 |
return Provider1SP.PARAMETERS;
|
259 |
} |
260 |
} |
261 |
|
262 |
/**
|
263 |
* Constructs a new map projection from the supplied parameters.
|
264 |
*
|
265 |
* @param parameters The parameter values in standard units.
|
266 |
* @param expected The expected parameter descriptors.
|
267 |
* @throws ParameterNotFoundException if a mandatory parameter is missing.
|
268 |
*/
|
269 |
Mercator(final ParameterValueGroup parameters, final Collection expected) |
270 |
throws ParameterNotFoundException
|
271 |
{ |
272 |
//Fetch parameters
|
273 |
super(parameters, expected);
|
274 |
if (expected.contains(Provider2SP.STANDARD_PARALLEL)) {
|
275 |
// scaleFactor is not a parameter in the Mercator_2SP case and is computed from
|
276 |
// the standard parallel. The super-class constructor should have initialized
|
277 |
// 'scaleFactor' to 1. We still use the '*=' operator rather than '=' in case a
|
278 |
// user implementation still provides a scale factor for its custom projections.
|
279 |
standardParallel = Math.abs(doubleValue(expected,
|
280 |
Provider2SP.STANDARD_PARALLEL, parameters)); |
281 |
ensureLatitudeInRange(Provider2SP.STANDARD_PARALLEL, standardParallel, false);
|
282 |
if (isSpherical) {
|
283 |
scaleFactor *= Math.cos(standardParallel);
|
284 |
} else {
|
285 |
scaleFactor *= msfn(Math.sin(standardParallel),
|
286 |
Math.cos(standardParallel));
|
287 |
} |
288 |
globalScale = scaleFactor*semiMajor; |
289 |
} else {
|
290 |
// No standard parallel. Instead, uses the scale factor explicitely provided.
|
291 |
standardParallel = Double.NaN;
|
292 |
} |
293 |
//assert latitudeOfOrigin == 0 : latitudeOfOrigin;
|
294 |
} |
295 |
|
296 |
/**
|
297 |
* {@inheritDoc}
|
298 |
*/
|
299 |
public ParameterDescriptorGroup getParameterDescriptors() {
|
300 |
return Double.isNaN(standardParallel) ? Provider1SP.PARAMETERS |
301 |
: Provider2SP.PARAMETERS; |
302 |
} |
303 |
|
304 |
/**
|
305 |
* {@inheritDoc}
|
306 |
*/
|
307 |
public ParameterValueGroup getParameterValues() {
|
308 |
final ParameterValueGroup values = super.getParameterValues(); |
309 |
if (!Double.isNaN(standardParallel)) { |
310 |
final Collection expected = getParameterDescriptors().descriptors(); |
311 |
set(expected, Provider2SP.STANDARD_PARALLEL, values, standardParallel); |
312 |
} |
313 |
return values;
|
314 |
} |
315 |
|
316 |
/**
|
317 |
* Transforms the specified (<var>x</var>,<var>y</var>) coordinate (units in radians)
|
318 |
* and stores the result in <code>ptDst</code> (linear distance on a unit sphere).
|
319 |
*/
|
320 |
protected Point2D transformNormalized(double x, double y, final Point2D ptDst) |
321 |
throws ProjectionException
|
322 |
{ |
323 |
if (Math.abs(y) > (Math.PI/2 - EPS)) { |
324 |
throw new ProjectionException(Resources.format( |
325 |
ResourceKeys.ERROR_POLE_PROJECTION_$1, new Latitude(Math.toDegrees(y)))); |
326 |
} |
327 |
|
328 |
y = - Math.log(tsfn(y, Math.sin(y))); |
329 |
|
330 |
if (ptDst != null) { |
331 |
ptDst.setLocation(x,y); |
332 |
return ptDst;
|
333 |
} |
334 |
return new Point2D.Double(x,y); |
335 |
} |
336 |
|
337 |
/**
|
338 |
* Transforms the specified (<var>x</var>,<var>y</var>) coordinate
|
339 |
* and stores the result in <code>ptDst</code>.
|
340 |
*/
|
341 |
protected Point2D inverseTransformNormalized(double x, double y, final Point2D ptDst) |
342 |
throws ProjectionException
|
343 |
{ |
344 |
y = Math.exp(-y);
|
345 |
y = cphi2(y); |
346 |
|
347 |
if (ptDst != null) { |
348 |
ptDst.setLocation(x,y); |
349 |
return ptDst;
|
350 |
} |
351 |
return new Point2D.Double(x,y); |
352 |
} |
353 |
|
354 |
|
355 |
/**
|
356 |
* Provides the transform equations for the spherical case of the Mercator projection.
|
357 |
*
|
358 |
* @version $Id: Mercator.java 12260 2007-06-21 09:48:20Z jlgomez $
|
359 |
* @author Martin Desruisseaux
|
360 |
* @author Rueben Schulz
|
361 |
*/
|
362 |
private static final class Spherical extends Mercator { |
363 |
/**
|
364 |
* Constructs a new map projection from the suplied parameters.
|
365 |
*
|
366 |
* @param parameters The parameter values in standard units.
|
367 |
* @param expected The expected parameter descriptors.
|
368 |
* @throws ParameterNotFoundException if a mandatory parameter is missing.
|
369 |
*/
|
370 |
protected Spherical(final ParameterValueGroup parameters, final Collection expected) |
371 |
throws ParameterNotFoundException
|
372 |
{ |
373 |
super(parameters, expected);
|
374 |
//assert isSpherical;
|
375 |
} |
376 |
|
377 |
/**
|
378 |
* Transforms the specified (<var>x</var>,<var>y</var>) coordinate
|
379 |
* and stores the result in <code>ptDst</code> using equations for a Sphere.
|
380 |
*/
|
381 |
protected Point2D transformNormalized(double x, double y, Point2D ptDst) |
382 |
throws ProjectionException
|
383 |
{ |
384 |
if (Math.abs(y) > (Math.PI/2 - EPS)) { |
385 |
throw new ProjectionException(Resources.format( |
386 |
ResourceKeys.ERROR_POLE_PROJECTION_$1, new Latitude(Math.toDegrees(y)))); |
387 |
} |
388 |
// Compute using ellipsoidal formulas, for comparaison later.
|
389 |
//assert (ptDst = super.transformNormalized(x, y, ptDst)) != null;
|
390 |
|
391 |
y = Math.log(Math.tan((Math.PI/4) + 0.5*y)); |
392 |
|
393 |
//assert Math.abs(ptDst.getX()-x) <= EPS*globalScale : x;
|
394 |
//assert Math.abs(ptDst.getY()-y) <= EPS*globalScale : y;
|
395 |
if (ptDst != null) { |
396 |
ptDst.setLocation(x,y); |
397 |
return ptDst;
|
398 |
} |
399 |
return new Point2D.Double(x,y); |
400 |
} |
401 |
|
402 |
/**
|
403 |
* Transforms the specified (<var>x</var>,<var>y</var>) coordinate
|
404 |
* and stores the result in <code>ptDst</code> using equations for a sphere.
|
405 |
*/
|
406 |
protected Point2D inverseTransformNormalized(double x, double y, Point2D ptDst) |
407 |
throws ProjectionException
|
408 |
{ |
409 |
// Compute using ellipsoidal formulas, for comparaison later.
|
410 |
//assert (ptDst = super.inverseTransformNormalized(x, y, ptDst)) != null;
|
411 |
|
412 |
y = (Math.PI/2) - 2.0*Math.atan(Math.exp(-y)); |
413 |
|
414 |
//assert Math.abs(ptDst.getX()-x) <= EPS : x;
|
415 |
//assert Math.abs(ptDst.getY()-y) <= EPS : y;
|
416 |
if (ptDst != null) { |
417 |
ptDst.setLocation(x,y); |
418 |
return ptDst;
|
419 |
} |
420 |
return new Point2D.Double(x,y); |
421 |
} |
422 |
} |
423 |
|
424 |
|
425 |
/**
|
426 |
* Returns a hash value for this projection.
|
427 |
*/
|
428 |
public int hashCode() { |
429 |
final long code = Double.doubleToLongBits(standardParallel); |
430 |
return ((int)code ^ (int)(code >>> 32)) + 37*super.hashCode(); |
431 |
} |
432 |
|
433 |
/**
|
434 |
* Compares the specified object with this map projection for equality.
|
435 |
*/
|
436 |
public boolean equals(final Object object) { |
437 |
if (object == this) { |
438 |
// Slight optimization
|
439 |
return true; |
440 |
} |
441 |
if (super.equals(object)) { |
442 |
final Mercator that = (Mercator) object;
|
443 |
return equals(this.standardParallel, that.standardParallel); |
444 |
} |
445 |
return false; |
446 |
} |
447 |
} |