Source Code Cross Referenced for TransverseMercator.java in  » GIS » GeoTools-2.4.1 » org » geotools » referencing » operation » projection » Java Source Code / Java DocumentationJava Source Code and Java Documentation

Java Source Code / Java Documentation
1. 6.0 JDK Core
2. 6.0 JDK Modules
3. 6.0 JDK Modules com.sun
4. 6.0 JDK Modules com.sun.java
5. 6.0 JDK Modules sun
6. 6.0 JDK Platform
7. Ajax
8. Apache Harmony Java SE
9. Aspect oriented
10. Authentication Authorization
11. Blogger System
12. Build
13. Byte Code
14. Cache
15. Chart
16. Chat
17. Code Analyzer
18. Collaboration
19. Content Management System
20. Database Client
21. Database DBMS
22. Database JDBC Connection Pool
23. Database ORM
24. Development
25. EJB Server geronimo
26. EJB Server GlassFish
27. EJB Server JBoss 4.2.1
28. EJB Server resin 3.1.5
29. ERP CRM Financial
30. ESB
31. Forum
32. GIS
33. Graphic Library
34. Groupware
35. HTML Parser
36. IDE
37. IDE Eclipse
38. IDE Netbeans
39. Installer
40. Internationalization Localization
41. Inversion of Control
42. Issue Tracking
43. J2EE
44. JBoss
45. JMS
46. JMX
47. Library
48. Mail Clients
49. Net
50. Parser
51. PDF
52. Portal
53. Profiler
54. Project Management
55. Report
56. RSS RDF
57. Rule Engine
58. Science
59. Scripting
60. Search Engine
61. Security
62. Sevlet Container
63. Source Control
64. Swing Library
65. Template Engine
66. Test Coverage
67. Testing
68. UML
69. Web Crawler
70. Web Framework
71. Web Mail
72. Web Server
73. Web Services
74. Web Services apache cxf 2.0.1
75. Web Services AXIS2
76. Wiki Engine
77. Workflow Engines
78. XML
79. XML UI
Java
Java Tutorial
Java Open Source
Jar File Download
Java Articles
Java Products
Java by API
Photoshop Tutorials
Maya Tutorials
Flash Tutorials
3ds-Max Tutorials
Illustrator Tutorials
GIMP Tutorials
C# / C Sharp
C# / CSharp Tutorial
C# / CSharp Open Source
ASP.Net
ASP.NET Tutorial
JavaScript DHTML
JavaScript Tutorial
JavaScript Reference
HTML / CSS
HTML CSS Reference
C / ANSI-C
C Tutorial
C++
C++ Tutorial
Ruby
PHP
Python
Python Tutorial
Python Open Source
SQL Server / T-SQL
SQL Server / T-SQL Tutorial
Oracle PL / SQL
Oracle PL/SQL Tutorial
PostgreSQL
SQL / MySQL
MySQL Tutorial
VB.Net
VB.Net Tutorial
Flash / Flex / ActionScript
VBA / Excel / Access / Word
XML
XML Tutorial
Microsoft Office PowerPoint 2007 Tutorial
Microsoft Office Excel 2007 Tutorial
Microsoft Office Word 2007 Tutorial
Java Source Code / Java Documentation » GIS » GeoTools 2.4.1 » org.geotools.referencing.operation.projection 
Source Cross Referenced  Class Diagram Java Document (Java Doc) 


001:        /*
002:         *    GeoTools - OpenSource mapping toolkit
003:         *    http://geotools.org
004:         *
005:         *   (C) 2003-2006, Geotools Project Managment Committee (PMC)
006:         *   (C) 2002, Centre for Computational Geography
007:         *   (C) 2001, Institut de Recherche pour le Développement
008:         *   (C) 2000, Frank Warmerdam
009:         *   (C) 1999, Fisheries and Oceans Canada
010:         *
011:         *    This library is free software; you can redistribute it and/or
012:         *    modify it under the terms of the GNU Lesser General Public
013:         *    License as published by the Free Software Foundation; either
014:         *    version 2.1 of the License, or (at your option) any later version.
015:         *
016:         *    This library is distributed in the hope that it will be useful,
017:         *    but WITHOUT ANY WARRANTY; without even the implied warranty of
018:         *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
019:         *    Lesser General Public License for more details.
020:         *
021:         *    This package contains formulas from the PROJ package of USGS.
022:         *    USGS's work is fully acknowledged here. This derived work has
023:         *    been relicensed under LGPL with Frank Warmerdam's permission.
024:         */
025:        package org.geotools.referencing.operation.projection;
026:
027:        // J2SE dependencies and extensions
028:        import java.awt.geom.AffineTransform;
029:        import java.awt.geom.Point2D;
030:        import java.util.Collection;
031:
032:        // OpenGIS dependencies
033:        import org.opengis.parameter.ParameterDescriptor;
034:        import org.opengis.parameter.ParameterDescriptorGroup;
035:        import org.opengis.parameter.ParameterNotFoundException;
036:        import org.opengis.parameter.ParameterValueGroup;
037:        import org.opengis.referencing.operation.CylindricalProjection;
038:        import org.opengis.referencing.operation.MathTransform;
039:        import org.opengis.referencing.ReferenceIdentifier;
040:
041:        // Geotools dependencies
042:        import org.geotools.metadata.iso.citation.Citations;
043:        import org.geotools.referencing.NamedIdentifier;
044:        import org.geotools.referencing.operation.transform.ConcatenatedTransform;
045:        import org.geotools.referencing.operation.transform.ProjectiveTransform;
046:        import org.geotools.resources.i18n.VocabularyKeys;
047:        import org.geotools.resources.i18n.Vocabulary;
048:        import org.geotools.resources.i18n.ErrorKeys;
049:        import org.geotools.resources.i18n.Errors;
050:
051:        /**
052:         * Transverse Mercator Projection (EPSG code 9807). This
053:         * is a cylindrical projection, in which the cylinder has been rotated 90°.
054:         * Instead of being tangent to the equator (or to an other standard latitude),
055:         * it is tangent to a central meridian. Deformation are more important as we
056:         * are going futher from the central meridian. The Transverse Mercator
057:         * projection is appropriate for region wich have a greater extent north-south
058:         * than east-west.
059:         * <p>
060:         *
061:         * The elliptical equations used here are series approximations, and their accuracy
062:         * decreases as points move farther from the central meridian of the projection.
063:         * The forward equations here are accurate to a less than a mm &plusmn;10 degrees from 
064:         * the central meridian, a few mm &plusmn;15 degrees from the 
065:         * central meridian and a few cm &plusmn;20 degrees from the central meridian.
066:         * The spherical equations are not approximations and should always give the 
067:         * correct values.
068:         * <p>
069:         *
070:         * There are a number of versions of the transverse mercator projection 
071:         * including the Universal (UTM) and Modified (MTM) Transverses Mercator 
072:         * projections. In these cases the earth is divided into zones. For the UTM
073:         * the zones are 6 degrees wide, numbered from 1 to 60 proceeding east from 
074:         * 180 degrees longitude, and between lats 84 degrees North and 80 
075:         * degrees South. The central meridian is taken as the center of the zone
076:         * and the latitude of origin is the equator. A scale factor of 0.9996 and 
077:         * false easting of 500000m is used for all zones and a false northing of 10000000m
078:         * is used for zones in the southern hemisphere.
079:         * <p>
080:         *
081:         * NOTE: formulas used below are not those of Snyder, but rather those
082:         *       from the {@code proj4} package of the USGS survey, which
083:         *       have been reproduced verbatim. USGS work is acknowledged here.
084:         * <p>
085:         *
086:         * <strong>References:</strong><ul>
087:         *   <li> Proj-4.4.6 available at <A HREF="http://www.remotesensing.org/proj">www.remotesensing.org/proj</A><br>
088:         *        Relevent files are: {@code PJ_tmerc.c}, {@code pj_mlfn.c}, {@code pj_fwd.c} and {@code pj_inv.c}.</li>
089:         *   <li> John P. Snyder (Map Projections - A Working Manual,
090:         *        U.S. Geological Survey Professional Paper 1395, 1987).</li>
091:         *   <li> "Coordinate Conversions and Transformations including Formulas",
092:         *        EPSG Guidence Note Number 7, Version 19.</li>
093:         * </ul>
094:         *
095:         * @see <A HREF="http://mathworld.wolfram.com/MercatorProjection.html">Transverse Mercator projection on MathWorld</A>
096:         * @see <A HREF="http://www.remotesensing.org/geotiff/proj_list/transverse_mercator.html">"Transverse_Mercator" on RemoteSensing.org</A>
097:         *
098:         * @since 2.1
099:         * @source $URL: http://svn.geotools.org/geotools/tags/2.4.1/modules/library/referencing/src/main/java/org/geotools/referencing/operation/projection/TransverseMercator.java $
100:         * @version $Id: TransverseMercator.java 25262 2007-04-23 21:11:16Z desruisseaux $
101:         * @author André Gosselin
102:         * @author Martin Desruisseaux
103:         * @author Rueben Schulz
104:         */
105:        public class TransverseMercator extends MapProjection {
106:            /**
107:             * Maximum number of iterations for iterative computations.
108:             */
109:            private static final int MAXIMUM_ITERATIONS = 15;
110:
111:            /**
112:             * Relative iteration precision used in the {@code mlfn} method.
113:             * This overrides the value in the {@link MapProjection} class.
114:             */
115:            private static final double ITERATION_TOLERANCE = 1E-11;
116:
117:            /**
118:             * Maximum difference allowed when comparing real numbers.
119:             */
120:            private static final double EPSILON = 1E-6;
121:
122:            /**
123:             * Maximum difference allowed when comparing latitudes.
124:             */
125:            private static final double EPSILON_LATITUDE = 1E-10;
126:
127:            /**
128:             * A derived quantity of excentricity, computed by <code>e'˛ = (a˛-b˛)/b˛ = es/(1-es)</code>
129:             * where <var>a</var> is the semi-major axis length and <var>b</bar> is the semi-minor axis
130:             * length.
131:             */
132:            private final double esp;
133:
134:            /**
135:             * Meridian distance at the {@code latitudeOfOrigin}.
136:             * Used for calculations for the ellipsoid.
137:             */
138:            private final double ml0;
139:
140:            /**
141:             * Constant needed for the <code>mlfn<code> method.
142:             * Setup at construction time.
143:             */
144:            private final double en0, en1, en2, en3, en4;
145:
146:            /**
147:             * Constants used to calculate {@link #en0}, {@link #en1},
148:             * {@link #en2}, {@link #en3}, {@link #en4}.
149:             */
150:            private static final double C00 = 1.0, C02 = 0.25, C04 = 0.046875,
151:                    C06 = 0.01953125, C08 = 0.01068115234375, C22 = 0.75,
152:                    C44 = 0.46875, C46 = 0.01302083333333333333,
153:                    C48 = 0.00712076822916666666, C66 = 0.36458333333333333333,
154:                    C68 = 0.00569661458333333333, C88 = 0.3076171875;
155:
156:            /**
157:             * Contants used for the forward and inverse transform for the eliptical
158:             * case of the Transverse Mercator.
159:             */
160:            private static final double FC1 = 1.00000000000000000000000, // 1/1
161:                    FC2 = 0.50000000000000000000000, // 1/2
162:                    FC3 = 0.16666666666666666666666, // 1/6
163:                    FC4 = 0.08333333333333333333333, // 1/12
164:                    FC5 = 0.05000000000000000000000, // 1/20
165:                    FC6 = 0.03333333333333333333333, // 1/30
166:                    FC7 = 0.02380952380952380952380, // 1/42
167:                    FC8 = 0.01785714285714285714285; // 1/56
168:
169:            /**
170:             * Constructs a new map projection from the supplied parameters.
171:             *
172:             * @param  parameters The parameter values in standard units.
173:             * @throws ParameterNotFoundException if a mandatory parameter is missing.
174:             */
175:            protected TransverseMercator(final ParameterValueGroup parameters)
176:                    throws ParameterNotFoundException {
177:                //Fetch parameters 
178:                super (parameters);
179:
180:                //  Compute constants
181:                esp = excentricitySquared / (1.0 - excentricitySquared);
182:
183:                double t;
184:                en0 = C00
185:                        - excentricitySquared
186:                        * (C02 + excentricitySquared
187:                                * (C04 + excentricitySquared
188:                                        * (C06 + excentricitySquared * C08)));
189:                en1 = excentricitySquared
190:                        * (C22 - excentricitySquared
191:                                * (C04 + excentricitySquared
192:                                        * (C06 + excentricitySquared * C08)));
193:                en2 = (t = excentricitySquared * excentricitySquared)
194:                        * (C44 - excentricitySquared
195:                                * (C46 + excentricitySquared * C48));
196:                en3 = (t *= excentricitySquared)
197:                        * (C66 - excentricitySquared * C68);
198:                en4 = t * excentricitySquared * C88;
199:                ml0 = mlfn(latitudeOfOrigin, Math.sin(latitudeOfOrigin), Math
200:                        .cos(latitudeOfOrigin));
201:            }
202:
203:            /**
204:             * {@inheritDoc}
205:             */
206:            public ParameterDescriptorGroup getParameterDescriptors() {
207:                return Provider.PARAMETERS;
208:            }
209:
210:            /**
211:             * Transforms the specified (<var>&lambda;</var>,<var>&phi;</var>) coordinates
212:             * (units in radians) and stores the result in {@code ptDst} (linear distance
213:             * on a unit sphere).
214:             */
215:            protected Point2D transformNormalized(double x, double y,
216:                    Point2D ptDst) throws ProjectionException {
217:                double sinphi = Math.sin(y);
218:                double cosphi = Math.cos(y);
219:
220:                double t = (Math.abs(cosphi) > EPSILON) ? sinphi / cosphi : 0;
221:                t *= t;
222:                double al = cosphi * x;
223:                double als = al * al;
224:                al /= Math.sqrt(1.0 - excentricitySquared * sinphi * sinphi);
225:                double n = esp * cosphi * cosphi;
226:
227:                /* NOTE: meridinal distance at latitudeOfOrigin is always 0 */
228:                y = (mlfn(y, sinphi, cosphi) - ml0 + sinphi
229:                        * al
230:                        * x
231:                        * FC2
232:                        * (1.0 + FC4
233:                                * als
234:                                * (5.0 - t + n * (9.0 + 4.0 * n) + FC6
235:                                        * als
236:                                        * (61.0 + t * (t - 58.0) + n
237:                                                * (270.0 - 330.0 * t) + FC8
238:                                                * als
239:                                                * (1385.0 + t
240:                                                        * (t * (543.0 - t) - 3111.0))))));
241:
242:                x = al
243:                        * (FC1 + FC3
244:                                * als
245:                                * (1.0 - t + n + FC5
246:                                        * als
247:                                        * (5.0 + t * (t - 18.0) + n
248:                                                * (14.0 - 58.0 * t) + FC7
249:                                                * als
250:                                                * (61.0 + t
251:                                                        * (t * (179.0 - t) - 479.0)))));
252:
253:                if (ptDst != null) {
254:                    ptDst.setLocation(x, y);
255:                    return ptDst;
256:                }
257:                return new Point2D.Double(x, y);
258:            }
259:
260:            /**
261:             * Transforms the specified (<var>x</var>,<var>y</var>) coordinates
262:             * and stores the result in {@code ptDst}.
263:             */
264:            protected Point2D inverseTransformNormalized(double x, double y,
265:                    Point2D ptDst) throws ProjectionException {
266:                double phi = inv_mlfn(ml0 + y);
267:
268:                if (Math.abs(phi) >= (Math.PI / 2)) {
269:                    y = y < 0.0 ? -(Math.PI / 2) : (Math.PI / 2);
270:                    x = 0.0;
271:                } else {
272:                    double sinphi = Math.sin(phi);
273:                    double cosphi = Math.cos(phi);
274:                    double t = (Math.abs(cosphi) > EPSILON) ? sinphi / cosphi
275:                            : 0.0;
276:                    double n = esp * cosphi * cosphi;
277:                    double con = 1.0 - excentricitySquared * sinphi * sinphi;
278:                    double d = x * Math.sqrt(con);
279:                    con *= t;
280:                    t *= t;
281:                    double ds = d * d;
282:
283:                    y = phi
284:                            - (con * ds / (1.0 - excentricitySquared))
285:                            * FC2
286:                            * (1.0 - ds
287:                                    * FC4
288:                                    * (5.0 + t * (3.0 - 9.0 * n) + n
289:                                            * (1.0 - 4 * n) - ds
290:                                            * FC6
291:                                            * (61.0
292:                                                    + t
293:                                                    * (90.0 - 252.0 * n + 45.0 * t)
294:                                                    + 46.0 * n - ds
295:                                                    * FC8
296:                                                    * (1385.0 + t
297:                                                            * (3633.0 + t
298:                                                                    * (4095.0 + 1574.0 * t))))));
299:
300:                    x = d
301:                            * (FC1 - ds
302:                                    * FC3
303:                                    * (1.0 + 2.0 * t + n - ds
304:                                            * FC5
305:                                            * (5.0 + t
306:                                                    * (28.0 + 24 * t + 8.0 * n)
307:                                                    + 6.0 * n - ds
308:                                                    * FC7
309:                                                    * (61.0 + t
310:                                                            * (662.0 + t
311:                                                                    * (1320.0 + 720.0 * t))))))
312:                            / cosphi;
313:                }
314:
315:                if (ptDst != null) {
316:                    ptDst.setLocation(x, y);
317:                    return ptDst;
318:                }
319:                return new Point2D.Double(x, y);
320:            }
321:
322:            /**
323:             * {@inheritDoc}
324:             */
325:            protected double getToleranceForAssertions(final double longitude,
326:                    final double latitude) {
327:                if (Math.abs(longitude - centralMeridian) > 0.26) { // 15 degrees
328:                    // When far from the valid area, use a larger tolerance.
329:                    return 2.5;
330:                } else if (Math.abs(longitude - centralMeridian) > 0.22) { // 12.5 degrees
331:                    return 1.0;
332:                } else if (Math.abs(longitude - centralMeridian) > 0.17) { // 10 degrees
333:                    return 0.5;
334:                }
335:                // a normal tolerance
336:                return 1E-6;
337:            }
338:
339:            /**
340:             * Provides the transform equations for the spherical case of the
341:             * TransverseMercator projection.
342:             * 
343:             * @version $Id: TransverseMercator.java 25262 2007-04-23 21:11:16Z desruisseaux $
344:             * @author André Gosselin
345:             * @author Martin Desruisseaux
346:             * @author Rueben Schulz
347:             */
348:            private static final class Spherical extends TransverseMercator {
349:                /**
350:                 * Constructs a new map projection from the suplied parameters.
351:                 *
352:                 * @param  parameters The parameter values in standard units.
353:                 * @throws ParameterNotFoundException if a mandatory parameter is missing.
354:                 */
355:                protected Spherical(final ParameterValueGroup parameters)
356:                        throws ParameterNotFoundException {
357:                    super (parameters);
358:                    ensureSpherical();
359:                }
360:
361:                /**
362:                 * {@inheritDoc}
363:                 */
364:                protected Point2D transformNormalized(double x, double y,
365:                        Point2D ptDst) throws ProjectionException {
366:                    // Compute using ellipsoidal formulas, for comparaison later.
367:                    final double normalizedLongitude = x;
368:                    assert (ptDst = super .transformNormalized(x, y, ptDst)) != null;
369:
370:                    double cosphi = Math.cos(y);
371:                    double b = cosphi * Math.sin(x);
372:                    if (Math.abs(Math.abs(b) - 1.0) <= EPSILON) {
373:                        throw new ProjectionException(Errors
374:                                .format(ErrorKeys.VALUE_TEND_TOWARD_INFINITY));
375:                    }
376:
377:                    //Using Snyder's equation for calculating y, instead of the one used in Proj4
378:                    //poential problems when y and x = 90 degrees, but behaves ok in tests   
379:                    y = Math.atan2(Math.tan(y), Math.cos(x)) - latitudeOfOrigin; /* Snyder 8-3 */
380:                    x = 0.5 * Math.log((1.0 + b) / (1.0 - b)); /* Snyder 8-1 */
381:
382:                    assert checkTransform(
383:                            x,
384:                            y,
385:                            ptDst,
386:                            getToleranceForSphereAssertions(normalizedLongitude));
387:                    if (ptDst != null) {
388:                        ptDst.setLocation(x, y);
389:                        return ptDst;
390:                    }
391:                    return new Point2D.Double(x, y);
392:                }
393:
394:                /**
395:                 * {@inheritDoc}
396:                 */
397:                protected Point2D inverseTransformNormalized(double x,
398:                        double y, Point2D ptDst) throws ProjectionException {
399:                    // Compute using ellipsoidal formulas, for comparaison later.
400:                    assert (ptDst = super .inverseTransformNormalized(x, y,
401:                            ptDst)) != null;
402:
403:                    double t = Math.exp(x);
404:                    double sinhX = 0.5 * (t - 1.0 / t); // sinh(x)
405:                    double cosD = Math.cos(latitudeOfOrigin + y);
406:                    double phi = Math.asin(Math.sqrt((1.0 - cosD * cosD)
407:                            / (1.0 + sinhX * sinhX)));
408:                    // correct for the fact that we made everything positive using sqrt(x*x)
409:                    y = ((y + latitudeOfOrigin) < 0.0) ? -phi : phi;
410:                    x = (Math.abs(sinhX) <= EPSILON && Math.abs(cosD) <= EPSILON) ? 0.0
411:                            : Math.atan2(sinhX, cosD);
412:
413:                    assert checkInverseTransform(x, y, ptDst,
414:                            getToleranceForSphereAssertions(x));
415:                    if (ptDst != null) {
416:                        ptDst.setLocation(x, y);
417:                        return ptDst;
418:                    }
419:                    return new Point2D.Double(x, y);
420:                }
421:
422:                /**
423:                 * Maximal error tolerated for assertions in the spherical case. When assertions
424:                 * are enabled, every projection using spherical formulas is followed by a projection
425:                 * using the ellipsical formulas, and the results are compared. If a distance greater
426:                 * than the tolerance level is found, then an {@link AssertionError} will be thrown.
427:                 * Subclasses may override this method if they need to relax the tolerance level.
428:                 * <p>
429:                 * The default implementation allows a larger tolerance when comparing spherical to
430:                 * elliptical calculations when we are far from the central meridian (elliptical
431:                 * calculations are not valid here).
432:                 *
433:                 * @param  longitude The longitude standardized (in radians) with
434:                 *         {@linkplain #centralMeridian} already removed from it.
435:                 * @return The tolerance level for assertions, in meters.
436:                 */
437:                protected double getToleranceForSphereAssertions(
438:                        final double longitude) {
439:                    if (Math.abs(Math.abs(longitude) - Math.PI / 2.0) < EPSILON_LATITUDE) { // 90 degrees
440:                        // elliptical equations are at their worst here
441:                        return 1E+18;
442:                    }
443:                    if (Math.abs(longitude) > 0.26) { // 15 degrees
444:                        // When far from the valid area, use a very larger tolerance.          
445:                        return 1E+6;
446:                    }
447:                    // a normal tolerance
448:                    return 1E-6;
449:                }
450:            }
451:
452:            /**
453:             * Calculates the meridian distance. This is the distance along the central 
454:             * meridian from the equator to {@code phi}. Accurate to < 1e-5 meters 
455:             * when used in conjuction with typical major axis values.
456:             *
457:             * @param phi latitude to calculate meridian distance for.
458:             * @param sphi sin(phi).
459:             * @param cphi cos(phi).
460:             * @return meridian distance for the given latitude.
461:             */
462:            private final double mlfn(final double phi, double sphi, double cphi) {
463:                cphi *= sphi;
464:                sphi *= sphi;
465:                return en0 * phi - cphi
466:                        * (en1 + sphi * (en2 + sphi * (en3 + sphi * (en4))));
467:            }
468:
469:            /**
470:             * Calculates the latitude ({@code phi}) from a meridian distance.
471:             * Determines phi to {@link #ITERATION_TOLERANCE} radians, about
472:             * 1e-6 seconds.
473:             * 
474:             * @param  arg meridian distance to calulate latitude for.
475:             * @return the latitude of the meridian distance.
476:             * @throws ProjectionException if the iteration does not converge.
477:             */
478:            private final double inv_mlfn(double arg)
479:                    throws ProjectionException {
480:                double s, t, phi, k = 1.0 / (1.0 - excentricitySquared);
481:                int i;
482:                phi = arg;
483:                for (i = MAXIMUM_ITERATIONS; true;) { // rarely goes over 5 iterations
484:                    if (--i < 0) {
485:                        throw new ProjectionException(Errors
486:                                .format(ErrorKeys.NO_CONVERGENCE));
487:                    }
488:                    s = Math.sin(phi);
489:                    t = 1.0 - excentricitySquared * s * s;
490:                    t = (mlfn(phi, s, Math.cos(phi)) - arg)
491:                            * (t * Math.sqrt(t)) * k;
492:                    phi -= t;
493:                    if (Math.abs(t) < ITERATION_TOLERANCE) {
494:                        return phi;
495:                    }
496:                }
497:            }
498:
499:            /**
500:             * Convenience method computing the zone code from the central meridian.
501:             * Information about zones convention must be specified in argument. Two
502:             * widely set of arguments are of Universal Transverse Mercator (UTM) and
503:             * Modified Transverse Mercator (MTM) projections:<br>
504:             * <p>
505:             *
506:             * UTM projection (zones numbered from 1 to 60):<br>
507:             * <p>
508:             *        {@code getZone(-177, 6);}<br>
509:             * <p>
510:             * MTM projection (zones numbered from 1 to 120):<br>
511:             * <p>
512:             *        {@code getZone(-52.5, -3);}<br>
513:             *
514:             * @param  centralLongitudeZone1 Longitude in the middle of zone 1, in decimal degrees
515:             *         relative to Greenwich. Positive longitudes are toward east, and negative
516:             *         longitudes toward west.
517:             * @param  zoneWidth Number of degrees of longitudes in one zone. A positive value
518:             *         means that zones are numbered from west to east (i.e. in the direction of
519:             *         positive longitudes). A negative value means that zones are numbered from
520:             *         east to west.
521:             * @return The zone number. First zone is numbered 1.
522:             */
523:            private int getZone(final double centralLongitudeZone1,
524:                    final double zoneWidth) {
525:                final double zoneCount = Math.abs(360 / zoneWidth);
526:                double t;
527:                t = centralLongitudeZone1 - 0.5 * zoneWidth; // Longitude at the beginning of the first zone.
528:                t = Math.toDegrees(centralMeridian) - t; // Degrees of longitude between the central longitude and longitude 1.
529:                t = Math.floor(t / zoneWidth + EPSILON); // Number of zones between the central longitude and longitude 1.
530:                t -= zoneCount * Math.floor(t / zoneCount); // If negative, bring back to the interval 0 to (zoneCount-1).
531:                return ((int) t) + 1;
532:            }
533:
534:            /**
535:             * Convenience method returning the meridian in the middle of current zone. This meridian is
536:             * typically the central meridian. This method may be invoked to make sure that the central
537:             * meridian is correctly set.
538:             *
539:             * @param  centralLongitudeZone1 Longitude in the middle of zone 1, in decimal degrees
540:             *         relative to Greenwich. Positive longitudes are toward east, and negative
541:             *         longitudes toward west.
542:             * @param  zoneWidth Number of degrees of longitudes in one zone. A positive value
543:             *         means that zones are numbered from west to east (i.e. in the direction of
544:             *         positive longitudes). A negative value means that zones are numbered from
545:             *         east to west.
546:             * @return The central meridian.
547:             */
548:            private double getCentralMedirian(
549:                    final double centralLongitudeZone1, final double zoneWidth) {
550:                double t;
551:                t = centralLongitudeZone1
552:                        + (getZone(centralLongitudeZone1, zoneWidth) - 1)
553:                        * zoneWidth;
554:                t -= 360 * Math.floor((t + 180) / 360); // Bring back into [-180..+180] range.
555:                return t;
556:            }
557:
558:            /**
559:             * Convenience method computing the zone code from the central meridian.
560:             *
561:             * @return The zone number, using the scalefactor and false easting 
562:             *         to decide if this is a UTM or MTM case. Returns 0 if the 
563:             *         case of the projection cannot be determined.
564:             */
565:            public int getZone() {
566:                // UTM
567:                if (scaleFactor == 0.9996 && falseEasting == 500000.0) {
568:                    return getZone(-177, 6);
569:                }
570:                // MTM
571:                if (scaleFactor == 0.9999 && falseEasting == 304800.0) {
572:                    return getZone(-52.5, -3);
573:                }
574:                // unknown
575:                throw new IllegalStateException(Errors
576:                        .format(ErrorKeys.UNKNOW_PROJECTION_TYPE));
577:            }
578:
579:            /**
580:             * Convenience method returning the meridian in the middle of current zone. This meridian is
581:             * typically the central meridian. This method may be invoked to make sure that the central
582:             * meridian is correctly set.
583:             *
584:             * @return The central meridian, using the scalefactor and false easting 
585:             *         to decide if this is a UTM or MTM case. Returns {@link Double#NaN}
586:             *         if the case of the projection cannot be determined.
587:             */
588:            public double getCentralMeridian() {
589:                // UTM
590:                if (scaleFactor == 0.9996 && falseEasting == 500000.0) {
591:                    return getCentralMedirian(-177, 6);
592:                }
593:                // MTM
594:                if (scaleFactor == 0.9999 && falseEasting == 304800.0) {
595:                    return getCentralMedirian(-52.5, -3);
596:                }
597:                // unknown
598:                throw new IllegalStateException(Errors
599:                        .format(ErrorKeys.UNKNOW_PROJECTION_TYPE));
600:            }
601:
602:            /**
603:             * Returns a hash value for this projection.
604:             */
605:            public int hashCode() {
606:                final long code = Double.doubleToLongBits(ml0);
607:                return ((int) code ^ (int) (code >>> 32)) + 37
608:                        * super .hashCode();
609:            }
610:
611:            /**
612:             * Compares the specified object with this map projection for equality.
613:             */
614:            public boolean equals(final Object object) {
615:                if (object == this ) {
616:                    // Slight optimization
617:                    return true;
618:                }
619:                // Relevant parameters are already compared in MapProjection
620:                return super .equals(object);
621:            }
622:
623:            //////////////////////////////////////////////////////////////////////////////////////////
624:            //////////////////////////////////////////////////////////////////////////////////////////
625:            ////////                                                                          ////////
626:            ////////                                 PROVIDERS                                ////////
627:            ////////                                                                          ////////
628:            //////////////////////////////////////////////////////////////////////////////////////////
629:            //////////////////////////////////////////////////////////////////////////////////////////
630:
631:            /**
632:             * The {@linkplain org.geotools.referencing.operation.MathTransformProvider math transform
633:             * provider} for a {@linkplain TransverseMercator Transverse Mercator} projection (EPSG code
634:             * 9807).
635:             *
636:             * @since 2.1
637:             * @version $Id: TransverseMercator.java 25262 2007-04-23 21:11:16Z desruisseaux $
638:             * @author Martin Desruisseaux
639:             * @author Rueben Schulz
640:             *
641:             * @see org.geotools.referencing.operation.DefaultMathTransformFactory
642:             */
643:            public static class Provider extends AbstractProvider {
644:                /**
645:                 * Returns a descriptor group for the specified parameters.
646:                 */
647:                static ParameterDescriptorGroup createDescriptorGroup(
648:                        final ReferenceIdentifier[] identifiers) {
649:                    return createDescriptorGroup(
650:                            identifiers,
651:                            new ParameterDescriptor[] { SEMI_MAJOR, SEMI_MINOR,
652:                                    CENTRAL_MERIDIAN, LATITUDE_OF_ORIGIN,
653:                                    SCALE_FACTOR, FALSE_EASTING, FALSE_NORTHING });
654:                }
655:
656:                /**
657:                 * The parameters group.
658:                 */
659:                static final ParameterDescriptorGroup PARAMETERS = createDescriptorGroup(new NamedIdentifier[] {
660:                        new NamedIdentifier(Citations.OGC,
661:                                "Transverse_Mercator"),
662:                        new NamedIdentifier(Citations.EPSG,
663:                                "Transverse Mercator"),
664:                        new NamedIdentifier(Citations.EPSG, "Gauss-Kruger"),
665:                        new NamedIdentifier(Citations.EPSG, "9807"),
666:                        new NamedIdentifier(Citations.GEOTIFF,
667:                                "CT_TransverseMercator"),
668:                        new NamedIdentifier(Citations.ESRI,
669:                                "Transverse_Mercator"),
670:                        new NamedIdentifier(Citations.ESRI, "Gauss_Kruger"),
671:                        new NamedIdentifier(
672:                                Citations.GEOTOOLS,
673:                                Vocabulary
674:                                        .formatInternational(VocabularyKeys.TRANSVERSE_MERCATOR_PROJECTION)) });
675:
676:                /**
677:                 * Constructs a new provider.
678:                 */
679:                public Provider() {
680:                    super (PARAMETERS);
681:                }
682:
683:                /**
684:                 * Constructs a new provider with the specified parameters.
685:                 */
686:                Provider(final ParameterDescriptorGroup descriptor) {
687:                    super (descriptor);
688:                }
689:
690:                /**
691:                 * Returns the operation type for this map projection.
692:                 */
693:                public Class getOperationType() {
694:                    return CylindricalProjection.class;
695:                }
696:
697:                /**
698:                 * Creates a transform from the specified group of parameter values.
699:                 *
700:                 * @param  parameters The group of parameter values.
701:                 * @return The created math transform.
702:                 * @throws ParameterNotFoundException if a required parameter was not found.
703:                 */
704:                protected MathTransform createMathTransform(
705:                        final ParameterValueGroup parameters)
706:                        throws ParameterNotFoundException {
707:                    if (isSpherical(parameters)) {
708:                        return new Spherical(parameters);
709:                    } else {
710:                        return new TransverseMercator(parameters);
711:                    }
712:                }
713:            }
714:
715:            /**
716:             * The {@linkplain org.geotools.referencing.operation.MathTransformProvider math transform
717:             * provider} for a South Orientated {@linkplain TransverseMercator Transverse Mercator}
718:             * projection (EPSG code 9808). Note that at the contrary of what this class name suggest,
719:             * the projected coordinates are still increasing toward North. This is because all
720:             * {@linkplain MapProjection map projections} in Geotools must complies with standard axis
721:             * orientations. The real axis flip is performed by the CRS framework outside this package.
722:             * See "<cite>Axis units and orientation</cite>" in
723:             * {@linkplain org.geotools.referencing.operation.projection package description} for details.
724:             * <p>
725:             * The usual Transverse Mercator formulas are:
726:             * <p>
727:             * <ul>
728:             *   <li><var>easting</var> = <var>{@linkplain #falseEasting false easting}</var> + <var>px</var></li>
729:             *   <li><var>northing</var> = <var>{@linkplain #falseNorthing false northing}</var> + <var>py</var></li>
730:             * </ul>
731:             * <p>
732:             * The Transverse Mercator South Orientated Projection formulas are:
733:             * <p>
734:             * <ul>
735:             *   <li><var>westing</var> = <var>{@linkplain #falseEasting false easting}</var> - <var>px</var></li>
736:             *   <li><var>southing</var> = <var>{@linkplain #falseNorthing false northing}</var> - <var>py</var></li>
737:             * </ul>
738:             * <p>
739:             * Where the <var>px</var> and <var>py</var> terms are the same in both cases.
740:             * Transforms created by this provider actually computes
741:             * (<var>easting</var>,<var>northing</var>) = (-<var>westing</var>,-<var>southing</var>).
742:             * This is equivalent to a {@link TransverseMercator} projection with
743:             * {@link #falseEasting falseEasting} and {@link #falseNorthing falseNorthing} sign reverted.
744:             * This operation is implemented as a concatenation of a North-orientated transverse mercator
745:             * projection with an affine transform for (<var>false easting</var>,<var>false northing</var>)
746:             * correction.
747:             *
748:             * @since 2.2
749:             * @version $Id: TransverseMercator.java 25262 2007-04-23 21:11:16Z desruisseaux $
750:             * @author Martin Desruisseaux
751:             *
752:             * @see org.geotools.referencing.operation.DefaultMathTransformFactory
753:             */
754:            public static class Provider_SouthOrientated extends Provider {
755:                /**
756:                 * Constructs a new provider.
757:                 */
758:                public Provider_SouthOrientated() {
759:                    super (createDescriptorGroup(new NamedIdentifier[] {
760:                            new NamedIdentifier(Citations.EPSG,
761:                                    "Transverse Mercator (South Orientated)"),
762:                            new NamedIdentifier(Citations.EPSG, "9808") }));
763:                }
764:
765:                /**
766:                 * Creates a transform from the specified group of parameter values.
767:                 *
768:                 * @param  parameters The group of parameter values.
769:                 * @return The created math transform.
770:                 * @throws ParameterNotFoundException if a required parameter was not found.
771:                 */
772:                public MathTransform createMathTransform(
773:                        final ParameterValueGroup parameters)
774:                        throws ParameterNotFoundException {
775:                    final MapProjection projection = (MapProjection) super 
776:                            .createMathTransform(parameters);
777:                    if (projection.falseEasting == 0
778:                            && projection.falseNorthing == 0) {
779:                        return projection;
780:                    }
781:                    final AffineTransform step = AffineTransform
782:                            .getTranslateInstance(-2 * projection.falseEasting,
783:                                    -2 * projection.falseNorthing);
784:                    return ConcatenatedTransform.create(ProjectiveTransform
785:                            .create(step), projection);
786:                }
787:            }
788:        }
www.java2java.com | Contact Us
Copyright 2009 - 12 Demo Source and Support. All rights reserved.
All other trademarks are property of their respective owners.