001    //$HeadURL: $
002    /*----------------    FILE HEADER  ------------------------------------------
003     This file is part of deegree.
004     Copyright (C) 2001-2008 by:
005     Department of Geography, University of Bonn
006     http://www.giub.uni-bonn.de/deegree/
007     lat/lon GmbH
008     http://www.lat-lon.de
009    
010     This library is free software; you can redistribute it and/or
011     modify it under the terms of the GNU Lesser General Public
012     License as published by the Free Software Foundation; either
013     version 2.1 of the License, or (at your option) any later version.
014     This library is distributed in the hope that it will be useful,
015     but WITHOUT ANY WARRANTY; without even the implied warranty of
016     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
017     Lesser General Public License for more details.
018     You should have received a copy of the GNU Lesser General Public
019     License along with this library; if not, write to the Free Software
020     Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
021     Contact:
022    
023     Andreas Poth
024     lat/lon GmbH
025     Aennchenstr. 19
026     53177 Bonn
027     Germany
028     E-Mail: poth@lat-lon.de
029    
030     Prof. Dr. Klaus Greve
031     Department of Geography
032     University of Bonn
033     Meckenheimer Allee 166
034     53115 Bonn
035     Germany
036     E-Mail: greve@giub.uni-bonn.de
037     ---------------------------------------------------------------------------*/
038    
039    package org.deegree.crs.projections;
040    
041    import javax.vecmath.Point2d;
042    
043    import org.deegree.crs.Identifiable;
044    import org.deegree.crs.components.Datum;
045    import org.deegree.crs.components.Ellipsoid;
046    import org.deegree.crs.components.PrimeMeridian;
047    import org.deegree.crs.components.Unit;
048    import org.deegree.crs.coordinatesystems.GeographicCRS;
049    import org.deegree.crs.exceptions.CRSException;
050    import org.deegree.i18n.Messages;
051    
052    /**
053     * Map <code>conversion</code> is the process of changing the map grid coordinates (usually, but
054     * not always, Eastings & Northings) of a Projected Coordinate Reference System to its corresponding
055     * geographical coordinates (Latitude & Longitude) or vice versa.
056     * <p>
057     * A projection is conformal if an infinitesimal small perfect circle on the earth's surface results
058     * in an infinitesimal small projected perfect circle (an ellisoid with no eccentrity). In other
059     * words, the relative local angles about every point on the map are shown correctly.
060     * </p>
061     * <p>
062     * An equal area projection can be best explained with a coin (Snyder), a coin (of any size) covers
063     * exactly the same area of the actual earth as the same coin on any other part of the map. This can
064     * only be done by distorting shape, scale and and angles of the original earth's layout.
065     * </p>
066     * 
067     * @author <a href="mailto:bezema@lat-lon.de">Rutger Bezema</a>
068     * 
069     * @author last edited by: $Author:$
070     * 
071     * @version $Revision:$, $Date:$
072     * 
073     */
074    
075    public abstract class Projection extends Identifiable {
076    
077        private final boolean conformal;
078    
079        private boolean equalArea;
080    
081        private double scale;
082    
083        /**
084         * the scale*semimajor-axis, often revered to as R*k_0 in Snyder.
085         */
086        private double scaleFactor;
087    
088        private final double falseNorthing;
089    
090        private double falseEasting;
091    
092        private final Point2d naturalOrigin;
093    
094        // Values gotten from the natural origin
095        private double projectionLatitude;
096    
097        private double projectionLongitude;
098    
099        // the sin of the projection latitude
100        private double sinphi0;
101    
102        // the cos of the projection latitude
103        private double cosphi0;
104    
105        private final Unit units;
106    
107        private final GeographicCRS geographicCRS;
108    
109        private boolean isSpherical;
110    
111        /**
112         * Creates a Projection. Causian, the given natural origin should be given in radians rather
113         * then degrees.
114         * 
115         * @param geographicCRS
116         *            which this projection uses.
117         * @param falseNorthing
118         *            in given units
119         * @param falseEasting
120         *            in given units
121         * @param naturalOrigin
122         *            in radians longitude, latitude.
123         * @param units
124         *            of the map projection
125         * @param scale
126         *            at the prime meridian (eg. 0.9996 for UTM)
127         * @param conformal
128         *            if the projection is conformal
129         * @param equalArea
130         *            if the projection result in an equal area map
131         * @param identifiers
132         * @param names
133         * @param versions
134         * @param descriptions
135         * @param areasOfUse
136         */
137        public Projection( GeographicCRS geographicCRS, double falseNorthing, double falseEasting, Point2d naturalOrigin,
138                           Unit units, double scale, boolean conformal, boolean equalArea, String[] identifiers,
139                           String[] names, String[] versions, String[] descriptions, String[] areasOfUse ) {
140            super( identifiers, names, versions, descriptions, areasOfUse );
141            this.scale = scale;
142            this.conformal = conformal;
143            this.equalArea = equalArea;
144            this.geographicCRS = geographicCRS;
145            this.falseNorthing = falseNorthing;
146            this.falseEasting = falseEasting;
147            this.units = units;
148    
149            checkForNullObject( geographicCRS, Messages.getMessage( "CRS_PARAMETER_NOT_NULL", "Projection", "geographic" ) );
150            checkForNullObject( geographicCRS.getGeodeticDatum(), Messages.getMessage( "CRS_PARAMETER_NOT_NULL",
151                                                                                       "Projection", "geographic.datum" ) );
152            checkForNullObject( geographicCRS.getGeodeticDatum().getEllipsoid(),
153                                Messages.getMessage( "CRS_PARAMETER_NOT_NULL", "Projection", "geographic.datum.ellipsoid" ) );
154            checkForNullObject( naturalOrigin,
155                                Messages.getMessage( "CRS_PARAMETER_NOT_NULL", "Projection", "naturalOrigin" ) );
156            checkForNullObject( units, Messages.getMessage( "CRS_PARAMETER_NOT_NULL", "Projection", "units" ) );
157    
158            this.scaleFactor = scale * getSemiMajorAxis();
159    
160            this.naturalOrigin = new Point2d( ProjectionUtils.normalizeLongitude( naturalOrigin.x ),
161                                              ProjectionUtils.normalizeLatitude( naturalOrigin.y ) );
162    
163            // uses different library
164            // this.projectionLongitude = this.naturalOrigin.getX();
165            // this.projectionLatitude = this.naturalOrigin.getY();
166            this.projectionLongitude = this.naturalOrigin.x;
167            this.projectionLatitude = this.naturalOrigin.y;
168    
169            sinphi0 = Math.sin( projectionLatitude );
170            cosphi0 = Math.cos( projectionLatitude );
171    
172            isSpherical = geographicCRS.getGeodeticDatum().getEllipsoid().getEccentricity() < 0.0000001;
173        }
174    
175        /**
176         * The actual transform method doing a projection from geographic coordinates to map
177         * coordinates.
178         * 
179         * @param lambda
180         *            the longitude
181         * @param phi
182         *            the latitude
183         * @return the projected Point or Point(Double.NAN, Double.NAN) if an error occurred.
184         * @throws CRSException
185         */
186        public abstract Point2d doProjection( double lambda, double phi )
187                                throws CRSException;
188    
189        /**
190         * Do an inverse projection from projected (map) coordinates to geogpraphic coordinates.
191         * 
192         * @param x
193         *            coordinate on the map
194         * @param y
195         *            coordinate on the map
196         * @return the projected Point with x = lambda and y = phi, or Point(Double.NAN, Double.NAN) if
197         *         an error occurred.
198         * @throws CRSException
199         */
200        public abstract Point2d doInverseProjection( double x, double y )
201                                throws CRSException;
202    
203        /**
204         * @return A deegree specific name which will be used for exportation of a projection.
205         */
206        public abstract String getDeegreeSpecificName();
207    
208        /**
209         * @return true if the projection projects conformal.
210         */
211        public final boolean isConformal() {
212            return conformal;
213        }
214    
215        /**
216         * @return true if the projection is projects equal Area.
217         */
218        public final boolean isEqualArea() {
219            return equalArea;
220        }
221    
222        /**
223         * @return the scale.
224         */
225        public final double getScale() {
226            return scale;
227        }
228    
229        /**
230         * Sets the old scale to the given scale, also adjusts the scaleFactor.
231         * 
232         * @param scale
233         *            the new scale
234         */
235        public void setScale( double scale ) {
236            this.scale = scale;
237            this.scaleFactor = scale * getSemiMajorAxis();
238        }
239    
240        /**
241         * @return the scale*semimajor-axis, often revered to as R*k_0 in Snyder.
242         */
243        public final double getScaleFactor() {
244            return scaleFactor;
245        }
246    
247        /**
248         * @return the datum.
249         */
250        public final Datum getDatum() {
251            return geographicCRS.getGeodeticDatum();
252        }
253    
254        /**
255         * @return the falseEasting.
256         */
257        public final double getFalseEasting() {
258            return falseEasting;
259        }
260    
261        /**
262         * sets the false easting to given value. (Used in for example transverse mercator, while
263         * setting the utm zone).
264         * 
265         * @param newFalseEasting
266         *            the new false easting parameter.
267         */
268        public void setFalseEasting( double newFalseEasting ) {
269            this.falseEasting = newFalseEasting;
270        }
271    
272        /**
273         * @return the falseNorthing.
274         */
275        public final double getFalseNorthing() {
276            return falseNorthing;
277        }
278    
279        /**
280         * @return the naturalOrigin.
281         */
282        public final Point2d getNaturalOrigin() {
283            return naturalOrigin;
284        }
285    
286        /**
287         * @return the units.
288         */
289        public final Unit getUnits() {
290            return units;
291        }
292    
293        /**
294         * @return the primeMeridian of the datum.
295         */
296        public final PrimeMeridian getPrimeMeridian() {
297            return geographicCRS.getGeodeticDatum().getPrimeMeridian();
298        }
299    
300        /**
301         * @return the ellipsoid of the datum.
302         */
303        public final Ellipsoid getEllipsoid() {
304            return geographicCRS.getGeodeticDatum().getEllipsoid();
305        }
306    
307        /**
308         * @return the eccentricity of the ellipsoid of the datum.
309         */
310        public final double getEccentricity() {
311            return geographicCRS.getGeodeticDatum().getEllipsoid().getEccentricity();
312        }
313    
314        /**
315         * @return the eccentricity of the ellipsoid of the datum.
316         */
317        public final double getSquaredEccentricity() {
318            return geographicCRS.getGeodeticDatum().getEllipsoid().getSquaredEccentricity();
319        }
320    
321        /**
322         * @return the semiMajorAxis (a) of the ellipsoid of the datum.
323         */
324        public final double getSemiMajorAxis() {
325            return geographicCRS.getGeodeticDatum().getEllipsoid().getSemiMajorAxis();
326        }
327    
328        /**
329         * @return the semiMinorAxis (a) of the ellipsoid of the datum.
330         */
331        public final double getSemiMinorAxis() {
332            return geographicCRS.getGeodeticDatum().getEllipsoid().getSemiMinorAxis();
333        }
334    
335        /**
336         * @return true if the ellipsoid of the datum is a sphere and not an ellipse.
337         */
338        public final boolean isSpherical() {
339            return isSpherical;
340        }
341    
342        /**
343         * Set the new projection latitude and calculate the sinphi0 and cosPhi0 accordingly
344         * 
345         * @param projectionLatitude
346         *            the new latitude.
347         */
348        public final void setProjectionLatitude( double projectionLatitude ) {
349            this.projectionLatitude = projectionLatitude;
350            sinphi0 = Math.sin( projectionLatitude );
351            cosphi0 = Math.cos( projectionLatitude );
352        }
353    
354        /**
355         * @param projectionLongitude
356         *            the new projection longitude
357         */
358        public final void setProjectionLongitude( double projectionLongitude ) {
359            this.projectionLongitude = projectionLongitude;
360        }
361    
362        /**
363         * @return the projectionLatitude also known as central-latitude or latitude-of-origin, in
364         *         Snyder referenced as phi_1 for azimuthal, phi_0 for other projections.
365         */
366        public final double getProjectionLatitude() {
367            return projectionLatitude;
368        }
369    
370        /**
371         * @return the projectionLongitude also known as projection-meridian or central-meridian, in
372         *         Snyder referenced as lambda_0
373         */
374        public final double getProjectionLongitude() {
375            return projectionLongitude;
376        }
377    
378        /**
379         * @return the sinphi0, the sine of the projection latitude
380         */
381        public final double getSinphi0() {
382            return sinphi0;
383        }
384    
385        /**
386         * @return the cosphi0, the cosine of the projection latitude
387         */
388        public final double getCosphi0() {
389            return cosphi0;
390        }
391    
392        /**
393         * @return the geographicCRS.
394         */
395        public final GeographicCRS getGeographicCRS() {
396            return geographicCRS;
397        }
398    
399        @Override
400        public boolean equals( Object other ) {
401            if ( other != null && other instanceof Projection ) {
402                final Projection that = (Projection) other;
403                return this.units.equals( that.units )
404                       && Math.abs( ( this.projectionLatitude - that.projectionLatitude ) ) < ProjectionUtils.EPS11
405                       && Math.abs( ( this.projectionLongitude - that.projectionLongitude ) ) < ProjectionUtils.EPS11
406                       && Math.abs( ( this.falseNorthing - that.falseNorthing ) ) < ProjectionUtils.EPS11
407                       && Math.abs( ( this.falseEasting - that.falseEasting ) ) < ProjectionUtils.EPS11
408                       && Math.abs( ( this.scale - that.scale ) ) < ProjectionUtils.EPS11
409                       && ( this.conformal == that.conformal ) && ( this.equalArea == that.equalArea );
410    
411                // && Math.abs( (this.projectionLatitude - that.projectionLatitude) ) >
412                // ProjectionUtils.EPS11 &&
413                // Math.abs( (this.projectionLongitude - that.projectionLongitude) ) >
414                // ProjectionUtils.EPS11 &&
415                // Math.abs( (this.falseNorthing - that.falseNorthing) ) > ProjectionUtils.EPS11 &&
416                // Math.abs( (this.falseEasting - that.falseEasting) ) > ProjectionUtils.EPS11 &&
417                // Math.abs( (this.scale- that.scale) ) > ProjectionUtils.EPS11 &&
418                // ( this.conformal && that.conformal) && (this.equalArea && that.equalArea );
419            }
420            return false;
421        }
422    
423        @Override
424        public String toString() {
425            StringBuilder sb = new StringBuilder( super.toString() );
426            sb.append( "\n - underlying-geographic-CRS: " ).append( geographicCRS );
427            sb.append( "\n - units: " ).append( units );
428            sb.append( "\n - projection-longitude: " ).append( projectionLongitude );
429            sb.append( "\n - projection-latitude: " ).append( projectionLatitude );
430            sb.append( "\n - is-spherical: " ).append( isSpherical() );
431            sb.append( "\n - is-conformal: " ).append( isConformal() );
432            sb.append( "\n - natural-origin: " ).append( getNaturalOrigin() );
433            sb.append( "\n - false-easting: " ).append( getFalseEasting() );
434            sb.append( "\n - false-northing: " ).append( getFalseNorthing() );
435            sb.append( "\n - scale: " ).append( getScale() );
436            return sb.toString();
437    
438        }
439        
440        /**
441         * Implementation as proposed by Joshua Block in Effective Java (Addison-Wesley 2001), which supplies an even
442         * distribution and is relatively fast. It is created from field <b>f</b> as follows:
443         * <ul>
444         * <li>boolean -- code = (f ? 0 : 1)</li>
445         * <li>byte, char, short, int -- code = (int)f </li>
446         * <li>long -- code = (int)(f ^ (f &gt;&gt;&gt;32))</li>
447         * <li>float -- code = Float.floatToIntBits(f);</li>
448         * <li>double -- long l = Double.doubleToLongBits(f); code = (int)(l ^ (l &gt;&gt;&gt; 32))</li>
449         * <li>all Objects, (where equals(&nbsp;) calls equals(&nbsp;) for this field) -- code = f.hashCode(&nbsp;)</li>
450         * <li>Array -- Apply above rules to each element</li>
451         * </ul>
452         * <p>
453         * Combining the hash code(s) computed above: result = 37 * result + code;
454         * </p>
455         * 
456         * @return (int) ( result >>> 32 ) ^ (int) result;
457         * 
458         * @see java.lang.Object#hashCode()
459         */
460        @Override
461        public int hashCode() {
462            // the 2nd millionth prime, :-)
463            long code = 32452843;
464            if( units != null ){
465                code = code * 37 + units.hashCode();
466            }
467    
468            long tmp = Double.doubleToLongBits( projectionLatitude );
469            code = code * 37 + (int)( tmp ^ (tmp >>> 32 ));
470            
471            tmp = Double.doubleToLongBits( projectionLongitude );
472            code = code * 37 + (int)( tmp ^ (tmp >>> 32 ));
473            
474            tmp = Double.doubleToLongBits( falseNorthing );
475            code = code * 37 + (int)( tmp ^ (tmp >>> 32 ));
476            
477            tmp = Double.doubleToLongBits( falseEasting );
478            code = code * 37 + (int)( tmp ^ (tmp >>> 32 ));
479            
480            tmp = Double.doubleToLongBits( scale );
481            code = code * 37 + (int)( tmp ^ (tmp >>> 32 ));
482            
483            code = code * 37 + (conformal? 0:1);
484            code = code * 37 + (equalArea? 0:1);
485            
486            return (int) ( code >>> 32 ) ^ (int) code;
487        }
488    }