blob: 840f6774071c35631b402bf2a75f8d733a9b1512 [file] [log] [blame]
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.sis.internal.referencing.provider;
import org.opengis.util.FactoryException;
import org.opengis.parameter.ParameterValueGroup;
import org.opengis.referencing.datum.Ellipsoid;
import org.opengis.referencing.operation.MathTransform;
import org.opengis.referencing.operation.MathTransformFactory;
import org.opengis.referencing.operation.NoninvertibleTransformException;
import org.opengis.referencing.operation.TransformException;
import org.apache.sis.internal.system.DefaultFactories;
import org.apache.sis.internal.referencing.Formulas;
import org.apache.sis.parameter.Parameters;
import org.apache.sis.referencing.CommonCRS;
import org.apache.sis.referencing.operation.matrix.Matrix4;
import org.apache.sis.referencing.operation.transform.CoordinateDomain;
import org.apache.sis.referencing.operation.transform.EllipsoidToCentricTransform;
import org.apache.sis.referencing.operation.transform.LinearTransform;
import org.apache.sis.referencing.operation.transform.MathTransformTestCase;
import org.apache.sis.test.DependsOnMethod;
import org.apache.sis.test.DependsOn;
import org.junit.Test;
import static java.lang.StrictMath.toRadians;
import static org.opengis.test.Assert.*;
/**
* Tests {@link GeocentricTranslation} and {@link GeocentricTranslation3D}.
*
* @author Martin Desruisseaux (Geomatys)
* @version 0.7
* @since 0.7
* @module
*/
@DependsOn({
AffineTest.class,
org.apache.sis.referencing.operation.transform.ProjectiveTransformTest.class,
org.apache.sis.referencing.operation.transform.ConcatenatedTransformTest.class,
org.apache.sis.referencing.operation.transform.EllipsoidToCentricTransformTest.class,
org.apache.sis.referencing.datum.BursaWolfParametersTest.class
})
public final strictfp class GeocentricTranslationTest extends MathTransformTestCase {
/**
* Geocentric translation parameters for transforming a point in the North Sea from WGS84 to ED50.
* They are the parameters to use for transformation of the point given by {@link #samplePoint(int)}.
*/
public static final double TX = 84.87, TY = 96.49, TZ = 116.95;
/**
* Returns the sample point for a step in the example given by the EPSG guidance note.
*
* <blockquote><b>Source:</b>
* §2.4.3.5 <cite>Three-parameter geocentric translations</cite> and
* §2.4.4.2 <cite>Abridged Molodensky transformation</cite> in
* IOGP Publication 373-7-2 – Geomatics Guidance Note number 7, part 2 – April 2015
* </blockquote>
*
* The example transforms a point in the North Sea from WGS84 to ED50.
* The transform can be decomposed in 4 steps:
* <pre>
* <b>Step 1:</b> Source point in WGS84: 53°48'33.820"N, 02°07'46.380"E, 73.00 metres.
* <b>Step 2:</b> Geocentric (X,Y,Z): 3771793.968, 140253.342, 5124304.349 metres.
* <b>Step 3:</b> (X,Y,Z) after shift: +84.87, +96.49, +116.95 metres
* <b>Step 4:</b> Target point in ED50: 53°48'36.565"N, 02'07"51.477"E, 28.02 metres.
* </pre>
*
* @param step the step as a value from 1 to 5 inclusive.
* @return the sample point at the given step.
*/
public static double[] samplePoint(final int step) {
switch (step) {
case 1: return new double[] {
2 + ( 7 + 46.380/60)/60, // λ: Longitude
53 + (48 + 33.820/60)/60, // φ: Latitude
73.00 // h: Height
};
case 2: return new double[] {
3771793.968, // X: Toward prime meridian
140253.342, // Y: Toward 90° east
5124304.349 // Z: Toward north pole
};
case 3: return new double[] {
3771878.84, // X: Toward prime meridian
140349.83, // Y: Toward 90° east
5124421.30 // Z: Toward north pole
};
case 4: return new double[] { // Result according geocentric translation
2 + ( 7 + 51.477/60)/60, // λ: Longitude
53 + (48 + 36.565/60)/60, // φ: Latitude
28.02 // h: Height
};
case 5: return new double[] { // Result according abridged Molodensky
2 + ( 7 + 51.477/60)/60, // λ: Longitude
53 + (48 + 36.563/60)/60, // φ: Latitude
28.091 // h: Height
};
default: throw new AssertionError(step);
}
}
/**
* Returns the sample point precision for the given step.
*
* @param step the step as a value from 1 to 4 inclusive.
* @return the sample point precision for the given step.
*/
public static double precision(final int step) {
switch (step) {
case 1:
case 4:
case 5: return 0.001 / 60 / 60 / 2; // Half the precision of (λ,φ) values given by EPSG
case 2: return 0.001 / 2; // Half the precision for (X,Y,Z) values given by EPSG
case 3: return 0.01 / 2; // Half the precision for (X,Y,Z) values given by EPSG
default: throw new AssertionError(step);
}
}
/**
* Creates a transformation for the given method.
*
* @param method the method to test.
*/
private void create(final GeocentricAffine method) throws FactoryException {
final ParameterValueGroup values = method.getParameters().createValue();
setTranslation(values);
if (method instanceof GeocentricAffineBetweenGeographic) {
setEllipsoids(values, CommonCRS.WGS84.ellipsoid(), CommonCRS.ED50.ellipsoid());
}
transform = method.createMathTransform(DefaultFactories.forBuildin(MathTransformFactory.class), values);
}
/**
* Creates a datum shift operation without using the provider. This method is equivalent to the
* <cite>"Geocentric translations (geog3D domain)"</cite> method (EPSG:1035), but the transform
* steps are created without the use of any {@link ParameterValueGroup}. This way to create the
* datum shift is provided for better separation of aspects being tested.
*
* @param factory the factory to use for creating the transforms.
* @param source the source ellipsoid.
* @param target the target ellipsoid.
* @param tX geocentric translation on the X axis.
* @param tY geocentric translation on the Y axis.
* @param tZ geocentric translation on the Z axis.
* @return transform performing the datum shift.
* @throws FactoryException if an error occurred while creating the transform.
* @throws NoninvertibleTransformException if an error occurred while creating the transform.
*
* @see #testGeographicDomain()
*/
public static MathTransform createDatumShiftForGeographic3D(
final MathTransformFactory factory,
final Ellipsoid source, final Ellipsoid target,
final double tX, final double tY, final double tZ)
throws FactoryException, NoninvertibleTransformException
{
final Matrix4 translation = new Matrix4();
translation.m03 = tX;
translation.m13 = tY;
translation.m23 = tZ;
final MathTransform step1 = EllipsoidToCentricTransform.createGeodeticConversion(factory, source, true);
final MathTransform step3 = EllipsoidToCentricTransform.createGeodeticConversion(factory, target, true).inverse();
final MathTransform step2 = factory.createAffineTransform(translation);
return factory.createConcatenatedTransform(step1,
factory.createConcatenatedTransform(step2, step3));
}
/**
* Creates a "Geographic 2D to 3D → Geocentric → Affine → Geographic → Geographic 3D to 2D" chain
* using EPSG or OGC standard operation methods and parameters. This is used for integration tests.
*
* @param factory the math transform factory to use for creating and concatenating the transform.
* @return the chain of transforms.
* @throws FactoryException if an error occurred while creating a transform.
*/
public static MathTransform createDatumShiftForGeographic2D(final MathTransformFactory factory) throws FactoryException {
final Parameters values = Parameters.castOrWrap(factory.getDefaultParameters("Geocentric translations (geog2D domain)"));
setTranslation(values);
setEllipsoids(values, CommonCRS.WGS84.ellipsoid(), CommonCRS.ED50.ellipsoid());
final MathTransform gt = new GeocentricTranslation().createMathTransform(factory, values);
assertFalse("isIdentity", gt.isIdentity());
assertEquals("sourceDimensions", 3, gt.getSourceDimensions());
assertEquals("targetDimensions", 3, gt.getTargetDimensions());
assertInstanceOf("Geocentric translation", LinearTransform.class, gt);
return Geographic3Dto2DTest.createDatumShiftForGeographic2D(factory, gt, values);
}
/**
* Sets the translation parameters in the given parameter value group.
*/
private static void setTranslation(final ParameterValueGroup values) {
values.parameter("X-axis translation").setValue(TX);
values.parameter("Y-axis translation").setValue(TY);
values.parameter("Z-axis translation").setValue(TZ);
}
/**
* Sets the source and target ellipsoid axes in the given parameter value group.
*/
static void setEllipsoids(final ParameterValueGroup values, final Ellipsoid source, final Ellipsoid target) {
values.parameter("src_semi_major").setValue(source.getSemiMajorAxis());
values.parameter("src_semi_minor").setValue(source.getSemiMinorAxis());
values.parameter("tgt_semi_major").setValue(target.getSemiMajorAxis());
values.parameter("tgt_semi_minor").setValue(target.getSemiMinorAxis());
}
/**
* Tests transformation of the sample point from WGS84 to ED50.
*
* @param sourceStep the {@link #samplePoint(int)} to use as the source coordinate.
* @param targetStep the {@link #samplePoint(int)} expected as a result of the transformation.
*/
private void datumShift(final int sourceStep, final int targetStep) throws TransformException {
tolerance = precision(targetStep);
tolerance = Formulas.LINEAR_TOLERANCE; // Other SIS branches use a stricter threshold.
verifyTransform(samplePoint(sourceStep), samplePoint(targetStep));
validate();
}
/**
* Tests <cite>"Geocentric translations (geocentric domain)"</cite> (EPSG:1031).
*
* @throws FactoryException if an error occurred while creating the transform.
* @throws TransformException if transformation of a point failed.
*/
@Test
public void testGeocentricDomain() throws FactoryException, TransformException {
create(new GeocentricTranslation());
assertTrue(transform instanceof LinearTransform);
derivativeDeltas = new double[] {100, 100, 100}; // In metres
datumShift(2, 3);
}
/**
* Tests <cite>"Geocentric translations (geog3D domain)"</cite> (EPSG:1035).
*
* @throws FactoryException if an error occurred while creating the transform.
* @throws TransformException if transformation of a point failed.
*/
@Test
@DependsOnMethod("testGeocentricDomain")
public void testGeographicDomain() throws FactoryException, TransformException {
create(new GeocentricTranslation3D());
assertFalse(transform instanceof LinearTransform);
final double delta = toRadians(100.0 / 60) / 1852; // Approximately 100 metres
derivativeDeltas = new double[] {delta, delta, 100}; // (Δλ, Δφ, Δh)
zTolerance = Formulas.LINEAR_TOLERANCE / 2; // Half the precision of h value given by EPSG
zDimension = new int[] {2}; // Dimension of h where to apply zTolerance
tolerance = Formulas.LINEAR_TOLERANCE; // Other SIS branches use a stricter threshold.
datumShift(1, 4);
}
/**
* Tests conversion of random points.
*
* @throws FactoryException if an error occurred while creating the transform.
* @throws TransformException if transformation of a point failed.
*/
@Test
@DependsOnMethod("testGeographicDomain")
public void testRandomPoints() throws FactoryException, TransformException {
testGeographicDomain(); // For creating the transform.
tolerance = Formulas.LINEAR_TOLERANCE;
// toleranceModifier = ToleranceModifier.GEOGRAPHIC;
verifyInDomain(CoordinateDomain.GEOGRAPHIC, 831342815);
}
/**
* Tests Well Known Text formatting of a two-dimensional transform.
* The main point of this test is to verify that Geographic 2D/3D conversions have been
* inserted before and after the transform, and that Geographic/Geocentric conversions
* have been replaced by Bursa-Wolf parameters for formatting purpose.
*
* @throws FactoryException if an error occurred while creating the transform.
*/
@Test
@DependsOnMethod("testWKT3D")
public void testWKT2D() throws FactoryException {
create(new GeocentricTranslation2D());
verifyWKT2D();
}
/**
* Verifies the WKT formatting of current transform.
*/
private void verifyWKT2D() {
assertWktEquals("CONCAT_MT[\n" +
" INVERSE_MT[PARAM_MT[“Geographic3D to 2D conversion”]],\n" +
" PARAM_MT[“Ellipsoid_To_Geocentric”,\n" +
" PARAMETER[“semi_major”, 6378137.0],\n" +
" PARAMETER[“semi_minor”, 6356752.314245179]],\n" +
" PARAM_MT[“Geocentric translations (geocentric domain)”,\n" +
" PARAMETER[“dx”, 84.87],\n" +
" PARAMETER[“dy”, 96.49],\n" +
" PARAMETER[“dz”, 116.95]],\n" +
" PARAM_MT[“Geocentric_To_Ellipsoid”,\n" +
" PARAMETER[“semi_major”, 6378388.0],\n" +
" PARAMETER[“semi_minor”, 6356911.9461279465]],\n" +
" PARAM_MT[“Geographic3D to 2D conversion”]]");
}
/**
* Tests Well Known Text formatting.
* The main point of this test is to verify that the affine transform between the two
* Geographic/Geocentric conversions have been replaced by Bursa-Wolf parameters for
* formatting purpose.
*
* @throws FactoryException if an error occurred while creating the transform.
*/
@Test
@DependsOnMethod("testGeographicDomain")
public void testWKT3D() throws FactoryException {
create(new GeocentricTranslation3D());
assertWktEquals("CONCAT_MT[\n" +
" PARAM_MT[“Ellipsoid_To_Geocentric”,\n" +
" PARAMETER[“semi_major”, 6378137.0],\n" +
" PARAMETER[“semi_minor”, 6356752.314245179]],\n" +
" PARAM_MT[“Geocentric translations (geocentric domain)”,\n" +
" PARAMETER[“dx”, 84.87],\n" +
" PARAMETER[“dy”, 96.49],\n" +
" PARAMETER[“dz”, 116.95]],\n" +
" PARAM_MT[“Geocentric_To_Ellipsoid”,\n" +
" PARAMETER[“semi_major”, 6378388.0],\n" +
" PARAMETER[“semi_minor”, 6356911.9461279465]]]");
/*
* The above WKT what was we show to users. But what we have in memory is different:
* affines before, after and between the two Geographic/Geocentric conversions.
*/
assertInternalWktEquals(
"Concat_MT[\n" +
" Param_MT[“Affine”,\n" +
" Parameter[“num_row”, 4],\n" +
" Parameter[“num_col”, 4],\n" +
" Parameter[“elt_0_0”, 0.017453292519943295],\n" +
" Parameter[“elt_1_1”, 0.017453292519943295],\n" +
" Parameter[“elt_2_2”, 1.567855942887398E-7]],\n" +
" Param_MT[“Ellipsoid (radians domain) to centric”,\n" +
" Parameter[“eccentricity”, 0.08181919084262157],\n" +
" Parameter[“target”, “CARTESIAN”],\n" +
" Parameter[“dim”, 3]],\n" +
" Param_MT[“Affine”,\n" +
" Parameter[“num_row”, 4],\n" +
" Parameter[“num_col”, 4],\n" +
" Parameter[“elt_0_0”, 0.9999606483644456],\n" +
" Parameter[“elt_0_3”, 1.3305869758942228E-5],\n" +
" Parameter[“elt_1_1”, 0.9999606483644456],\n" +
" Parameter[“elt_1_3”, 1.512764667185502E-5],\n" +
" Parameter[“elt_2_2”, 0.9999606483644456],\n" +
" Parameter[“elt_2_3”, 1.8335353697517302E-5]],\n" +
" Param_MT[“Centric to ellipsoid (radians domain)”,\n" +
" Parameter[“eccentricity”, 0.08199188997902956],\n" +
" Parameter[“target”, “CARTESIAN”],\n" +
" Parameter[“dim”, 3]],\n" +
" Param_MT[“Affine”,\n" +
" Parameter[“num_row”, 4],\n" +
" Parameter[“num_col”, 4],\n" +
" Parameter[“elt_0_0”, 57.29577951308232],\n" +
" Parameter[“elt_1_1”, 57.29577951308232],\n" +
" Parameter[“elt_2_2”, 6378388.0]]]");
}
/**
* Tests "Geographic 2D to 3D → Geocentric → Affine → Geographic → Geographic 3D to 2D" chain
* created from a {@link MathTransformFactory}. Because this test involves a lot of steps,
* this is more an integration test than a unit test: a failure here may not be easy to debug.
*
* @throws FactoryException if an error occurred while creating the transform.
*/
@Test
@DependsOnMethod("testWKT2D")
public void testIntegration() throws FactoryException {
transform = createDatumShiftForGeographic2D(DefaultFactories.forBuildin(MathTransformFactory.class));
verifyWKT2D();
}
}