Statistics
| Revision:

root / branches / F2 / libraries / libJCRS / src / org / geotools / referencing / operation / projection / TransverseMercator.java @ 12186

History | View | Annotate | Download (29.8 KB)

1
/*
2
 *    GeoTools - OpenSource mapping toolkit
3
 *    http://geotools.org
4
 *
5
 *   (C) 2003-2006, Geotools Project Managment Committee (PMC)
6
 *   (C) 2002, Centre for Computational Geography
7
 *   (C) 2001, Institut de Recherche pour le Dveloppement
8
 *   (C) 2000, Frank Warmerdam
9
 *   (C) 1999, Fisheries and Oceans Canada
10
 *
11
 *    This library is free software; you can redistribute it and/or
12
 *    modify it under the terms of the GNU Lesser General Public
13
 *    License as published by the Free Software Foundation; either
14
 *    version 2.1 of the License, or (at your option) any later version.
15
 *
16
 *    This library is distributed in the hope that it will be useful,
17
 *    but WITHOUT ANY WARRANTY; without even the implied warranty of
18
 *    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
19
 *    Lesser General Public License for more details.
20
 *
21
 *    This package contains formulas from the PROJ package of USGS.
22
 *    USGS's work is fully acknowledged here. This derived work has
23
 *    been relicensed under LGPL with Frank Warmerdam's permission.
24
 */
25
package org.geotools.referencing.operation.projection;
26

    
27
// J2SE dependencies and extensions
28
import java.awt.geom.AffineTransform;
29
import java.awt.geom.Point2D;
30
import java.util.Collection;
31

    
32
// OpenGIS dependencies
33
import org.opengis.metadata.Identifier;
34
import org.opengis.parameter.ParameterDescriptor;
35
import org.opengis.parameter.ParameterDescriptorGroup;
36
import org.opengis.parameter.ParameterNotFoundException;
37
import org.opengis.parameter.ParameterValueGroup;
38
import org.opengis.referencing.operation.CylindricalProjection;
39
import org.opengis.referencing.operation.MathTransform;
40

    
41
// Geotools dependencies
42
import org.geotools.metadata.iso.citation.CitationImpl;
43
import org.geotools.referencing.NamedIdentifier;
44
import org.geotools.referencing.operation.transform.ConcatenatedTransform;
45
import org.geotools.referencing.operation.transform.ProjectiveTransform;
46
//import org.geotools.resources.i18n.VocabularyKeys;
47
//import org.geotools.resources.i18n.Vocabulary;
48
import org.geotools.resources.cts.ResourceKeys;
49
import org.geotools.resources.cts.Resources;
50

    
51

    
52
/**
53
 * Transverse Mercator Projection (EPSG code 9807). This
54
 * is a cylindrical projection, in which the cylinder has been rotated 90.
55
 * Instead of being tangent to the equator (or to an other standard latitude),
56
 * it is tangent to a central meridian. Deformation are more important as we
57
 * are going futher from the central meridian. The Transverse Mercator
58
 * projection is appropriate for region wich have a greater extent north-south
59
 * than east-west.
60
 * <p>
61
 *
62
 * The elliptical equations used here are series approximations, and their accuracy
63
 * decreases as points move farther from the central meridian of the projection.
64
 * The forward equations here are accurate to a less than a mm &plusmn;10 degrees from 
65
 * the central meridian, a few mm &plusmn;15 degrees from the 
66
 * central meridian and a few cm &plusmn;20 degrees from the central meridian.
67
 * The spherical equations are not approximations and should always give the 
68
 * correct values.
69
 * <p>
70
 *
71
 * There are a number of versions of the transverse mercator projection 
72
 * including the Universal (UTM) and Modified (MTM) Transverses Mercator 
73
 * projections. In these cases the earth is divided into zones. For the UTM
74
 * the zones are 6 degrees wide, numbered from 1 to 60 proceeding east from 
75
 * 180 degrees longitude, and between lats 84 degrees North and 80 
76
 * degrees South. The central meridian is taken as the center of the zone
77
 * and the latitude of origin is the equator. A scale factor of 0.9996 and 
78
 * false easting of 500000m is used for all zones and a false northing of 10000000m
79
 * is used for zones in the southern hemisphere.
80
 * <p>
81
 *
82
 * NOTE: formulas used below are not those of Snyder, but rather those
83
 *       from the {@code proj4} package of the USGS survey, which
84
 *       have been reproduced verbatim. USGS work is acknowledged here.
85
 * <p>
86
 *
87
 * <strong>References:</strong><ul>
88
 *   <li> Proj-4.4.6 available at <A HREF="http://www.remotesensing.org/proj">www.remotesensing.org/proj</A><br>
89
 *        Relevent files are: PJ_tmerc.c, pj_mlfn.c, pj_fwd.c and pj_inv.c </li>
90
 *   <li> John P. Snyder (Map Projections - A Working Manual,
91
 *        U.S. Geological Survey Professional Paper 1395, 1987)</li>
92
 *   <li> "Coordinate Conversions and Transformations including Formulas",
93
 *        EPSG Guidence Note Number 7, Version 19.</li>
94
 * </ul>
95
 *
96
 * @see <A HREF="http://mathworld.wolfram.com/MercatorProjection.html">Transverse Mercator projection on MathWorld</A>
97
 * @see <A HREF="http://www.remotesensing.org/geotiff/proj_list/transverse_mercator.html">"Transverse_Mercator" on Remote Sensing</A>
98
 *
99
 * @since 2.1
100
 * @source $URL: http://svn.geotools.org/geotools/tags/2.3.0/module/referencing/src/org/geotools/referencing/operation/projection/TransverseMercator.java $
101
 * @version $Id: TransverseMercator.java 12186 2007-06-18 08:46:09Z jlgomez $
102
 * @author Andr Gosselin
103
 * @author Martin Desruisseaux
104
 * @author Rueben Schulz
105
 */
106
public class TransverseMercator extends MapProjection {
107
    /**
108
     * A derived quantity of excentricity, computed by <code>e' = (a鲲-b)/b = es/(1-es)</code>
109
     * where <var>a</var> is the semi-major axis length and <var>b</bar> is the semi-minor axis
110
     * length.
111
     */
112
    private final double esp;
113
    
114
    /**
115
     * Meridian distance at the {@code latitudeOfOrigin}.
116
     * Used for calculations for the ellipsoid.
117
     */
118
    private final double ml0;
119

    
120
     /**
121
     * Constant needed for the <code>mlfn<code> method.
122
     * Setup at construction time.
123
     */
124
    private final double en0,en1,en2,en3,en4;
125
    
126
    /**
127
     * Constants used to calculate {@link #en0}, {@link #en1},
128
     * {@link #en2}, {@link #en3}, {@link #en4}.
129
     */
130
    private static final double C00= 1.0,
131
                                C02= 0.25,
132
                                C04= 0.046875,
133
                                C06= 0.01953125,
134
                                C08= 0.01068115234375,
135
                                C22= 0.75,
136
                                C44= 0.46875,
137
                                C46= 0.01302083333333333333,
138
                                C48= 0.00712076822916666666,
139
                                C66= 0.36458333333333333333,
140
                                C68= 0.00569661458333333333,
141
                                C88= 0.3076171875;
142

    
143
    /**
144
     * Contants used for the forward and inverse transform for the eliptical
145
     * case of the Transverse Mercator.
146
     */
147
    private static final double FC1= 1.00000000000000000000000,  // 1/1
148
                                FC2= 0.50000000000000000000000,  // 1/2
149
                                FC3= 0.16666666666666666666666,  // 1/6
150
                                FC4= 0.08333333333333333333333,  // 1/12
151
                                FC5= 0.05000000000000000000000,  // 1/20
152
                                FC6= 0.03333333333333333333333,  // 1/30
153
                                FC7= 0.02380952380952380952380,  // 1/42
154
                                FC8= 0.01785714285714285714285;  // 1/56
155

    
156
    /**
157
     * Relative iteration precision used in the <code>mlfn<code> method. This 
158
     * overrides the value in the MapProjection class.
159
     */
160
    private static final double TOL = 1E-11;
161

    
162
    /**
163
     * Constructs a new map projection from the supplied parameters.
164
     *
165
     * @param  parameters The parameter values in standard units.
166
     * @throws ParameterNotFoundException if a mandatory parameter is missing.
167
     */
168
    protected TransverseMercator(final ParameterValueGroup parameters)
169
            throws ParameterNotFoundException
170
    {
171
        //Fetch parameters 
172
        super(parameters);
173
        
174
        //  Compute constants
175
        esp = excentricitySquared / (1.0 - excentricitySquared);
176
          
177
        double t;
178
        en0 = C00 - excentricitySquared  *  (C02 + excentricitySquared  * 
179
             (C04 + excentricitySquared  *  (C06 + excentricitySquared  * C08)));
180
        en1 =       excentricitySquared  *  (C22 - excentricitySquared  *
181
             (C04 + excentricitySquared  *  (C06 + excentricitySquared  * C08)));
182
        en2 =  (t = excentricitySquared  *         excentricitySquared) * 
183
             (C44 - excentricitySquared  *  (C46 + excentricitySquared  * C48));
184
        en3 = (t *= excentricitySquared) *  (C66 - excentricitySquared  * C68);
185
        en4 =   t * excentricitySquared  *  C88;
186
        ml0 = mlfn(latitudeOfOrigin, Math.sin(latitudeOfOrigin), Math.cos(latitudeOfOrigin));
187
    }
188
    
189
    /**
190
     * {@inheritDoc}
191
     */
192
    public ParameterDescriptorGroup getParameterDescriptors() {
193
        return Provider.PARAMETERS;
194
    }
195
    
196
    /**
197
     * Transforms the specified (<var>x</var>,<var>y</var>) coordinate (units in radians)
198
     * and stores the result in {@code ptDst} (linear distance on a unit sphere).
199
     */
200
    protected Point2D transformNormalized(double x, double y, Point2D ptDst) 
201
            throws ProjectionException 
202
    {
203
        double sinphi = Math.sin(y);
204
        double cosphi = Math.cos(y);
205
        
206
        double t = (Math.abs(cosphi)>EPS) ? sinphi/cosphi : 0;
207
        t *= t;
208
        double al = cosphi*x;
209
        double als = al*al;
210
        al /= Math.sqrt(1.0 - excentricitySquared * sinphi*sinphi);
211
        double n = esp * cosphi*cosphi;
212

    
213
        /* NOTE: meridinal distance at latitudeOfOrigin is always 0 */
214
        y = (mlfn(y, sinphi, cosphi) - ml0 + 
215
            sinphi*al*x*
216
            FC2 * ( 1.0 +
217
            FC4 * als * (5.0 - t + n*(9.0 + 4.0*n) +
218
            FC6 * als * (61.0 + t * (t - 58.0) + n*(270.0 - 330.0*t) +
219
            FC8 * als * (1385.0 + t * ( t*(543.0 - t) - 3111.0))))));
220
        
221
        x = al*(FC1 + FC3 * als*(1.0 - t + n +
222
            FC5 * als * (5.0 + t*(t - 18.0) + n*(14.0 - 58.0*t) +
223
            FC7 * als * (61.0+ t*(t*(179.0 - t) - 479.0 )))));
224
               
225
        if (ptDst != null) {
226
            ptDst.setLocation(x,y);
227
            return ptDst;
228
        }
229
        return new Point2D.Double(x,y);          
230
    }
231
    
232
    /**
233
     * Transforms the specified (<var>x</var>,<var>y</var>) coordinate
234
     * and stores the result in {@code ptDst}.
235
     */
236
    protected Point2D inverseTransformNormalized(double x, double y, Point2D ptDst) 
237
            throws ProjectionException 
238
    {
239
        double phi = inv_mlfn(ml0 + y);
240
        
241
        if (Math.abs(phi) >= (Math.PI/2)) {
242
            y = y<0.0 ? -(Math.PI/2) : (Math.PI/2);
243
            x = 0.0;
244
        } else {
245
            double sinphi = Math.sin(phi);
246
            double cosphi = Math.cos(phi);
247
            double t = (Math.abs(cosphi) > EPS) ? sinphi/cosphi : 0.0;
248
            double n = esp * cosphi*cosphi;
249
            double con = 1.0 - excentricitySquared * sinphi*sinphi;
250
            double d = x*Math.sqrt(con);
251
            con *= t;
252
            t *= t;
253
            double ds = d*d;
254
            
255
            y = phi - (con*ds / (1.0 - excentricitySquared)) *
256
                FC2 * (1.0 - ds *
257
                FC4 * (5.0 + t*(3.0 - 9.0*n) + n*(1.0 - 4*n) - ds *
258
                FC6 * (61.0 + t*(90.0 - 252.0*n + 45.0*t) + 46.0*n - ds *
259
                FC8 * (1385.0 + t*(3633.0 + t*(4095.0 + 1574.0*t))))));
260
            
261
            x = d*(FC1 - ds * FC3 * (1.0 + 2.0*t + n -
262
                ds*FC5*(5.0 + t*(28.0 + 24* t + 8.0*n) + 6.0*n -
263
                ds*FC7*(61.0 + t*(662.0 + t*(1320.0 + 720.0*t))))))/cosphi;
264
        }
265
        
266
        if (ptDst != null) {
267
            ptDst.setLocation(x,y);
268
            return ptDst;
269
        }
270
        return new Point2D.Double(x,y);        
271
    }
272
    
273
    /**
274
     * {@inheritDoc}
275
     */
276
    protected double getToleranceForAssertions(final double longitude, final double latitude) {
277
        if (Math.abs(longitude - centralMeridian) > 0.26) {   //15 degrees
278
            // When far from the valid area, use a larger tolerance.
279
            return 2.5;
280
        } else if (Math.abs(longitude - centralMeridian) > 0.22) {  //12.5 degrees
281
            return 1.0;
282
        } else if (Math.abs(longitude - centralMeridian) > 0.17) {  //10 degrees
283
            return 0.5;
284
        }
285
        // a normal tolerance
286
        return 1E-6;
287
    }
288
    
289

    
290
    /**
291
     * Provides the transform equations for the spherical case of the
292
     * TransverseMercator projection.
293
     * 
294
     * @version $Id: TransverseMercator.java 12186 2007-06-18 08:46:09Z jlgomez $
295
     * @author Andr Gosselin
296
     * @author Martin Desruisseaux
297
     * @author Rueben Schulz
298
     */
299
    private static final class Spherical extends TransverseMercator {
300
        /**
301
         * Constructs a new map projection from the suplied parameters.
302
         *
303
         * @param  parameters The parameter values in standard units.
304
         * @throws ParameterNotFoundException if a mandatory parameter is missing.
305
         */
306
        protected Spherical(final ParameterValueGroup parameters)
307
                throws ParameterNotFoundException
308
        {
309
            super(parameters);
310
            //assert isSpherical;
311
        }
312
        
313
        /**
314
         * {@inheritDoc}
315
         */
316
        protected Point2D transformNormalized(double x, double y, Point2D ptDst)
317
                throws ProjectionException 
318
        {
319
            // Compute using ellipsoidal formulas, for comparaison later.
320
            double normalizedLongitude = x;
321
            //assert (ptDst = super.transformNormalized(x, y, ptDst)) != null;
322
                    
323
            double cosphi = Math.cos(y);
324
            double b = cosphi * Math.sin(x);
325
            if (Math.abs(Math.abs(b) - 1.0) <= EPS) {
326
                throw new ProjectionException(Resources.format(ResourceKeys.ERROR_VALUE_TEND_TOWARD_INFINITY));
327
            }
328
            
329
            //Using Snyder's equation for calculating y, instead of the one used in Proj4
330
            //poential problems when y and x = 90 degrees, but behaves ok in tests   
331
            y = Math.atan2(Math.tan(y),Math.cos(x)) - latitudeOfOrigin;   /* Snyder 8-3 */
332
            x = 0.5 * Math.log((1.0+b)/(1.0-b));    /* Snyder 8-1 */
333

    
334
            //assert Math.abs(ptDst.getX()-x) <= getToleranceForSphereAssertions(normalizedLongitude,0) : ptDst.getX()-x;
335
            //assert Math.abs(ptDst.getY()-y) <= getToleranceForSphereAssertions(normalizedLongitude,0) : ptDst.getY()-y;
336
            if (ptDst != null) {
337
                ptDst.setLocation(x,y);
338
                return ptDst;
339
            }
340
            return new Point2D.Double(x,y);
341
        }        
342
        
343
        /**
344
         * {@inheritDoc}
345
         */
346
        protected Point2D inverseTransformNormalized(double x, double y, Point2D ptDst) 
347
                throws ProjectionException 
348
        {
349
            // Compute using ellipsoidal formulas, for comparaison later.
350
            //assert (ptDst = super.inverseTransformNormalized(x, y, ptDst)) != null;
351
            
352
            double t = Math.exp(x);
353
            double sinhX = 0.5 * (t-1.0/t);                //sinh(x)
354
            double cosD = Math.cos(latitudeOfOrigin + y);
355
            double phi = Math.asin(Math.sqrt((1.0-cosD*cosD)/(1.0+sinhX*sinhX)));
356
            // correct for the fact that we made everything positive using sqrt(x*x)
357
            y = ((y + latitudeOfOrigin)<0.0) ? -phi : phi;   
358
            x = (Math.abs(sinhX)<=EPS && Math.abs(cosD)<=EPS) ? 
359
                    0.0 : Math.atan2(sinhX,cosD);
360
           
361
            //assert Math.abs(ptDst.getX()-x) <= getToleranceForSphereAssertions(x,0) : ptDst.getX()-x;
362
            //assert Math.abs(ptDst.getY()-y) <= getToleranceForSphereAssertions(x,0) : ptDst.getY()-y;
363
            if (ptDst != null) {
364
                ptDst.setLocation(x,y);
365
                return ptDst;
366
            }
367
            return new Point2D.Double(x,y);
368
        }
369
        
370
        /*
371
         *  Allow a larger tolerance when comparing spherical to elliptical calculations
372
         *  when we are far from the central meridian (elliptical calculations are
373
         *  not valid here).
374
         *
375
         *  longitude here is standardized (in radians) and centralMeridian has 
376
         *  already been removed from it.
377
         */
378
        protected double getToleranceForSphereAssertions(final double longitude, final double latitude) {
379
            if (Math.abs(Math.abs(longitude)- Math.PI/2.0) < TOL) {  //90 degrees
380
                // elliptical equations are at their worst here
381
                return 1E18;
382
            }
383
            if (Math.abs(longitude) > 0.26) {   //15 degrees
384
                // When far from the valid area, use a very larger tolerance.          
385
                return 1000000;
386
            }
387
            // a normal tolerance
388
            return 1E-6;
389
        }
390
    }
391
    
392
    
393
    /**
394
     * Calculates the meridian distance. This is the distance along the central 
395
     * meridian from the equator to {@code phi}. Accurate to < 1e-5 meters 
396
     * when used in conjuction with typical major axis values.
397
     *
398
     * @param phi latitude to calculate meridian distance for.
399
     * @param sphi sin(phi).
400
     * @param cphi cos(phi).
401
     * @return meridian distance for the given latitude.
402
     */
403
    private final double mlfn(final double phi, double sphi, double cphi) {        
404
        cphi *= sphi;
405
        sphi *= sphi;
406
        return en0 * phi - cphi *
407
              (en1 + sphi *
408
              (en2 + sphi *
409
              (en3 + sphi *
410
              (en4))));
411
    }
412
    
413
    /**
414
     * Calculates the latitude ({@code phi}) from a meridian distance.
415
     * Determines phi to TOL (1e-11) radians, about 1e-6 seconds.
416
     * 
417
     * @param arg meridian distance to calulate latitude for.
418
     * @return the latitude of the meridian distance.
419
     * @throws ProjectionException if the itteration does not converge.
420
     */
421
    private final double inv_mlfn(double arg) throws ProjectionException {
422
        double s, t, phi, k = 1.0/(1.0 - excentricitySquared);
423
        int i;
424
        phi = arg;
425
        for (i=MAX_ITER; true;) { // rarely goes over 5 iterations
426
            if (--i < 0) {
427
                throw new ProjectionException(Resources.format(ResourceKeys.ERROR_NO_CONVERGENCE));
428
            }
429
            s = Math.sin(phi);
430
            t = 1.0 - excentricitySquared * s * s;
431
            t = (mlfn(phi, s, Math.cos(phi)) - arg) * (t * Math.sqrt(t)) * k;
432
            phi -= t;
433
            if (Math.abs(t) < TOL) {
434
                return phi;
435
            }
436
        }
437
    }
438
    
439
    /**
440
     * Convenience method computing the zone code from the central meridian.
441
     * Information about zones convention must be specified in argument. Two
442
     * widely set of arguments are of Universal Transverse Mercator (UTM) and
443
     * Modified Transverse Mercator (MTM) projections:<br>
444
     * <p>
445
     *
446
     * UTM projection (zones numbered from 1 to 60):<br>
447
     * <p>
448
     *        {@code getZone(-177, 6);}<br>
449
     * <p>
450
     * MTM projection (zones numbered from 1 to 120):<br>
451
     * <p>
452
     *        {@code getZone(-52.5, -3);}<br>
453
     *
454
     * @param  centralLongitudeZone1 Longitude in the middle of zone 1, in decimal degrees
455
     *         relative to Greenwich. Positive longitudes are toward east, and negative
456
     *         longitudes toward west.
457
     * @param  zoneWidth Number of degrees of longitudes in one zone. A positive value
458
     *         means that zones are numbered from west to east (i.e. in the direction of
459
     *         positive longitudes). A negative value means that zones are numbered from
460
     *         east to west.
461
     * @return The zone number. First zone is numbered 1.
462
     */
463
    private int getZone(final double centralLongitudeZone1, final double zoneWidth) {
464
        final double zoneCount = Math.abs(360/zoneWidth);
465
        double t;
466
        t  = centralLongitudeZone1 - 0.5*zoneWidth; // Longitude at the beginning of the first zone.
467
        t  = Math.toDegrees(centralMeridian) - t;   // Degrees of longitude between the central longitude and longitude 1.
468
        t  = Math.floor(t/zoneWidth + EPS);         // Number of zones between the central longitude and longitude 1.
469
        t -= zoneCount*Math.floor(t/zoneCount);     // If negative, bring back to the interval 0 to (zoneCount-1).
470
        return ((int) t)+1;
471
    }
472
    
473
    /**
474
     * Convenience method returning the meridian in the middle of current zone. This meridian is
475
     * typically the central meridian. This method may be invoked to make sure that the central
476
     * meridian is correctly set.
477
     *
478
     * @param  centralLongitudeZone1 Longitude in the middle of zone 1, in decimal degrees
479
     *         relative to Greenwich. Positive longitudes are toward east, and negative
480
     *         longitudes toward west.
481
     * @param  zoneWidth Number of degrees of longitudes in one zone. A positive value
482
     *         means that zones are numbered from west to east (i.e. in the direction of
483
     *         positive longitudes). A negative value means that zones are numbered from
484
     *         east to west.
485
     * @return The central meridian.
486
     */
487
    private double getCentralMedirian(final double centralLongitudeZone1, final double zoneWidth) {
488
        double t;
489
        t  = centralLongitudeZone1 + (getZone(centralLongitudeZone1, zoneWidth)-1)*zoneWidth;
490
        t -= 360*Math.floor((t+180)/360); // Bring back into [-180..+180] range.
491
        return t;
492
    }
493

    
494
    /**
495
     * Convenience method computing the zone code from the central meridian.
496
     *
497
     * @return The zone number, using the scalefactor and false easting 
498
     *         to decide if this is a UTM or MTM case. Returns 0 if the 
499
     *         case of the projection cannot be determined.
500
     */
501
    public int getZone() {
502
        // UTM
503
        if (scaleFactor == 0.9996 && falseEasting == 500000.0) {
504
            return getZone(-177, 6);
505
        }
506
        // MTM
507
        if (scaleFactor == 0.9999 && falseEasting == 304800.0){
508
            return getZone(-52.5, -3);
509
        }
510
        // unknown
511
        throw new IllegalStateException(Resources.format(ResourceKeys.ERROR_UNKNOW_TYPE_$1));//.UNKNOW_PROJECTION_TYPE));
512
    }
513

    
514
    /**
515
     * Convenience method returning the meridian in the middle of current zone. This meridian is
516
     * typically the central meridian. This method may be invoked to make sure that the central
517
     * meridian is correctly set.
518
     *
519
     * @return The central meridian, using the scalefactor and false easting 
520
     *         to decide if this is a UTM or MTM case. Returns {@link Double#NaN}
521
     *         if the case of the projection cannot be determined.
522
     */
523
    public double getCentralMeridian() {
524
        // UTM
525
        if (scaleFactor == 0.9996 && falseEasting == 500000.0) {
526
            return getCentralMedirian(-177, 6);
527
        }
528
        // MTM
529
        if (scaleFactor == 0.9999 && falseEasting == 304800.0){
530
            return getCentralMedirian(-52.5, -3);
531
        }
532
        // unknown
533
        throw new IllegalStateException(Resources.format(ResourceKeys.ERROR_UNKNOW_TYPE_$1));//.UNKNOW_PROJECTION_TYPE));
534
    }
535

    
536
    /**
537
     * Returns a hash value for this projection.
538
     */
539
    public int hashCode() { 
540
        final long code = Double.doubleToLongBits(ml0);
541
        return ((int)code ^ (int)(code >>> 32)) + 37*super.hashCode();
542
    }
543

    
544
    /**
545
     * Compares the specified object with this map projection for equality.
546
     */
547
    public boolean equals(final Object object) {
548
        if (object == this) {
549
            // Slight optimization
550
            return true;
551
        }
552
        // Relevant parameters are already compared in MapProjection
553
        return super.equals(object);
554
    }
555
    
556
    
557
    
558
    
559
    //////////////////////////////////////////////////////////////////////////////////////////
560
    //////////////////////////////////////////////////////////////////////////////////////////
561
    ////////                                                                          ////////
562
    ////////                                 PROVIDER                                 ////////
563
    ////////                                                                          ////////
564
    //////////////////////////////////////////////////////////////////////////////////////////
565
    //////////////////////////////////////////////////////////////////////////////////////////
566

    
567
    /**
568
     * The {@link org.geotools.referencing.operation.MathTransformProvider}
569
     * for a {@link TransverseMercator} projection.
570
     *
571
     * @see org.geotools.referencing.operation.DefaultMathTransformFactory
572
     *
573
     * @since 2.1
574
     * @version $Id: TransverseMercator.java 12186 2007-06-18 08:46:09Z jlgomez $
575
     * @author Martin Desruisseaux
576
     * @author Rueben Schulz
577
     */
578
    public static class Provider extends AbstractProvider {
579
        /**
580
         * Returns a descriptor group for the specified parameters.
581
         */
582
        static ParameterDescriptorGroup createDescriptorGroup(final Identifier[] identifiers) {
583
            return createDescriptorGroup(identifiers, new ParameterDescriptor[] {
584
                SEMI_MAJOR,       SEMI_MINOR,
585
                CENTRAL_MERIDIAN, LATITUDE_OF_ORIGIN,
586
                SCALE_FACTOR,     FALSE_EASTING,
587
                FALSE_NORTHING
588
            });
589
        }
590

    
591
        /**
592
         * The parameters group.
593
         */
594
        static final ParameterDescriptorGroup PARAMETERS = createDescriptorGroup(new NamedIdentifier[] {
595
                new NamedIdentifier(CitationImpl.OGC,      "Transverse_Mercator"),
596
                new NamedIdentifier(CitationImpl.EPSG,     "Transverse Mercator"),
597
                new NamedIdentifier(CitationImpl.EPSG,     "Gauss-Kruger"),
598
                new NamedIdentifier(CitationImpl.EPSG,     "9807"),
599
                new NamedIdentifier(CitationImpl.GEOTIFF,  "CT_TransverseMercator"),
600
                new NamedIdentifier(CitationImpl.ESRI,     "Transverse_Mercator"),
601
                new NamedIdentifier(CitationImpl.ESRI,     "Gauss_Kruger")//,
602
                //new NamedIdentifier(CitationImpl.GEOTOOLS, Vocabulary.formatInternational(
603
                //                    VocabularyKeys.TRANSVERSE_MERCATOR_PROJECTION))
604
            });
605

    
606
        /**
607
         * Constructs a new provider.
608
         */
609
        public Provider() {
610
            super(PARAMETERS);
611
        }
612

    
613
        /**
614
         * Constructs a new provider with the specified parameters.
615
         */
616
        Provider(final ParameterDescriptorGroup descriptor) {
617
            super(descriptor);
618
        }
619

    
620
        /**
621
         * Returns the operation type for this map projection.
622
         */
623
        protected Class getOperationType() {
624
            return CylindricalProjection.class;
625
        }
626

    
627
        /**
628
         * Creates a transform from the specified group of parameter values.
629
         *
630
         * @param  parameters The group of parameter values.
631
         * @return The created math transform.
632
         * @throws ParameterNotFoundException if a required parameter was not found.
633
         */
634
        public MathTransform createMathTransform(final ParameterValueGroup parameters)
635
                throws ParameterNotFoundException
636
        {
637
            if (isSpherical(parameters)) {
638
                return new Spherical(parameters);
639
            } else {
640
                return new TransverseMercator(parameters);
641
            }
642
        }
643
    }
644

    
645
    /**
646
     * The {@link org.geotools.referencing.operation.MathTransformProvider} for a South Orientated
647
     * {@link TransverseMercator} projection (EPSG code 9808). Note that at the contrary of what
648
     * this class name suggest, the projected coordinates are still increasing toward North. This
649
     * is because all {@link MapProjection}s must complies with standard axis orientations. The
650
     * real axis flip is performed by the CRS framework outside this package.
651
     * See "<cite>Axis units and orientation</cite>" in
652
     * {@linkplain org.geotools.referencing.operation.projection package description} for details.
653
     * <p>
654
     * The usual Transverse Mercator formulas are:
655
     * <p>
656
     * <ul>
657
     *   <li><var>easting</var> = <var>{@linkplain #falseEasting false easting}</var> + <var>px</var></li>
658
     *   <li><var>northing</var> = <var>{@linkplain #falseNorthing false northing}</var> + <var>py</var></li>
659
     * </ul>
660
     * <p>
661
     * The Transverse Mercator South Orientated Projection formulas are:
662
     * <p>
663
     * <ul>
664
     *   <li><var>westing</var> = <var>{@linkplain #falseEasting false easting}</var> - <var>px</var></li>
665
     *   <li><var>southing</var> = <var>{@linkplain #falseNorthing false northing}</var> - <var>py</var></li>
666
     * </ul>
667
     * <p>
668
     * Where the <var>px</var> and <var>py</var> terms are the same in both cases.
669
     * Transforms created by this provider actually computes
670
     * (<var>easting</var>,<var>northing</var>) = (-<var>westing</var>,-<var>southing</var>).
671
     * This is equivalent to a {@link TransverseMercator} projection with
672
     * {@link #falseEasting falseEasting} and {@link #falseNorthing falseNorthing} sign reverted.
673
     * This operation is implemented as a concatenation of a North-orientated transverse mercator
674
     * projection with an affine transform for (<var>false easting</var>,<var>false northing</var>)
675
     * correction.
676
     *
677
     * @since 2.2
678
     * @version $Id: TransverseMercator.java 12186 2007-06-18 08:46:09Z jlgomez $
679
     * @author Martin Desruisseaux
680
     */
681
    public static class Provider_SouthOrientated extends Provider {
682
        /**
683
         * Constructs a new provider.
684
         */
685
        public Provider_SouthOrientated() {
686
            super(createDescriptorGroup(new NamedIdentifier[] {
687
                new NamedIdentifier(CitationImpl.EPSG, "Transverse Mercator (South Orientated)"),
688
                new NamedIdentifier(CitationImpl.EPSG, "9808")
689
            }));
690
        }
691

    
692
        /**
693
         * Creates a transform from the specified group of parameter values.
694
         *
695
         * @param  parameters The group of parameter values.
696
         * @return The created math transform.
697
         * @throws ParameterNotFoundException if a required parameter was not found.
698
         */
699
        public MathTransform createMathTransform(final ParameterValueGroup parameters)
700
                throws ParameterNotFoundException
701
        {
702
            final MapProjection projection = (MapProjection) super.createMathTransform(parameters);
703
            if (projection.falseEasting==0 && projection.falseNorthing==0) {
704
                return projection;
705
            }
706
            final AffineTransform step = AffineTransform.getTranslateInstance(
707
                    -2*projection.falseEasting, -2*projection.falseNorthing);
708
            return ConcatenatedTransform.create(ProjectiveTransform.create(step), projection);
709
        }
710
    }
711
}