/*
* Geotools 2 - OpenSource mapping toolkit
* (C) 2003, Geotools Project Managment Committee (PMC)
* (C) 2001, Institut de Recherche pour le Développement
*
* This library is free software; you can redistribute it and/or
* modify it under the terms of the GNU Lesser General Public
* License as published by the Free Software Foundation; either
* version 2.1 of the License, or (at your option) any later version.
*
* This library is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
* Lesser General Public License for more details.
*
* You should have received a copy of the GNU Lesser General Public
* License along with this library; if not, write to the Free Software
* Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
*
*
* This package contains documentation from OpenGIS specifications.
* OpenGIS consortium's work is fully acknowledged here.
*/
package org.geotools.ct;
// J2SE and JAI dependencies
import java.io.Serializable;
import javax.media.jai.ParameterList;
import javax.media.jai.util.Range;
import org.geotools.cs.Ellipsoid;
import org.geotools.cs.HorizontalDatum;
import org.geotools.cs.WGS84ConversionInfo;
import org.geotools.resources.i18n.VocabularyKeys;
import org.geotools.resources.i18n.ErrorKeys;
import org.geotools.resources.i18n.Errors;
/**
* The Molodensky transformation (EPSG code 9604) transforms three dimensional
* geographic points from one geographic coordinate reference system to another
* (a datum shift), using three shift parameters (delta X, delta Y, delta Z) and
* the difference between the semi-major axis and flattenings of the two ellipsoids.
*
*
* Unlike the Bursa-Wolf 3 parameter method (which acts on geocentric coordinates),
* this transformation can be performed directly on geographic coordinates.
*
*
* References:
true
for a 3D transformation, or
* false
for a 2D transformation.
*/
private final boolean source3D, target3D;
/**
* X,Y,Z shift in meters.
*/
private final double dx, dy, dz;
/**
* Semi-major (a) semi-minor (b/) radius in meters.
*/
private final double a, b;
/**
* Difference in the semi-major (da = target a - source a
) and semi-minor
* (db=target b - source b
) axes of the target and source ellipsoids.
*/
private final double da, db;
/**
* Difference between the flattenings (df = target f - source f
)
* of the target and source ellipsoids.
*/
private final double df;
/**
* Ratio of the Semi-major (a) semi-minor (b/) axis
* values (a_b = a/b
and b_a = b/a
).
*/
private final double b_a, a_b;
/**
* Some more constants (daa = da*a
and da_a = da/a
).
*/
private final double daa, da_a;
/**
* The square of excentricity of the ellipsoid: e˛ = (a˛-b˛)/a˛ where
* a is the semi-major axis length and
* b is the semi-minor axis length.
*/
private final double e2;
/**
* Construct a MolodenskiTransform from the specified datums.
*
* @param source source horizontal datum you are transforming from.
* @param target target horizontal datum you are transforming to.
* @param source3D true
if the source geographic CRS has a Z-axis (3 dimentional)
* @param target3D true
if the target geographic CRS has a Z-axis (3 dimentional)
*/
protected MolodenskiTransform(final HorizontalDatum source,
final HorizontalDatum target,
final boolean source3D, final boolean target3D)
{
double f;
final WGS84ConversionInfo srcInfo = source.getWGS84Parameters();
final WGS84ConversionInfo tgtInfo = target.getWGS84Parameters();
final Ellipsoid srcEllipsoid = source.getEllipsoid();
final Ellipsoid tgtEllipsoid = target.getEllipsoid();
dx = srcInfo.dx - tgtInfo.dx;
dy = srcInfo.dy - tgtInfo.dy;
dz = srcInfo.dz - tgtInfo.dz;
a = srcEllipsoid.getSemiMajorAxis();
b = srcEllipsoid.getSemiMinorAxis();
da = tgtEllipsoid.getSemiMajorAxis() - a;
db = tgtEllipsoid.getSemiMinorAxis() - b;
a_b = a/b;
b_a = b/a;
daa = da*a;
da_a = da/a;
f = 1 / srcEllipsoid.getInverseFlattening();
df = 1/tgtEllipsoid.getInverseFlattening() - f;
e2 = 1 - (b*b)/(a*a);
this.source3D = source3D;
this.target3D = target3D;
}
/**
* Construct a MolodenskiTransform from the specified parameters.
*
* @param parameters The parameter values in standard units.
*/
protected MolodenskiTransform(final ParameterList parameters) {
final int dim = parameters.getIntParameter("dim");
switch (dim) {
case 2: source3D=target3D=false; break;
case 3: source3D=target3D=true; break;
default: throw new IllegalArgumentException(Errors.format(
ErrorKeys.ILLEGAL_ARGUMENT_$2, "dim", new Integer(dim)));
}
final double ta, tb, f;
dx = parameters.getDoubleParameter("dx");
dy = parameters.getDoubleParameter("dy");
dz = parameters.getDoubleParameter("dz");
a = parameters.getDoubleParameter("src_semi_major");
b = parameters.getDoubleParameter("src_semi_minor");
ta = parameters.getDoubleParameter("tgt_semi_major");
tb = parameters.getDoubleParameter("tgt_semi_minor");
da = ta - a;
db = tb - b;
a_b = a/b;
b_a = b/a;
daa = da*a;
da_a = da/a;
f = (a-b)/a;
df = (ta-tb)/ta - f;
e2 = 1 - (b*b)/(a*a);
}
/**
* Transforms a list of coordinate point ordinal values.
* This method is provided for efficiently transforming many points.
* The supplied array of ordinal values will contain packed ordinal
* values. For example, if the source dimension is 3, then the ordinals
* will be packed in this order:
*
* (x0,y0,z0,
* x1,y1,z1 ...).
*
* @param srcPts the array containing the source point coordinates.
* @param srcOff the offset to the first point to be transformed
* in the source array.
* @param dstPts the array into which the transformed point
* coordinates are returned. May be the same
* than srcPts
.
* @param dstOff the offset to the location of the first
* transformed point that is stored in the
* destination array.
* @param numPts the number of point objects to be transformed.
*/
public void transform(final double[] srcPts, int srcOff,
final double[] dstPts, int dstOff, int numPts)
{
int step = 0;
if (srcPts==dstPts && srcOff