001 //$HeadURL: svn+ssh://rbezema@svn.wald.intevation.org/deegree/base/tags/2.1/src/org/deegree/model/csct/cs/Ellipsoid.java $
002 /*---------------- FILE HEADER ------------------------------------------
003
004 This file is part of deegree.
005 Copyright (C) 2001 by:
006 EXSE, Department of Geography, University of Bonn
007 http://www.giub.uni-bonn.de/exse/
008 lat/lon GmbH
009 http://www.lat-lon.de
010
011 It has been implemented within SEAGIS - An OpenSource implementation of OpenGIS specification
012 (C) 2001, Institut de Recherche pour le D�veloppement (http://sourceforge.net/projects/seagis/)
013 SEAGIS Contacts: Surveillance de l'Environnement Assist�e par Satellite
014 Institut de Recherche pour le D�veloppement / US-Espace
015 mailto:seasnet@teledetection.fr
016
017
018 This library is free software; you can redistribute it and/or
019 modify it under the terms of the GNU Lesser General Public
020 License as published by the Free Software Foundation; either
021 version 2.1 of the License, or (at your option) any later version.
022
023 This library is distributed in the hope that it will be useful,
024 but WITHOUT ANY WARRANTY; without even the implied warranty of
025 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
026 Lesser General Public License for more details.
027
028 You should have received a copy of the GNU Lesser General Public
029 License along with this library; if not, write to the Free Software
030 Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
031
032 Contact:
033
034 Andreas Poth
035 lat/lon GmbH
036 Aennchenstr. 19
037 53115 Bonn
038 Germany
039 E-Mail: poth@lat-lon.de
040
041 Klaus Greve
042 Department of Geography
043 University of Bonn
044 Meckenheimer Allee 166
045 53115 Bonn
046 Germany
047 E-Mail: klaus.greve@uni-bonn.de
048
049
050 ---------------------------------------------------------------------------*/
051 package org.deegree.model.csct.cs;
052
053 // OpenGIS dependencies
054 import java.awt.geom.Point2D;
055 import java.util.Map;
056
057 import org.deegree.model.csct.resources.Utilities;
058 import org.deegree.model.csct.resources.XMath;
059 import org.deegree.model.csct.resources.css.ResourceKeys;
060 import org.deegree.model.csct.resources.css.Resources;
061 import org.deegree.model.csct.units.Unit;
062
063 /**
064 * The figure formed by the rotation of an ellipse about an axis. In this context, the axis of
065 * rotation is always the minor axis. It is named geodetic ellipsoid if the parameters are derived
066 * by the measurement of the shape and the size of the Earth to approximate the geoid as close as
067 * possible.
068 *
069 * @version 1.00
070 * @author OpenGIS (www.opengis.org)
071 * @author Martin Desruisseaux
072 *
073 * @author last edited by: $Author: bezema $
074 *
075 * @version $Revision: 6259 $, $Date: 2007-03-20 10:15:15 +0100 (Di, 20 Mär 2007) $
076 *
077 * @see "org.opengis.cs.CS_Ellipsoid"
078 */
079 public class Ellipsoid extends Info {
080 /**
081 * Serial number for interoperability with different versions.
082 */
083 private static final long serialVersionUID = -1047804526105439230L;
084
085 /**
086 * WGS 1984 ellipsoid. This ellipsoid is used in GPS system and is the default for most
087 * <code>org.deegree.model</code> packages.
088 */
089 public static final Ellipsoid WGS84 = (Ellipsoid) pool.intern( createFlattenedSphere(
090 "WGS84",
091 6378137.0,
092 298.257223563,
093 Unit.METRE ) );
094
095 /**
096 * The equatorial radius.
097 */
098 private final double semiMajorAxis;
099
100 /**
101 * The polar radius.
102 */
103 private final double semiMinorAxis;
104
105 /**
106 * The inverse of the flattening value, or {@link Double#POSITIVE_INFINITY} if the ellipsoid is
107 * a sphere.
108 *
109 */
110 private final double inverseFlattening;
111
112 /**
113 * Is the Inverse Flattening definitive for this ellipsoid?
114 *
115 */
116 private final boolean ivfDefinitive;
117
118 /**
119 * The units of the semi-major and semi-minor axis values.
120 */
121 private final Unit unit;
122
123 /**
124 * Construct a new sphere using the specified radius.
125 *
126 * @param name
127 * Name of this sphere.
128 * @param radius
129 * The equatorial and polar radius.
130 * @param unit
131 * The units of the semi-major and semi-minor axis values.
132 */
133 public Ellipsoid( final String name, final double radius, final Unit unit ) {
134 this( name, check( "radius", radius ), radius, Double.POSITIVE_INFINITY, false, unit );
135 }
136
137 /**
138 * Construct a new ellipsoid using the specified axis length.
139 *
140 * @param name
141 * Name of this ellipsoid.
142 * @param semiMajorAxis
143 * The equatorial radius.
144 * @param semiMinorAxis
145 * The polar radius.
146 * @param unit
147 * The units of the semi-major and semi-minor axis values.
148 *
149 */
150 public Ellipsoid( final String name, final double semiMajorAxis, final double semiMinorAxis,
151 final Unit unit ) {
152 this( name, semiMajorAxis, semiMinorAxis,
153 semiMajorAxis / ( semiMajorAxis - semiMinorAxis ), false, unit );
154 }
155
156 /**
157 * Construct a new ellipsoid using the specified axis length.
158 *
159 * @param name
160 * Name of this ellipsoid.
161 * @param semiMajorAxis
162 * The equatorial radius.
163 * @param semiMinorAxis
164 * The polar radius.
165 * @param inverseFlattening
166 * The inverse of the flattening value.
167 * @param ivfDefinitive
168 * Is the Inverse Flattening definitive for this ellipsoid?
169 * @param unit
170 * The units of the semi-major and semi-minor axis values.
171 */
172 private Ellipsoid( final String name, final double semiMajorAxis, final double semiMinorAxis,
173 final double inverseFlattening, final boolean ivfDefinitive, final Unit unit ) {
174 super( name );
175 this.unit = unit;
176 this.semiMajorAxis = check( "semiMajorAxis", semiMajorAxis );
177 this.semiMinorAxis = check( "semiMinorAxis", semiMinorAxis );
178 this.inverseFlattening = check( "inverseFlattening", inverseFlattening );
179 this.ivfDefinitive = ivfDefinitive;
180 ensureNonNull( "unit", unit );
181 ensureLinearUnit( unit );
182 }
183
184 /**
185 * Construct a new ellipsoid using the specified axis length.
186 *
187 * @param properties
188 * The set of properties (see {@link Info}).
189 * @param semiMajorAxis
190 * The equatorial radius.
191 * @param semiMinorAxis
192 * The polar radius.
193 * @param inverseFlattening
194 * The inverse of the flattening value.
195 * @param ivfDefinitive
196 * Is the Inverse Flattening definitive for this ellipsoid?
197 * @param unit
198 * The units of the semi-major and semi-minor axis values.
199 */
200 Ellipsoid( final Map properties, final double semiMajorAxis, final double semiMinorAxis,
201 final double inverseFlattening, final boolean ivfDefinitive, final Unit unit ) {
202 super( properties );
203 this.unit = unit;
204 this.semiMajorAxis = semiMajorAxis;
205 this.semiMinorAxis = semiMinorAxis;
206 this.inverseFlattening = inverseFlattening;
207 this.ivfDefinitive = ivfDefinitive;
208 // Accept null values.
209 }
210
211 /**
212 * Construct a new ellipsoid using the specified axis length and inverse flattening value.
213 *
214 * @param name
215 * Name of this ellipsoid.
216 * @param semiMajorAxis
217 * The equatorial radius.
218 * @param inverseFlattening
219 * The inverse flattening value.
220 * @param unit
221 * The units of the semi-major and semi-minor axis values.
222 * @return
223 *
224 */
225 public static Ellipsoid createFlattenedSphere( final String name, final double semiMajorAxis,
226 final double inverseFlattening, final Unit unit ) {
227 return new Ellipsoid( name, semiMajorAxis, semiMajorAxis * ( 1 - 1 / inverseFlattening ),
228 inverseFlattening, true, unit );
229 }
230
231 /**
232 * Check the argument validity. Argument <code>value</code> should be greater than zero.
233 *
234 * @param name
235 * Argument name.
236 * @param value
237 * Argument value.
238 * @return <code>value</code>.
239 * @throws IllegalArgumentException
240 * if <code>value</code> is not greater than 0.
241 */
242 private static double check( final String name, final double value )
243 throws IllegalArgumentException {
244 if ( value > 0 )
245 return value;
246 throw new IllegalArgumentException(
247 Resources.format(
248 ResourceKeys.ERROR_ILLEGAL_ARGUMENT_$2,
249 name, new Double( value ) ) );
250 }
251
252 /**
253 * Gets the equatorial radius. The returned length is expressed in this object's axis units.
254 *
255 * @return the equatorial radius. The returned length is expressed in this object's axis units.
256 *
257 * @see "org.opengis.cs.CS_Ellipsoid#getSemiMajorAxis()"
258 */
259 public double getSemiMajorAxis() {
260 return semiMajorAxis;
261 }
262
263 /**
264 * Gets the polar radius. The returned length is expressed in this object's axis units.
265 *
266 * @return the polar radius. The returned length is expressed in this object's axis units.
267 *
268 * @see "org.opengis.cs.CS_Ellipsoid#getSemiMinorAxis()"
269 */
270 public double getSemiMinorAxis() {
271 return semiMinorAxis;
272 }
273
274 /**
275 * The ratio of the distance between the center and a focus of the ellipse to the length of its
276 * semimajor axis. The eccentricity can alternately be computed from the equation:
277 * <code>e=sqrt(2f-f�)</code>.
278 *
279 * @return
280 */
281 public double getEccentricity() {
282 final double f = 1 - getSemiMinorAxis() / getSemiMajorAxis();
283 return Math.sqrt( 2 * f - f * f );
284 }
285
286 /**
287 * Returns the value of the inverse of the flattening constant. Flattening is a value used to
288 * indicate how closely an ellipsoid approaches a spherical shape. The inverse flattening is
289 * related to the equatorial/polar radius (<var>r<sub>e</sub></var> and <var>r<sub>p</sub></var>
290 * respectively) by the formula <code>ivf=r<sub>e</sub>/(r<sub>e</sub>-r<sub>p</sub>)</code>.
291 * For perfect spheres, this method returns {@link Double#POSITIVE_INFINITY} (which is the
292 * correct value).
293 *
294 * @return the value of the inverse of the flattening constant.
295 *
296 * @see "org.opengis.cs.CS_Ellipsoid#getInverseFlattening()"
297 */
298 public double getInverseFlattening() {
299 return inverseFlattening;
300 }
301
302 /**
303 * Is the Inverse Flattening definitive for this ellipsoid? Some ellipsoids use the IVF as the
304 * defining value, and calculate the polar radius whenever asked. Other ellipsoids use the polar
305 * radius to calculate the IVF whenever asked. This distinction can be important to avoid
306 * floating-point rounding errors.
307 *
308 * @return
309 */
310 public boolean isIvfDefinitive() {
311 return ivfDefinitive;
312 }
313
314 /**
315 * Returns an <em>estimation</em> of orthodromic distance between two geographic coordinates.
316 * The orthodromic distance is the shortest distance between two points on a sphere's surface.
317 * The orthodromic path is always on a great circle. An other possible distance measurement is
318 * the loxodromic distance, which is a longer distance on a path with a constant direction on
319 * the compas.
320 *
321 * @param P1
322 * Longitude and latitude of first point (in degrees).
323 * @param P2
324 * Longitude and latitude of second point (in degrees).
325 * @return The orthodromic distance (in the units of this ellipsoid).
326 */
327 public double orthodromicDistance( final Point2D P1, final Point2D P2 ) {
328 return orthodromicDistance( P1.getX(), P1.getY(), P2.getX(), P2.getY() );
329 }
330
331 /**
332 * Returns an <em>estimation</em> of orthodromic distance between two geographic coordinates.
333 * The orthodromic distance is the shortest distance between two points on a sphere's surface.
334 * The orthodromic path is always on a great circle. An other possible distance measurement is
335 * the loxodromic distance, which is a longer distance on a path with a constant direction on
336 * the compas.
337 *
338 * @param x1
339 * Longitude of first point (in degrees).
340 * @param y1
341 * Latitude of first point (in degrees).
342 * @param x2
343 * Longitude of second point (in degrees).
344 * @param y2
345 * Latitude of second point (in degrees).
346 * @return The orthodromic distance (in the units of this ellipsoid).
347 */
348 public double orthodromicDistance( double x1, double y1, double x2, double y2 ) {
349 /*
350 * Le calcul de la distance orthodromique sur une surface ellipso�dale est complexe, sujet �
351 * des erreurs d'arrondissements et sans solution � proximit� des p�les. Nous utilisont
352 * plut�t un calcul bas� sur une forme sph�rique de la terre. Un programme en Fortran
353 * calculant les distances orthodromiques sur une surface ellipso�dale peut �tre t�l�charg� �
354 * partir du site de NOAA:
355 *
356 * ftp://ftp.ngs.noaa.gov/pub/pcsoft/for_inv.3d/source/
357 */
358 y1 = Math.toRadians( y1 );
359 y2 = Math.toRadians( y2 );
360 final double y = 0.5 * ( y1 + y2 );
361 final double dx = Math.toRadians( Math.abs( x2 - x1 ) % 360 );
362 double rho = Math.sin( y1 ) * Math.sin( y2 ) + Math.cos( y1 ) * Math.cos( y2 )
363 * Math.cos( dx );
364 if ( rho > +1 )
365 rho = +1; // Catch rounding error.
366 if ( rho < -1 )
367 rho = -1; // Catch rounging error.
368 return Math.acos( rho )
369 / XMath.hypot( Math.sin( y ) / getSemiMajorAxis(), Math.cos( y )
370 / getSemiMinorAxis() );
371 // 'hypot' calcule l'inverse du rayon **apparent** de la terre � la latitude 'y'.
372 }
373
374 /**
375 * Returns the units of the semi-major and semi-minor axis values.
376 *
377 * @return he units of the semi-major and semi-minor axis values.
378 *
379 * @see "org.opengis.cs.CS_Ellipsoid#getAxisUnit()"
380 */
381 public Unit getAxisUnit() {
382 return unit;
383 }
384
385 /**
386 * Compares the specified object with this ellipsoid for equality.
387 *
388 * @param object
389 * @return
390 */
391 public boolean equals( final Object object ) {
392 if ( super.equals( object ) ) {
393 final Ellipsoid that = (Ellipsoid) object;
394 return this.ivfDefinitive == that.ivfDefinitive
395 && Double.doubleToLongBits( this.semiMajorAxis ) == Double.doubleToLongBits( that.semiMajorAxis )
396 && Double.doubleToLongBits( this.semiMinorAxis ) == Double.doubleToLongBits( that.semiMinorAxis )
397 && Double.doubleToLongBits( this.inverseFlattening ) == Double.doubleToLongBits( that.inverseFlattening )
398 && Utilities.equals( this.unit, that.unit );
399 }
400 return false;
401 }
402
403 /**
404 * Returns a hash value for this ellipsoid.
405 *
406 * @return a hash value for this ellipsoid.
407 */
408 public int hashCode() {
409 final long longCode = Double.doubleToLongBits( getSemiMajorAxis() );
410 return ( ( (int) ( longCode >>> 32 ) ) ^ (int) longCode ) + 37 * super.hashCode();
411 }
412
413 /**
414 * Fill the part inside "[...]". Used for formatting Well Know Text (WKT).
415 *
416 * @param buffer
417 * @return
418 */
419 String addString( final StringBuffer buffer ) {
420 buffer.append( ", " );
421 buffer.append( semiMajorAxis );
422 buffer.append( ", " );
423 buffer.append( Double.isInfinite( inverseFlattening ) ? 0 : inverseFlattening );
424 return "SPHEROID";
425 }
426
427 }