use of org.opengis.referencing.crs.CoordinateReferenceSystem in project sis by apache.
the class CoordinateOperationTest method testGeocentricTransform.
/**
* Tests a "geographic to geocentric" conversion.
*
* @throws FactoryException if an error occurred while creating a test CRS.
* @throws TransformException if an error occurred while testing a coordinate conversion.
*/
@Test
public void testGeocentricTransform() throws FactoryException, TransformException {
final Random random = new Random(661597560);
/*
* Gets the math transform from WGS84 to a geocentric transform.
*/
final Ellipsoid ellipsoid = CommonCRS.WGS84.ellipsoid();
final CoordinateReferenceSystem sourceCRS = AbstractCRS.castOrCopy(CommonCRS.WGS84.geographic3D()).forConvention(AxesConvention.RIGHT_HANDED);
final CoordinateReferenceSystem targetCRS = CommonCRS.WGS84.geocentric();
final CoordinateOperation operation = opFactory.createOperation(sourceCRS, targetCRS);
transform = operation.getMathTransform();
final int dimension = transform.getSourceDimensions();
assertEquals("Source dimension", 3, dimension);
assertEquals("Target dimension", 3, transform.getTargetDimensions());
assertSame("Inverse transform", transform, transform.inverse().inverse());
validate();
/*
* Constructs an array of random points. The first 8 points
* are initialized to know values. Other points are left random.
*/
final double[] cartesianDistance = new double[4];
final double[] orthodromicDistance = new double[4];
// Must be divisible by 3.
final double[] array0 = new double[900];
for (int i = 0; i < array0.length; i++) {
final int range;
switch(i % 3) {
// Longitude
case 0:
range = 360;
break;
// Latitidue
case 1:
range = 180;
break;
// Altitude
case 2:
range = 10000;
break;
// Should not happen
default:
range = 0;
break;
}
array0[i] = random.nextDouble() * range - (range / 2);
}
// 24°N 35°E 8km
array0[0] = 35.0;
// 24°N 35°E 8km
array0[1] = 24.0;
// 24°N 35°E 8km
array0[2] = 8000;
// … about 80 km away
array0[3] = 34.8;
// … about 80 km away
array0[4] = 24.7;
// … about 80 km away
array0[5] = 5000;
cartesianDistance[0] = 80284.00;
// Not really exact.
orthodromicDistance[0] = 80302.99;
array0[6] = 0;
array0[7] = 0.0;
array0[8] = 0;
// Antipodes; distance should be 2*6378.137 km
array0[9] = 180;
// Antipodes; distance should be 2*6378.137 km
array0[10] = 0.0;
// Antipodes; distance should be 2*6378.137 km
array0[11] = 0;
cartesianDistance[1] = ellipsoid.getSemiMajorAxis() * 2;
orthodromicDistance[1] = ellipsoid.getSemiMajorAxis() * PI;
array0[12] = 0;
array0[13] = -90;
array0[14] = 0;
// Antipodes; distance should be 2*6356.752 km
array0[15] = 180;
// Antipodes; distance should be 2*6356.752 km
array0[16] = +90;
// Antipodes; distance should be 2*6356.752 km
array0[17] = 0;
cartesianDistance[2] = ellipsoid.getSemiMinorAxis() * 2;
orthodromicDistance[2] = 20003931.46;
array0[18] = 95;
array0[19] = -38;
array0[20] = 0;
// Antipodes
array0[21] = -85;
// Antipodes
array0[22] = +38;
// Antipodes
array0[23] = 0;
cartesianDistance[3] = 12740147.19;
orthodromicDistance[3] = 20003867.86;
/*
* Transforms all points, and then inverse transform them. The resulting
* array2 should be equal to array0 except for rounding errors. We tolerate
* maximal error of 0.1 second in longitude or latitude and 1 cm in height.
*/
final double[] array1 = new double[array0.length];
final double[] array2 = new double[array0.length];
transform.transform(array0, 0, array1, 0, array0.length / dimension);
transform.inverse().transform(array1, 0, array2, 0, array1.length / dimension);
for (int i = 0; i < array0.length; ) {
assertEquals("Longitude", array2[i], array0[i], 0.1 / 3600);
i++;
assertEquals("Latitude", array2[i], array0[i], 0.1 / 3600);
i++;
assertEquals("Height", array2[i], array0[i], 0.01);
i++;
}
/*
* Compares the distances between "special" points with expected distances.
* This tests the ellipsoid orthodromic distance computation as well.
* We require a precision of 10 centimetres.
*/
for (int i = 0; i < array0.length / 6; i++) {
final int base = i * 6;
final double cartesian = MathFunctions.magnitude(array1[base + 0] - array1[base + 3], array1[base + 1] - array1[base + 4], array1[base + 2] - array1[base + 5]);
if (i < cartesianDistance.length) {
assertEquals("Cartesian distance", cartesianDistance[i], cartesian, 0.1);
}
/*
* Compares with orthodromic distance. Distance is computed using an ellipsoid
* at the maximal altitude (i.e. the length of semi-major axis is increased to
* fit the maximal altitude).
*/
try {
final double altitude = max(array0[base + 2], array0[base + 5]);
final DefaultEllipsoid ellip = DefaultEllipsoid.createFlattenedSphere(Collections.singletonMap(Ellipsoid.NAME_KEY, "Temporary"), ellipsoid.getSemiMajorAxis() + altitude, ellipsoid.getInverseFlattening(), ellipsoid.getAxisUnit());
double orthodromic = ellip.orthodromicDistance(array0[base + 0], array0[base + 1], array0[base + 3], array0[base + 4]);
orthodromic = hypot(orthodromic, array0[base + 2] - array0[base + 5]);
if (i < orthodromicDistance.length) {
assertEquals("Orthodromic distance", orthodromicDistance[i], orthodromic, 0.1);
}
assertTrue("Distance consistency", cartesian <= orthodromic);
} catch (ArithmeticException exception) {
// Orthodromic distance computation didn't converge. Ignore...
}
}
}
use of org.opengis.referencing.crs.CoordinateReferenceSystem in project sis by apache.
the class Angle method valueOf.
/**
* Returns the angular value of the axis having the given direction.
* This helper method is used for subclass constructors expecting a {@link DirectPosition} argument.
*
* @param position the position from which to get an angular value.
* @param positive axis direction of positive values.
* @param negative axis direction of negative values.
* @return angular value in degrees.
* @throws IllegalArgumentException if the given coordinate it not associated to a CRS,
* or if no axis oriented toward the given directions is found, or if that axis
* does not use {@linkplain Units#isAngular angular units}.
*/
static double valueOf(final DirectPosition position, final AxisDirection positive, final AxisDirection negative) {
final CoordinateReferenceSystem crs = position.getCoordinateReferenceSystem();
if (crs == null) {
throw new IllegalArgumentException(Errors.format(Errors.Keys.UnspecifiedCRS));
}
final CoordinateSystem cs = crs.getCoordinateSystem();
final int dimension = cs.getDimension();
IncommensurableException cause = null;
for (int i = 0; i < dimension; i++) {
final CoordinateSystemAxis axis = cs.getAxis(i);
final AxisDirection dir = axis.getDirection();
final boolean isPositive = dir.equals(positive);
if (isPositive || dir.equals(negative)) {
double value = position.getOrdinate(i);
if (!isPositive)
value = -value;
final Unit<?> unit = axis.getUnit();
if (unit != Units.DEGREE)
try {
value = unit.getConverterToAny(Units.DEGREE).convert(value);
} catch (IncommensurableException e) {
cause = e;
break;
}
return value;
}
}
throw new IllegalArgumentException(Errors.format(Errors.Keys.IllegalCRSType_1, Classes.getLeafInterfaces(crs.getClass(), CoordinateReferenceSystem.class)[0]), cause);
}
use of org.opengis.referencing.crs.CoordinateReferenceSystem in project sis by apache.
the class ServicesForMetadata method addElements.
/**
* Initializes a horizontal, vertical and temporal extent with the values inferred from the given envelope.
*
* @param envelope the source envelope.
* @param target the target extent where to store envelope information.
* @throws TransformException if a coordinate transformation was required and failed.
*/
@Override
public void addElements(final Envelope envelope, final DefaultExtent target) throws TransformException {
final CoordinateReferenceSystem crs = envelope.getCoordinateReferenceSystem();
final SingleCRS horizontalCRS = CRS.getHorizontalComponent(crs);
final VerticalCRS verticalCRS = CRS.getVerticalComponent(crs, true);
final TemporalCRS temporalCRS = CRS.getTemporalComponent(crs);
if (horizontalCRS == null && verticalCRS == null && temporalCRS == null) {
throw new TransformException(dimensionNotFound(Resources.Keys.MissingSpatioTemporalDimension_1, crs));
}
if (horizontalCRS != null) {
final DefaultGeographicBoundingBox extent = new DefaultGeographicBoundingBox();
extent.setInclusion(Boolean.TRUE);
setBounds(envelope, extent);
target.getGeographicElements().add(extent);
}
if (verticalCRS != null) {
final DefaultVerticalExtent extent = new DefaultVerticalExtent();
setVerticalExtent(envelope, extent, crs, verticalCRS);
target.getVerticalElements().add(extent);
}
if (temporalCRS != null) {
final DefaultTemporalExtent extent = new DefaultTemporalExtent();
setTemporalExtent(envelope, extent, crs, temporalCRS);
target.getTemporalElements().add(extent);
}
}
use of org.opengis.referencing.crs.CoordinateReferenceSystem in project sis by apache.
the class ServicesForMetadata method setBounds.
/**
* Sets a geographic bounding box from the specified envelope.
* If the envelope contains a CRS which is not geographic, then the bounding box will be transformed
* to a geographic CRS (without datum shift if possible). Otherwise, the envelope is assumed already
* in a geographic CRS using (<var>longitude</var>, <var>latitude</var>) axis order.
*
* @param envelope the source envelope.
* @param target the target bounding box where to store envelope information.
* @throws TransformException if the given envelope can not be transformed.
*/
@Override
public void setBounds(Envelope envelope, final DefaultGeographicBoundingBox target) throws TransformException {
final CoordinateReferenceSystem crs = envelope.getCoordinateReferenceSystem();
GeographicCRS normalizedCRS = ReferencingUtilities.toNormalizedGeographicCRS(crs);
if (normalizedCRS == null) {
if (crs != null) {
normalizedCRS = CommonCRS.defaultGeographic();
} else if (envelope.getDimension() != 2) {
throw new TransformException(dimensionNotFound(Resources.Keys.MissingHorizontalDimension_1, crs));
}
}
setGeographicExtent(envelope, target, crs, normalizedCRS);
}
use of org.opengis.referencing.crs.CoordinateReferenceSystem in project sis by apache.
the class ServicesForMetadata method setBounds.
/**
* Sets the geographic, vertical and temporal extents with the values inferred from the given envelope.
* If the given {@code target} has more geographic or vertical extents than needed (0 or 1), then the
* extraneous extents are removed.
*
* @param envelope the source envelope.
* @param target the target spatiotemporal extent where to store envelope information.
* @throws TransformException if no temporal component can be extracted from the given envelope.
*/
@Override
public void setBounds(final Envelope envelope, final DefaultSpatialTemporalExtent target) throws TransformException {
final CoordinateReferenceSystem crs = envelope.getCoordinateReferenceSystem();
final SingleCRS horizontalCRS = CRS.getHorizontalComponent(crs);
final VerticalCRS verticalCRS = CRS.getVerticalComponent(crs, true);
final TemporalCRS temporalCRS = CRS.getTemporalComponent(crs);
if (horizontalCRS == null && verticalCRS == null && temporalCRS == null) {
throw new TransformException(dimensionNotFound(Resources.Keys.MissingSpatioTemporalDimension_1, crs));
}
/*
* Try to set the geographic bounding box first, because this operation may fail with a
* TransformException while the other operations (vertical and temporal) should not fail.
* So doing the geographic part first help us to get a "all or nothing" behavior.
*/
DefaultGeographicBoundingBox box = null;
boolean useExistingBox = (horizontalCRS != null);
final Collection<GeographicExtent> spatialExtents = target.getSpatialExtent();
final Iterator<GeographicExtent> it = spatialExtents.iterator();
while (it.hasNext()) {
final GeographicExtent extent = it.next();
if (extent instanceof GeographicBoundingBox) {
if (useExistingBox && (extent instanceof DefaultGeographicBoundingBox)) {
box = (DefaultGeographicBoundingBox) extent;
useExistingBox = false;
} else {
it.remove();
}
}
}
if (horizontalCRS != null) {
if (box == null) {
box = new DefaultGeographicBoundingBox();
spatialExtents.add(box);
}
GeographicCRS normalizedCRS = ReferencingUtilities.toNormalizedGeographicCRS(crs);
if (normalizedCRS == null) {
normalizedCRS = CommonCRS.defaultGeographic();
}
setGeographicExtent(envelope, box, crs, normalizedCRS);
}
/*
* Other dimensions (vertical and temporal).
*/
if (verticalCRS != null) {
VerticalExtent e = target.getVerticalExtent();
if (!(e instanceof DefaultVerticalExtent)) {
e = new DefaultVerticalExtent();
target.setVerticalExtent(e);
}
setVerticalExtent(envelope, (DefaultVerticalExtent) e, crs, verticalCRS);
} else {
target.setVerticalExtent(null);
}
if (temporalCRS != null) {
setTemporalExtent(envelope, target, crs, temporalCRS);
} else {
target.setExtent(null);
}
}
Aggregations