001 //$HeadURL: svn+ssh://jwilden@svn.wald.intevation.org/deegree/base/branches/2.5_testing/src/org/deegree/crs/transformations/helmert/Helmert.java $
002 /*----------------------------------------------------------------------------
003 This file is part of deegree, http://deegree.org/
004 Copyright (C) 2001-2009 by:
005 Department of Geography, University of Bonn
006 and
007 lat/lon GmbH
008
009 This library is free software; you can redistribute it and/or modify it under
010 the terms of the GNU Lesser General Public License as published by the Free
011 Software Foundation; either version 2.1 of the License, or (at your option)
012 any later version.
013 This library is distributed in the hope that it will be useful, but WITHOUT
014 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
015 FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
016 details.
017 You should have received a copy of the GNU Lesser General Public License
018 along with this library; if not, write to the Free Software Foundation, Inc.,
019 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
020
021 Contact information:
022
023 lat/lon GmbH
024 Aennchenstr. 19, 53177 Bonn
025 Germany
026 http://lat-lon.de/
027
028 Department of Geography, University of Bonn
029 Prof. Dr. Klaus Greve
030 Postfach 1147, 53001 Bonn
031 Germany
032 http://www.geographie.uni-bonn.de/deegree/
033
034 e-mail: info@deegree.org
035 ----------------------------------------------------------------------------*/
036 package org.deegree.crs.transformations.helmert;
037
038 import static org.deegree.crs.projections.ProjectionUtils.EPS11;
039
040 import java.util.List;
041
042 import javax.vecmath.Matrix4d;
043 import javax.vecmath.Point3d;
044
045 import org.deegree.crs.Identifiable;
046 import org.deegree.crs.coordinatesystems.CoordinateSystem;
047 import org.deegree.crs.exceptions.TransformationException;
048 import org.deegree.crs.transformations.Transformation;
049
050 /**
051 * Parameters for a geographic transformation into another datum. The Bursa Wolf parameters should be applied to
052 * geocentric coordinates, where the X axis points towards the Greenwich Prime Meridian, the Y axis points East, and the
053 * Z axis points North.
054 *
055 *
056 * @author <a href="mailto:bezema@lat-lon.de">Rutger Bezema</a>
057 *
058 * @author last edited by: $Author: mschneider $
059 *
060 * @version $Revision: 18195 $, $Date: 2009-06-18 17:55:39 +0200 (Do, 18 Jun 2009) $
061 *
062 */
063 public class Helmert extends Transformation {
064
065 private static final long serialVersionUID = -1490341500541671907L;
066
067 /** Bursa Wolf shift in meters. */
068 public double dx;
069
070 /** Bursa Wolf shift in meters. */
071 public double dy;
072
073 /** Bursa Wolf shift in meters. */
074 public double dz;
075
076 /** Bursa Wolf rotation in arc seconds, which is 1/3600 of a degree. */
077 public double ex;
078
079 /** Bursa Wolf rotation in arc seconds. */
080 public double ey;
081
082 /** Bursa Wolf rotation in arc seconds. */
083 public double ez;
084
085 /** Bursa Wolf scaling in parts per million. */
086 public double ppm;
087
088 private Matrix4d transformMatrix = null;
089
090 private Matrix4d inverseMatrix = null;
091
092 private boolean rotationInRadians = false;
093
094 /**
095 * @param dx
096 * Bursa Wolf shift in meters.
097 * @param dy
098 * Bursa Wolf shift in meters.
099 * @param dz
100 * Bursa Wolf shift in meters.
101 * @param ex
102 * Bursa Wolf rotation in arc seconds or in radians (by setting the flag).
103 * @param ey
104 * Bursa Wolf rotation in arc seconds or in radians (by setting the flag).
105 * @param ez
106 * Bursa Wolf rotation in arc seconds or in radians (by setting the flag).
107 * @param ppm
108 * Bursa Wolf scaling in parts per million.
109 * @param sourceCRS
110 * of this helmert transformation
111 * @param targetCRS
112 * of this helmert transformation
113 * @param identifiable
114 * object containing all relevant id.
115 * @param inRadians
116 * true if the rotation parameters are in radians
117 */
118 public Helmert( double dx, double dy, double dz, double ex, double ey, double ez, double ppm,
119 CoordinateSystem sourceCRS, CoordinateSystem targetCRS, Identifiable identifiable, boolean inRadians ) {
120 super( sourceCRS, targetCRS, identifiable );
121 this.dx = dx;
122 this.dy = dy;
123 this.dz = dz;
124 this.ex = ex;
125 this.ey = ey;
126 this.ez = ez;
127 this.ppm = ppm;
128 rotationInRadians = inRadians;
129 }
130
131 /**
132 * @param dx
133 * Bursa Wolf shift in meters.
134 * @param dy
135 * Bursa Wolf shift in meters.
136 * @param dz
137 * Bursa Wolf shift in meters.
138 * @param ex
139 * Bursa Wolf rotation in arc seconds.
140 * @param ey
141 * Bursa Wolf rotation in arc seconds.
142 * @param ez
143 * Bursa Wolf rotation in arc seconds.
144 * @param ppm
145 * Bursa Wolf scaling in parts per million.
146 * @param sourceCRS
147 * of this helmert transformation
148 * @param targetCRS
149 * of this helmert transformation
150 * @param identifiable
151 * object containing all relevant id.
152 */
153 public Helmert( double dx, double dy, double dz, double ex, double ey, double ez, double ppm,
154 CoordinateSystem sourceCRS, CoordinateSystem targetCRS, Identifiable identifiable ) {
155 this( dx, dy, dz, ex, ey, ez, ppm, sourceCRS, targetCRS, identifiable, false );
156 }
157
158 /**
159 * Construct a conversion info with all parameters set to 0;
160 *
161 * @param sourceCRS
162 * of this helmert transformation
163 * @param targetCRS
164 * of this helmert transformation
165 *
166 * @param identifiers
167 * @param names
168 * @param versions
169 * @param descriptions
170 * @param areasOfUse
171 */
172 public Helmert( CoordinateSystem sourceCRS, CoordinateSystem targetCRS, String[] identifiers, String[] names,
173 String[] versions, String[] descriptions, String[] areasOfUse ) {
174 this( 0, 0, 0, 0, 0, 0, 0, sourceCRS, targetCRS, new Identifiable( identifiers, names, versions, descriptions,
175 areasOfUse ) );
176 }
177
178 /**
179 * Construct a conversion info with all parameters set to 0;
180 *
181 * @param sourceCRS
182 * of this helmert transformation
183 * @param targetCRS
184 * of this helmert transformation
185 * @param identifier
186 */
187 public Helmert( CoordinateSystem sourceCRS, CoordinateSystem targetCRS, String identifier ) {
188 this( sourceCRS, targetCRS, new String[] { identifier } );
189 }
190
191 /**
192 * Construct a conversion info with all parameters set to 0;
193 *
194 * @param sourceCRS
195 * of this helmert transformation
196 * @param targetCRS
197 * of this helmert transformation
198 * @param identifiers
199 */
200 public Helmert( CoordinateSystem sourceCRS, CoordinateSystem targetCRS, String[] identifiers ) {
201 this( sourceCRS, targetCRS, identifiers, null, null, null, null );
202 }
203
204 /**
205 * @param dx
206 * Bursa Wolf shift in meters.
207 * @param dy
208 * Bursa Wolf shift in meters.
209 * @param dz
210 * Bursa Wolf shift in meters.
211 * @param ex
212 * Bursa Wolf rotation in arc seconds.
213 * @param ey
214 * Bursa Wolf rotation in arc seconds.
215 * @param ez
216 * Bursa Wolf rotation in arc seconds.
217 * @param ppm
218 * Bursa Wolf scaling in parts per million.
219 * @param sourceCRS
220 * of this helmert transformation
221 * @param targetCRS
222 * of this helmert transformation
223 * @param identifiers
224 * @param names
225 * @param versions
226 * @param descriptions
227 * @param areaOfUses
228 */
229 public Helmert( double dx, double dy, double dz, double ex, double ey, double ez, double ppm,
230 CoordinateSystem sourceCRS, CoordinateSystem targetCRS, String[] identifiers, String[] names,
231 String[] versions, String[] descriptions, String[] areaOfUses ) {
232 this( dx, dy, dz, ex, ey, ez, ppm, sourceCRS, targetCRS, new Identifiable( identifiers, names, versions,
233 descriptions, areaOfUses ) );
234 }
235
236 /**
237 * @param dx
238 * Bursa Wolf shift in meters.
239 * @param dy
240 * Bursa Wolf shift in meters.
241 * @param dz
242 * Bursa Wolf shift in meters.
243 * @param ex
244 * Bursa Wolf rotation in arc seconds.
245 * @param ey
246 * Bursa Wolf rotation in arc seconds.
247 * @param ez
248 * Bursa Wolf rotation in arc seconds.
249 * @param ppm
250 * Bursa Wolf scaling in parts per million.
251 * @param sourceCRS
252 * of this helmert transformation
253 * @param targetCRS
254 * of this helmert transformation
255 * @param identifier
256 * @param name
257 * @param version
258 * @param description
259 * @param areaOfUse
260 */
261 public Helmert( double dx, double dy, double dz, double ex, double ey, double ez, double ppm,
262 CoordinateSystem sourceCRS, CoordinateSystem targetCRS, String identifier, String name,
263 String version, String description, String areaOfUse ) {
264 this( dx, dy, dz, ex, ey, ez, ppm, sourceCRS, targetCRS, new String[] { identifier }, new String[] { name },
265 new String[] { version }, new String[] { description }, new String[] { areaOfUse } );
266 }
267
268 /**
269 * @param dx
270 * Bursa Wolf shift in meters.
271 * @param dy
272 * Bursa Wolf shift in meters.
273 * @param dz
274 * Bursa Wolf shift in meters.
275 * @param ex
276 * Bursa Wolf rotation in arc seconds.
277 * @param ey
278 * Bursa Wolf rotation in arc seconds.
279 * @param ez
280 * Bursa Wolf rotation in arc seconds.
281 * @param ppm
282 * Bursa Wolf scaling in parts per million.
283 * @param sourceCRS
284 * of this helmert transformation
285 * @param targetCRS
286 * of this helmert transformation
287 * @param identifiers
288 */
289 public Helmert( double dx, double dy, double dz, double ex, double ey, double ez, double ppm,
290 CoordinateSystem sourceCRS, CoordinateSystem targetCRS, String[] identifiers ) {
291 this( dx, dy, dz, ex, ey, ez, ppm, sourceCRS, targetCRS, identifiers, null, null, null, null );
292 }
293
294 /**
295 * @param dx
296 * Bursa Wolf shift in meters.
297 * @param dy
298 * Bursa Wolf shift in meters.
299 * @param dz
300 * Bursa Wolf shift in meters.
301 * @param ex
302 * Bursa Wolf rotation in arc seconds.
303 * @param ey
304 * Bursa Wolf rotation in arc seconds.
305 * @param ez
306 * Bursa Wolf rotation in arc seconds.
307 * @param ppm
308 * Bursa Wolf scaling in parts per million.
309 * @param sourceCRS
310 * of this helmert transformation
311 * @param targetCRS
312 * of this helmert transformation
313 * @param identifier
314 */
315 public Helmert( double dx, double dy, double dz, double ex, double ey, double ez, double ppm,
316 CoordinateSystem sourceCRS, CoordinateSystem targetCRS, String identifier ) {
317 this( dx, dy, dz, ex, ey, ez, ppm, sourceCRS, targetCRS, new String[] { identifier } );
318 }
319
320 /**
321 * Returns an affine transformation also known as the "Helmert" transformation. The matrix representation of this
322 * transformation (also known as "Bursa Wolf" formula) is as follows:
323 *
324 * <blockquote>
325 *
326 * <pre>
327 * S = 1 + {@link #ppm}*1E-6
328 *
329 * [ X ] [ S -{@link #ez}*S +{@link #ey}*S {@link #dx} ] [ X ]
330 * [ Y ] = [ +{@link #ez}*S S -{@link #ex}*S {@link #dy} ] [ Y ]
331 * [ Z ] [ -{@link #ey}*S +{@link #ex}*S S {@link #dz} ] [ Z ]
332 * [ 1 ] [ 0 0 0 1 ] [ 1 ]
333 * </pre>
334 *
335 * </blockquote>
336 *
337 * This affine transform can be applied to transform <code>geocentric</code> coordinates from one datum into
338 * <code>geocentric</code> coordinates of an other datum. see <a
339 * href="http://www.posc.org/Epicentre.2_2/DataModel/ExamplesofUsage/eu_cs35.html#CS3523_helmert">
340 * http://www.posc.org/Epicentre.2_2/DataModel/ExamplesofUsage/eu_cs35.html</a> for more information.
341 *
342 * @return the affine "Helmert" transformation as a Matrix4d.
343 */
344 public Matrix4d getAsAffineTransform() {
345 if ( transformMatrix != null ) {
346 return new Matrix4d( transformMatrix );
347 }
348 // Note: (ex, ey, ez) is a rotation in arc seconds. We need to convert it into radians (the
349 // R factor in RS).
350 final double S = 1 + ( ppm * 1E-6 );
351 double arcToRad = .00000484813681109535; // Math.PI / ( 180. * 3600. )
352 if ( rotationInRadians ) {
353 arcToRad = 1;
354 }
355 final double RS = arcToRad * S;
356
357 transformMatrix = new Matrix4d( S, -ez * RS, +ey * RS, dx, +ez * RS, S, -ex * RS, dy, -ey * RS, +ex * RS, S,
358 dz, 0, 0, 0, 1. );
359 return new Matrix4d( transformMatrix );
360 }
361
362 /**
363 * @return true if any of the helmert parameters were set.
364 */
365 public boolean hasValues() {
366 return !( ex == 0 && ey == 0 && ez == 0 && dx == 0 && dy == 0 && dz == 0 && ppm == 0 );
367 }
368
369 @Override
370 public boolean equals( final Object other ) {
371 if ( other != null && other instanceof Helmert ) {
372 final Helmert that = (Helmert) other;
373 return ( Math.abs( this.dx - that.dx ) < EPS11 ) && ( Math.abs( this.dy - that.dy ) < EPS11 )
374 && ( Math.abs( this.dz - that.dz ) < EPS11 ) && ( Math.abs( this.ex - that.ex ) < EPS11 )
375 && ( Math.abs( this.ey - that.ey ) < EPS11 ) && ( Math.abs( this.ez - that.ez ) < EPS11 )
376 && ( Math.abs( this.ppm - that.ppm ) < EPS11 ) && super.equals( that );
377
378 }
379 return false;
380 }
381
382 /**
383 * Returns the Well Know Text (WKT) for this object. The WKT is part of OpenGIS's specification and looks like
384 * <code>TOWGS84[dx, dy, dz, ex, ey, ez, ppm]</code>.
385 *
386 * @return the Well Know Text (WKT) for this object.
387 */
388 @Override
389 public String toString() {
390 final StringBuffer buffer = new StringBuffer( super.getIdAndName() );
391 buffer.append( "\n[\"" );
392 buffer.append( dx );
393 buffer.append( ", " );
394 buffer.append( dy );
395 buffer.append( ", " );
396 buffer.append( dz );
397 buffer.append( ", " );
398 buffer.append( ex );
399 buffer.append( ", " );
400 buffer.append( ey );
401 buffer.append( ", " );
402 buffer.append( ez );
403 buffer.append( ", " );
404 buffer.append( ppm );
405 buffer.append( ']' );
406 return buffer.toString();
407 }
408
409 /**
410 * Implementation as proposed by Joshua Block in Effective Java (Addison-Wesley 2001), which supplies an even
411 * distribution and is relatively fast. It is created from field <b>f</b> as follows:
412 * <ul>
413 * <li>boolean -- code = (f ? 0 : 1)</li>
414 * <li>byte, char, short, int -- code = (int)f</li>
415 * <li>long -- code = (int)(f ^ (f >>>32))</li>
416 * <li>float -- code = Float.floatToIntBits(f);</li>
417 * <li>double -- long l = Double.doubleToLongBits(f); code = (int)(l ^ (l >>> 32))</li>
418 * <li>all Objects, (where equals( ) calls equals( ) for this field) -- code = f.hashCode( )</li>
419 * <li>Array -- Apply above rules to each element</li>
420 * </ul>
421 * <p>
422 * Combining the hash code(s) computed above: result = 37 * result + code;
423 * </p>
424 *
425 * @return (int) ( result >>> 32 ) ^ (int) result;
426 *
427 * @see java.lang.Object#hashCode()
428 */
429 @Override
430 public int hashCode() {
431 // the 2nd millionth prime, :-)
432 long code = 32452843;
433 long tmp = Double.doubleToLongBits( dx );
434 code = code * 37 + (int) ( tmp ^ ( tmp >>> 32 ) );
435
436 tmp = Double.doubleToLongBits( dy );
437 code = code * 37 + (int) ( tmp ^ ( tmp >>> 32 ) );
438
439 tmp = Double.doubleToLongBits( dz );
440 code = code * 37 + (int) ( tmp ^ ( tmp >>> 32 ) );
441
442 tmp = Double.doubleToLongBits( ex );
443 code = code * 37 + (int) ( tmp ^ ( tmp >>> 32 ) );
444
445 tmp = Double.doubleToLongBits( ey );
446 code = code * 37 + (int) ( tmp ^ ( tmp >>> 32 ) );
447
448 tmp = Double.doubleToLongBits( ez );
449 code = code * 37 + (int) ( tmp ^ ( tmp >>> 32 ) );
450
451 tmp = Double.doubleToLongBits( ppm );
452 code = code * 37 + (int) ( tmp ^ ( tmp >>> 32 ) );
453 return (int) ( code >>> 32 ) ^ (int) code;
454 }
455
456 @Override
457 public synchronized List<Point3d> doTransform( List<Point3d> srcPts )
458 throws TransformationException {
459
460 if ( srcPts == null || srcPts.size() == 0 ) {
461 return srcPts;
462 }
463
464 // lazy instantiation
465 if ( transformMatrix == null ) {
466 transformMatrix = getAsAffineTransform();
467 }
468 Matrix4d matrix = transformMatrix;
469
470 // create the inverse matrix
471 if ( isInverseTransform() ) {
472 if ( inverseMatrix == null ) {
473 transformMatrix.invert( inverseMatrix );
474 }
475 matrix = inverseMatrix;
476 }
477 for ( Point3d p : srcPts ) {
478 matrix.transform( p );
479 }
480
481 return srcPts;
482 }
483
484 @Override
485 public String getImplementationName() {
486 return "Helmert";
487 }
488
489 @Override
490 public boolean isIdentity() {
491 return hasValues();
492 }
493 }