This is an automated email from the ASF dual-hosted git repository.
asf-gitbox-commits pushed a commit to branch geoapi-4.0
in repository https://gitbox.apache.org/repos/asf/sis.git
The following commit(s) were added to refs/heads/geoapi-4.0 by this push:
new dc2ef55d6f Revisit the way that CRS are created from Grid Mapping
attributes in netCDF files. This is an attempt to resolve confusion about
expected axis order and number of dimensions.
dc2ef55d6f is described below
commit dc2ef55d6fa2f5f89c6b0f9404419351ee81d331
Author: Martin Desruisseaux <[email protected]>
AuthorDate: Fri Oct 23 03:31:38 2026 +0200
Revisit the way that CRS are created from Grid Mapping attributes in netCDF
files.
This is an attempt to resolve confusion about expected axis order and
number of dimensions.
---
.../sis/coverage/grid/ResampledGridCoverage.java | 2 +-
.../apache/sis/profile/japan/netcdf/GCOM_C.java | 20 +-
.../main/org/apache/sis/referencing/CRS.java | 22 +-
.../sis/referencing/operation/matrix/Matrices.java | 37 +-
.../referencing/operation/matrix/package-info.java | 2 +-
.../referencing/operation/projection/Mercator.java | 4 +-
.../operation/transform/MathTransforms.java | 4 +-
.../main/module-info.java | 2 +-
.../apache/sis/storage/netcdf/base/Convention.java | 33 +-
.../apache/sis/storage/netcdf/base/Decoder.java | 2 +
.../sis/storage/netcdf/base/GridMapping.java | 718 ++++++++++++---------
.../sis/storage/netcdf/internal/Resources.java | 6 +-
.../storage/netcdf/internal/Resources.properties | 2 +-
.../netcdf/internal/Resources_fr.properties | 2 +-
.../apache/sis/storage/netcdf/package-info.java | 2 +-
.../apache/sis/storage/netcdf/NetcdfStoreTest.java | 4 +-
.../org/apache/sis/storage/esri/WritableStore.java | 2 +-
.../main/org/apache/sis/math/Vector.java | 4 +
18 files changed, 522 insertions(+), 346 deletions(-)
diff --git
a/endorsed/src/org.apache.sis.feature/main/org/apache/sis/coverage/grid/ResampledGridCoverage.java
b/endorsed/src/org.apache.sis.feature/main/org/apache/sis/coverage/grid/ResampledGridCoverage.java
index 237401e492..6187164bc7 100644
---
a/endorsed/src/org.apache.sis.feature/main/org/apache/sis/coverage/grid/ResampledGridCoverage.java
+++
b/endorsed/src/org.apache.sis.feature/main/org/apache/sis/coverage/grid/ResampledGridCoverage.java
@@ -157,7 +157,7 @@ final class ResampledGridCoverage extends
DerivedGridCoverage {
*/
private static long[] getIntegerTranslation(final MathTransform toSource) {
final Matrix m = MathTransforms.getMatrix(toSource);
- if (m == null || !Matrices.isTranslation(m)) {
+ if (!Matrices.isTranslation(m)) {
return null;
}
final int tc = m.getNumCol() - 1;
diff --git
a/endorsed/src/org.apache.sis.profile.japan/main/org/apache/sis/profile/japan/netcdf/GCOM_C.java
b/endorsed/src/org.apache.sis.profile.japan/main/org/apache/sis/profile/japan/netcdf/GCOM_C.java
index abfeebc8a6..0c924ca985 100644
---
a/endorsed/src/org.apache.sis.profile.japan/main/org/apache/sis/profile/japan/netcdf/GCOM_C.java
+++
b/endorsed/src/org.apache.sis.profile.japan/main/org/apache/sis/profile/japan/netcdf/GCOM_C.java
@@ -373,7 +373,7 @@ public final class GCOM_C extends Convention {
* @return the map projection definition as a modifiable map, or {@code
null} if none.
*/
@Override
- public Map<String,Object> projection(final Node node) {
+ public Map<String, Object> projection(final Node node) {
final String name = node.getAttributeAsString("Image_projection");
if (name == null) {
return super.projection(node);
@@ -390,14 +390,14 @@ public final class GCOM_C extends Convention {
} else {
return super.projection(node);
}
- final var definition = new HashMap<String,Object>(4);
+ final var definition = new HashMap<String, Object>(4);
definition.put(CF.GRID_MAPPING_NAME, method);
definition.put(CONVERSION_NAME, name);
return definition;
}
/**
- * The attributes storing values of the 4 corners in degrees with
(latitude, longitude) axis order.
+ * The attributes storing values of the 4 corners in degrees with
(latitude,longitude) axis order.
* This is used by {@link #projection(Node)} for inferring a "grid to CRS"
transform.
*/
private static final String[] CORNERS = {
@@ -412,22 +412,22 @@ public final class GCOM_C extends Convention {
};
/**
- * Returns the <i>grid corners to CRS</i> transform for the given node.
+ * Returns the <i>grid corners to <abbr>CRS</abbr></i> transform for the
given node.
* This method is invoked after call to {@link #projection(Node)} resulted
in creation of a projected CRS.
- * The {@linkplain ProjectedCRS#getBaseCRS() base CRS} shall have
(latitude, longitude) axes in degrees.
+ * The {@linkplain ProjectedCRS#getBaseCRS() base CRS} shall have
(latitude,longitude) axes in degrees.
*
- * @param node the same node as the one given to {@link
#projection(Node)}.
- * @param baseToCRS conversion from (latitude, longitude) in degrees to
the projected CRS.
+ * @param node the same node as the one given to {@link
#projection(Node)}.
+ * @param geographicToCRS conversion from (latitude,longitude) in
degrees to the projected CRS.
* @return the "grid corner to CRS" transform, or {@code null} if none or
unknown.
* @throws TransformException if a coordinate operation was required but
failed.
*/
@Override
- public MathTransform gridToCRS(final Node node, final MathTransform
baseToCRS) throws TransformException {
+ public MathTransform gridToCRS(final Node node, final MathTransform
geographicToCRS) throws TransformException {
final double[] corners = new double[CORNERS.length];
for (int i=0; i<corners.length; i++) {
corners[i] = node.getAttributeAsDouble(CORNERS[i]);
}
- baseToCRS.transform(corners, 0, corners, 0, corners.length / 2);
+ geographicToCRS.transform(corners, 0, corners, 0, corners.length / 2);
/*
* Compute spans of data (typically in metres) as the average of the
spans on both sides
* (width as length of top and bottom edges, height as length of left
and right edges).
@@ -448,7 +448,7 @@ public final class GCOM_C extends Convention {
m.m12 = corners[1];
return MathTransforms.linear(m);
}
- return super.gridToCRS(node, baseToCRS);
+ return super.gridToCRS(node, geographicToCRS);
}
/**
diff --git
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/CRS.java
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/CRS.java
index 5d0a008e63..195baa5cb3 100644
---
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/CRS.java
+++
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/CRS.java
@@ -92,6 +92,7 @@ import
org.apache.sis.metadata.iso.extent.DefaultGeographicBoundingBox;
import org.apache.sis.metadata.iso.extent.Extents;
import org.apache.sis.util.ArgumentChecks;
import org.apache.sis.util.ComparisonMode;
+import org.apache.sis.util.LenientComparable;
import org.apache.sis.util.OptionalCandidate;
import org.apache.sis.util.Utilities;
import org.apache.sis.util.internal.shared.Constants;
@@ -151,7 +152,7 @@ import org.opengis.coordinate.CoordinateMetadata;
*
* @author Martin Desruisseaux (IRD, Geomatys)
* @author Alexis Manin (Geomatys)
- * @version 1.7
+ * @version 1.8
* @since 0.3
*/
public final class CRS {
@@ -515,7 +516,19 @@ public final class CRS {
* @since 1.5
*/
public static boolean equivalent(final CoordinateReferenceSystem object1,
final CoordinateReferenceSystem object2) {
- return Utilities.deepEquals(object1, object2,
ComparisonMode.COMPATIBILITY);
+ if (object1 == object2) {
+ return true;
+ }
+ if (object1 == null || object2 == null) {
+ return false;
+ }
+ if (object1 instanceof LenientComparable) {
+ return ((LenientComparable) object1).equals(object2,
ComparisonMode.COMPATIBILITY);
+ }
+ if (object2 instanceof LenientComparable) {
+ return ((LenientComparable) object2).equals(object1,
ComparisonMode.COMPATIBILITY);
+ }
+ return object1.equals(object2);
}
/**
@@ -1220,9 +1233,8 @@ public final class CRS {
* then this method returns a {@link BitSet} in which all bits are set to
{@code true} in the range from 0 to
* <var>n</var>−1, where <var>n</var> is the number of dimensions of
{@code crs}.
* Otherwise, if {@code crs} is an instance of {@link CompoundCRS}, then
the {@code crs} components are traversed
- * recursively and compared in the same way as described above, except
that the range of bits set to {@code true}
- * does not start at index 0. Instead, the range starts at the index of
the first dimension of the matched
- * {@code crs} component. If no match is found, then the returned {@link
BitSet} is empty.
+ * recursively and compared in the same way as described above, with bit
shifted at the position corresponding to
+ * each {@code crs} component. If no match is found, then the returned
{@link BitSet} is empty.
*
* <h4>Example</h4>
* The following snippet gets the dimensions of the temporal and vertical
components of a <abbr>CRS</abbr>:
diff --git
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/matrix/Matrices.java
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/matrix/Matrices.java
index 65b7abe5e4..adacb13018 100644
---
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/matrix/Matrices.java
+++
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/matrix/Matrices.java
@@ -74,7 +74,7 @@ import org.opengis.coordinate.MismatchedDimensionException;
* </ul>
*
* @author Martin Desruisseaux (IRD, Geomatys)
- * @version 1.6
+ * @version 1.8
*
* @see org.apache.sis.parameter.MatrixParameters
*
@@ -1059,13 +1059,15 @@ search: while (freeColumn < numCol) {
* A transform is affine if the matrix is square and its last row contains
* only zeros, except in the last column which contains 1.
*
- * @param matrix the matrix to test.
- * @return {@code true} if the matrix represents an affine transform.
+ * @param matrix the matrix to test, or {@code null}.
+ * @return {@code true} if the given matrix is non-null and represents an
affine transform.
*
* @see MatrixSIS#isAffine()
*/
public static boolean isAffine(final Matrix matrix) {
- if (matrix instanceof MatrixSIS) {
+ if (matrix == null) {
+ return false;
+ } else if (matrix instanceof MatrixSIS) {
return ((MatrixSIS) matrix).isAffine();
} else {
return MatrixSIS.isAffine(matrix);
@@ -1077,8 +1079,8 @@ search: while (freeColumn < numCol) {
* This method returns {@code true} if the given matrix {@linkplain
#isAffine(Matrix) is affine}
* and differs from the identity matrix only in the last column.
*
- * @param matrix the matrix to test.
- * @return {@code true} if the matrix represents a translation.
+ * @param matrix the matrix to test, or {@code null}.
+ * @return {@code true} if the given matrix is non-null and represents a
translation.
*
* @since 0.7
*/
@@ -1109,13 +1111,16 @@ search: while (freeColumn < numCol) {
* errors. In case of doubt, it is often safer to use the strict {@link
MatrixSIS#isIdentity()} method instead
* than this one.</p>
*
- * @param matrix the matrix to test for identity.
+ * @param matrix the matrix to test for identity, or {@code null}.
* @param tolerance the tolerance value, or 0 for a strict comparison.
- * @return {@code true} if this matrix is close to the identity matrix
given the tolerance threshold.
+ * @return {@code true} if this matrix is non-null and close to the
identity matrix given the tolerance threshold.
*
* @see MatrixSIS#isIdentity()
*/
public static boolean isIdentity(final Matrix matrix, final double
tolerance) {
+ if (matrix == null) {
+ return false;
+ }
final int numRow = matrix.getNumRow();
final int numCol = matrix.getNumCol();
if (numRow != numCol) {
@@ -1138,18 +1143,20 @@ search: while (freeColumn < numCol) {
/**
* Returns whether the given matrix has any NaN value.
*
- * @param matrix the matrix to test.
+ * @param matrix the matrix to test, or {@code null}.
* @return {@code true} if at least one matrix element is NaN.
*
* @since 1.6
*/
public static boolean hasNaN(final Matrix matrix) {
- final int numCol = matrix.getNumCol();
- final int numRow = matrix.getNumRow();
- for (int j=0; j<numRow; j++) {
- for (int i=0; i<numCol; i++) {
- if (Double.isNaN(matrix.getElement(j, i))) {
- return true;
+ if (matrix != null) {
+ final int numCol = matrix.getNumCol();
+ final int numRow = matrix.getNumRow();
+ for (int j=0; j<numRow; j++) {
+ for (int i=0; i<numCol; i++) {
+ if (Double.isNaN(matrix.getElement(j, i))) {
+ return true;
+ }
}
}
}
diff --git
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/matrix/package-info.java
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/matrix/package-info.java
index 1a6ed24747..0720d77c21 100644
---
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/matrix/package-info.java
+++
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/matrix/package-info.java
@@ -71,7 +71,7 @@
* Like SIS, Vecmath is optimized for small matrices of interest for 2D and 3D
graphics.</p>
*
* @author Martin Desruisseaux (IRD, Geomatys)
- * @version 1.7
+ * @version 1.8
* @since 0.4
*/
package org.apache.sis.referencing.operation.matrix;
diff --git
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/projection/Mercator.java
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/projection/Mercator.java
index 7259120223..968c7b1f68 100644
---
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/projection/Mercator.java
+++
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/projection/Mercator.java
@@ -185,6 +185,7 @@ public class Mercator extends ConformalProjection {
/**
* Workaround while waiting for JEP 513: Flexible Constructor Bodies.
*/
+ @SuppressWarnings("fallthrough")
@Workaround(library="JDK", version="11", fixed="25")
private static Initializer initializer(final OperationMethod method, final
Parameters parameters) {
final Variant variant = variant(method, Variant.values(),
Variant.TWO_PARALLELS);
@@ -622,8 +623,7 @@ subst: if (variant.spherical || eccentricity == 0) {
* (but no scale and no translation are allowed). Ignore the
longitude row because it just passes through
* this Mercator projection with no impact on any calculation.
*/
- return (m != null
- && Matrices.isAffine(m)
+ return (Matrices.isAffine(m)
&& m.getNumRow() == DIMENSION+1
&& m.getNumCol() == DIMENSION+1
&& m.getElement(1, DIMENSION) == 0
diff --git
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/transform/MathTransforms.java
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/transform/MathTransforms.java
index e15829a938..1a85580232 100644
---
a/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/transform/MathTransforms.java
+++
b/endorsed/src/org.apache.sis.referencing/main/org/apache/sis/referencing/operation/transform/MathTransforms.java
@@ -96,13 +96,13 @@ public final class MathTransforms {
* The most typical use case is for switching between
(<var>latitude</var>, <var>longitude</var>)
* and (<var>longitude</var>, <var>latitude</var>) axis order of
geographic <abbr>CRS</abbr>.
*
- * <h4>Warning: use in last resort only</h4>
+ * <p><strong>Warning: use in last resort only.</strong>
* While switching the axis order of geographic <abbr>CRS</abbr> seems a
very frequent problem,
* the need to invoke this method should be very rare. Apache
<abbr>SIS</abbr> already handles
* axis order changes when deriving Coordinate Operations between
Coordinate Reference Systems (<abbr>CRS</abbr>).
* If a user's code needs to change axis order manually, this is often the
symptom of a design problem such as
* data without {@link
org.opengis.referencing.crs.CoordinateReferenceSystem} or with an incorrect
<abbr>CRS</abbr>,
- * or computation that does not use another <abbr>CRS</abbr> for
describing their desired coordinate system.
+ * or computation that does not use another <abbr>CRS</abbr> for
describing their desired coordinate system.</p>
*
* @param dimension number of dimensions of the transform to be returned.
* @return a transform of the specified dimension which swaps the two
first axes.
diff --git a/endorsed/src/org.apache.sis.storage.netcdf/main/module-info.java
b/endorsed/src/org.apache.sis.storage.netcdf/main/module-info.java
index 1fcc4b3cd1..f745db414f 100644
--- a/endorsed/src/org.apache.sis.storage.netcdf/main/module-info.java
+++ b/endorsed/src/org.apache.sis.storage.netcdf/main/module-info.java
@@ -19,7 +19,7 @@
* NetCDF store.
*
* @author Martin Desruisseaux (IRD, Geomatys)
- * @version 1.4
+ * @version 1.8
* @since 0.3
*/
module org.apache.sis.storage.netcdf {
diff --git
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/Convention.java
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/Convention.java
index d45e630254..03f81923b2 100644
---
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/Convention.java
+++
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/Convention.java
@@ -492,7 +492,18 @@ public class Convention {
final String nameLC = name.toLowerCase(Decoder.DATA_LOCALE);
Object value;
switch (nameLC) {
- case CF.GRID_MAPPING_NAME: continue; // Already stored.
+ case CF.GRID_MAPPING_NAME: {
+ // Already stored.
+ continue;
+ }
+ case CF.SEMI_MAJOR_AXIS:
+ case CF.SEMI_MINOR_AXIS:
+ case CF.INVERSE_FLATTENING:
+ case LONGITUDE_OF_PRIME_MERIDIAN: {
+ value = node.getAttributeAsNumber(name);
+ if (value == null) continue;
+ break;
+ }
case TOWGS84: {
/*
* Conversion to WGS 84 datum may be specified as
Bursa-Wolf parameters. Encoding this information
@@ -509,18 +520,13 @@ public class Convention {
default: {
if (nameLC.endsWith(NAME_SUFFIX)) {
value = node.getAttributeAsString(name);
- if (value == null) continue;
} else {
- /*
- * Assume that all map projection parameters in netCDF
files are numbers or array of numbers.
- * If values are array, then they are converted to an
array of {@code double[]} type.
- */
value = node.getAttributeValue(name);
- if (value == null) continue;
if (value instanceof Vector) {
value = ((Vector) value).doubleValues();
}
}
+ if (value == null) continue;
break;
}
}
@@ -535,9 +541,10 @@ public class Convention {
}
/**
- * Returns the <i>grid corner to CRS</i> transform for the given node.
This method is invoked after call
- * to {@link #projection(Node)} method resulted in creation of a projected
coordinate reference system.
- * The {@linkplain ProjectedCRS#getBaseCRS() base CRS} is fixed to
(latitude, longitude) axes in degrees,
+ * Returns the <i>grid corner to <abbr>CRS</abbr></i> transform for the
given node.
+ * This method is invoked after a call to {@link #projection(Node)} method
resulted
+ * in the creation of a geographic or projected coordinate reference
system.
+ * The {@linkplain ProjectedCRS#getBaseCRS() base CRS} is fixed to
(latitude,longitude) axes in degrees,
* but the projected CRS axes may have any order and units. In the
particular case of "latitude_longitude"
* pseudo-projection, the "projected" CRS is actually a {@link
GeographicCRS} instance.
* The returned transform, if non-null, shall map cell corners.
@@ -550,12 +557,12 @@ public class Convention {
* <li>Base CRS axis order is (latitude, longitude) for increasing the
chances to have a CRS identified by EPSG.</li>
* </ul>
*
- * @param node the same node as the one given to {@link
#projection(Node)}.
- * @param baseToCRS conversion from (latitude, longitude) in degrees to
the projected CRS.
+ * @param node the same node as the one given to {@link
#projection(Node)}.
+ * @param geographicToCRS conversion from (latitude,longitude) in
degrees to the projected CRS.
* @return the <i>grid corner to CRS</i> transform, or {@code null} if
none or unknown.
* @throws TransformException if a coordinate operation was required but
failed.
*/
- public MathTransform gridToCRS(final Node node, final MathTransform
baseToCRS) throws TransformException {
+ public MathTransform gridToCRS(Node node, MathTransform geographicToCRS)
throws TransformException {
return null;
}
diff --git
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/Decoder.java
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/Decoder.java
index 18b8f64d17..33bd11a5e8 100644
---
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/Decoder.java
+++
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/Decoder.java
@@ -249,6 +249,8 @@ public abstract class Decoder extends
ReferencingFactoryContainer {
* <li>{@link MetadataBuilder#setPredefinedFormat(String,
StoreListeners, boolean)}</li>
* <li>{@link MetadataBuilder#addFormatReaderSIS(String)} (if
applicable)</li>
* </ul>
+ *
+ * @param builder the builder where to add metadata.
*/
public abstract void addFormatDescription(MetadataBuilder builder);
diff --git
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/GridMapping.java
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/GridMapping.java
index 2699f6a977..c4c6577c81 100644
---
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/GridMapping.java
+++
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/base/GridMapping.java
@@ -21,6 +21,7 @@ import java.util.List;
import java.util.ArrayList;
import java.util.HashMap;
import java.util.Iterator;
+import java.util.Optional;
import java.util.Locale;
import java.util.logging.Level;
import java.util.logging.LogRecord;
@@ -44,11 +45,11 @@ import org.opengis.referencing.crs.CRSFactory;
import org.opengis.referencing.crs.ProjectedCRS;
import org.opengis.referencing.crs.GeographicCRS;
import org.opengis.referencing.crs.CoordinateReferenceSystem;
+import org.opengis.referencing.operation.OperationNotFoundException;
import org.opengis.referencing.operation.TransformException;
import org.opengis.referencing.operation.OperationMethod;
import org.opengis.referencing.operation.MathTransform;
import org.opengis.referencing.operation.Conversion;
-import org.opengis.referencing.operation.Matrix;
import org.opengis.referencing.datum.DatumFactory;
import org.opengis.referencing.datum.GeodeticDatum;
import org.opengis.referencing.datum.PrimeMeridian;
@@ -80,6 +81,7 @@ import org.apache.sis.util.ArraysExt;
import org.apache.sis.util.Utilities;
import org.apache.sis.util.Exceptions;
import org.apache.sis.util.internal.shared.Strings;
+import org.apache.sis.util.internal.shared.Numerics;
import org.apache.sis.util.internal.shared.Constants;
import org.apache.sis.util.resources.Vocabulary;
import org.apache.sis.util.resources.Errors;
@@ -87,11 +89,9 @@ import org.apache.sis.util.resources.IndexedResourceBundle;
import org.apache.sis.io.wkt.WKTFormat;
import org.apache.sis.io.wkt.Warnings;
import org.apache.sis.math.NumberType;
+import org.apache.sis.math.Vector;
import org.apache.sis.measure.Units;
-// Specific to the geoapi-3.1 and geoapi-4.0 branches:
-import org.opengis.referencing.datum.DatumEnsemble;
-
// Specific to the geoapi-4.0 branch:
import org.opengis.referencing.crs.DerivedCRS;
@@ -120,13 +120,19 @@ final class GridMapping {
private static final String CRS_WKT = "crs_wkt", SPATIAL_REF =
"spatial_ref";
/**
- * Name of attributes where the "grid to <abbr>CRS</abbr> transform may be
encoded as an affine transform.
+ * Other attribute names, may be in mixed case.
+ */
+ private static final String EPSG_CODE = "EPSG_code";
+
+ /**
+ * Name of attributes where the "grid to <abbr>CRS</abbr>" transform may
be encoded as an affine transform.
* The {@code GeoTransform} attribute is specific to <abbr>GDAL</abbr>. It
uses pixel-corner convention and
* interprets data as if it was an image (as if the row shown on the top
had index 0), ignoring the netCDF
* cell indices (where row 0 is often in the bottom).
*
* @see #gridToCRS
* @see #SOURCE_AXIS_TO_FLIP
+ * @see #parseGeoTransform(Object)
*/
private static final String GEOTRANSFORM = "GeoTransform";
@@ -137,6 +143,11 @@ final class GridMapping {
*/
private static final int SOURCE_AXIS_TO_FLIP = 1;
+ /**
+ * Number of points to sample for consistency check. Used only for
checking if we should log a warning.
+ */
+ private static final int NUM_SAMPLING_POINTS = 5;
+
/**
* The variable on which projection parameters are defined as attributes.
* This is typically an empty variable referenced by the value of the
@@ -149,27 +160,34 @@ final class GridMapping {
/**
* The Coordinate Reference System inferred from grid mapping attribute
values, or {@code null} if none.
* This <abbr>CRS</abbr> may have been constructed from Well Known Text or
<abbr>EPSG</abbr> codes
- * declared in {@value #SPATIAL_REF}, {@code "ESRI_pe_string"} or {@code
"EPSG_code"} attributes.
- *
- * <h4>Usage note</h4>
+ * declared in {@value #SPATIAL_REF}, {@code "ESRI_pe_string"} or {@value
#EPSG_CODE} attributes.
* This is built from different information than the one used by {@link
CRSBuilder},
* which creates <abbr>CRS</abbr> by inspection of coordinate system axes.
*
+ * <p>If this <abbr>CRS</abbr> is geographic, its axis order may be
(latitude, longitude).
+ * This order may or may not be consistent with the axis order of netCDF
variables using this CRS.
+ * Each netCDF variable can have its own axis order and be associated to
the same <abbr>CRS</abbr>,
+ * so users of this <abbr>CRS</abbr> need to check anyway for change of
axis order.</p>
+ *
* @see #crs()
*/
private CoordinateReferenceSystem crs;
/**
- * The <i>grid corner to CRS</i> transform, or {@code null} if none.
+ * The two-dimensional <i>grid corner to <abbr>CRS</abbr></i> transform,
or {@code null} if none.
* This information is usually not specified except when using
<abbr>GDAL</abbr> conventions.
+ * The source coordinates are (column, row) indexes. The target
coordinates are in the CRS
+ * (including axis order and units) specified by the horizontal component
of {@link #crs}.
* If {@code null}, then the transform should be inferred by {@link Grid}.
*
* <h4>Image flip</h4>
* If the {@code gridToCRS} transform has been specified by the
<abbr>GDAL</abbr>'s {@value #GEOTRANSFORM}
* attribute, then it needs to have the <var>y</var> axis flipped. This is
because <abbr>GDAL</abbr> seems
* to ignore the netCDF cell coordinate system and to handle the data has
if it was an image with the last
- * row (in netCDF order) shown on the top.
+ * row (in netCDF order) shown on the top. Flipping cannot be done at
{@value #GEOTRANSFORM} reading time,
+ * because it depends on the shape of the variable associated with the
transform.
*
+ * @see #GEOTRANSFORM
* @see #SOURCE_AXIS_TO_FLIP
* @see #gridToCRS(Variable)
*/
@@ -182,9 +200,10 @@ final class GridMapping {
private PixelInCell anchorFromConvention;
/**
- * Whether the {@link #crs} was defined by a WKT string.
+ * The attribute that provided the {@link #gridToCRS}, or {@code null} if
none.
+ * Used for debugging purposes and for producing warning messages.
*/
- private boolean isWKT;
+ private String sourceOfGridToCRS;
/**
* Creates an initially empty instance.
@@ -244,19 +263,27 @@ final class GridMapping {
* This method tries to parse <abbr>CF</abbr>-compliant attributes,
potentially mixed with
* non-standard extensions (for example <abbr>GDAL</abbr> and
<abbr>ESRI</abbr>).
*
+ * <p>If the <abbr>CRS</abbr> is defined in many ways, an older version of
CF-conventions
+ * said that grid mapping attributes have precedence over definition by
<abbr>WKT</abbr>.
+ * A later revision of CF-conventions (around 2020) relaxed that
requirement.
+ * In case of inconsistent definitions, the choice of which definition
prevails
+ * is implementation dependent.</p>
+ *
* @param variable the variable from which to get the map projection
parameters.
* @return grid geometry information, or {@code null} if not found or if
an error occurred.
+ *
+ * @see <a
href="https://github.com/cf-convention/cf-conventions/issues/222";>Allow CRS WKT
+ * to represent the CRS without requiring reader to compare with grid
mapping parameters</a>
*/
@SuppressWarnings("UseSpecificCatch")
private static GridMapping tryAllConventions(final Node mapping) {
var gm = new GridMapping(mapping);
- // Tries CF-convention first, and if it doesn't work, try GDAL
convention.
- int i = 0;
+ int convention = 0;
boolean stop = false;
Exception warning = null;
do {
try {
- switch (i++) {
+ switch (convention++) {
case 0: stop = gm.parseProjectionParameters(); break;
case 1: stop = gm.parseGDAL(); break;
case 2: stop = gm.parseESRI(); break;
@@ -299,26 +326,33 @@ final class GridMapping {
*/
private boolean parseProjectionParameters() throws FactoryException,
TransformException {
final Decoder decoder = mapping.decoder;
- final Map<String, Object> definition =
decoder.convention().projection(mapping);
+ final Convention convention = decoder.convention();
+ final Map<String, Object> definition = convention.projection(mapping);
if (definition == null) {
return false;
}
/*
- * Search in advance for a CRS that we can use as a fallback. We need
only one such CRS for now,
- * the other CRS definitions will be parsed later in this method. We
need this fallback in order
- * to provide default names to geodetic objects in the common case
where the `definition` map
- * does not contain these names.
+ * Parse the WKT first because we will use some of their components
(in particular their names)
+ * as fallbacks if not specified in the grid mapping attributes. The
WKT are redundant with the
+ * grid mapping attributes, and the CF-conventions does not said which
one has precedence.
+ * This omission is intentional, so we have a bit of flexibility for
implementation convenience.
*/
- final var alreadyParsedWKT = new ArrayList<String>(2);
- if (crs == null) setOrVerifyWKT(alreadyParsedWKT, definition, CRS_WKT);
- if (crs == null) setOrVerifyWKT(alreadyParsedWKT, definition,
SPATIAL_REF);
+ { // For keeping variable in local scope.
+ final var alreadyParsedWKT = new ArrayList<String>(2);
+ setOrVerifyWKT(alreadyParsedWKT, definition, CRS_WKT);
+ setOrVerifyWKT(alreadyParsedWKT, definition, SPATIAL_REF);
+ }
final CoordinateReferenceSystem fromWKT = crs;
/*
- * Fetch now numerical values that are not map projection parameters.
- * This step needs to be done before to try to set parameter values.
+ * Fetch and remove the attributes that are not map projection
parameters.
+ * This step needs to be done before to iterate over remaining
parameters.
+ * Note that `gridMappingName` cannot be null as per `definition`
contract.
*/
+ parseGeoTransform(definition.remove("geotransform")); //
"GeoTransform" made lower-case.
final Object greenwichLongitude =
definition.remove(Convention.LONGITUDE_OF_PRIME_MERIDIAN);
- final String mappingName = (String)
definition.remove(CF.GRID_MAPPING_NAME);
+ final String gridMappingName =
definition.remove(CF.GRID_MAPPING_NAME).toString();
+ final Object crsIdentifier = definition.remove("epsg_code"); //
"EPSG_code" made lower-case.
+ final Object variableName = definition.remove(CF.LONG_NAME);
/*
* Prepare the group of projection parameters. The set of legal
parameter depends on the map projection.
* We assume that all numerical values are map projection parameters.
Character sequences (assumed to be
@@ -326,7 +360,8 @@ final class GridMapping {
* than OGC names, but Apache SIS implementations of map projections
know how to handle them, including
* the redundant parameters like "inverse_flattening" and
"earth_radius".
*/
- final OperationMethod method =
decoder.findOperationMethod(mappingName);
+ final var properties = new HashMap<String, Object>(4);
+ final OperationMethod method =
decoder.findOperationMethod(gridMappingName);
final ParameterValueGroup parameters =
method.getParameters().createValue();
for (final Iterator<Map.Entry<String, Object>> it =
definition.entrySet().iterator(); it.hasNext();) {
final Map.Entry<String, Object> entry = it.next();
@@ -334,27 +369,16 @@ final class GridMapping {
final Object value = entry.getValue();
try {
if (value instanceof Number || value instanceof double[] ||
value instanceof float[]) {
- it.remove();
parameters.parameter(name).setValue(value);
} else if (value instanceof String) {
final var text = (String) value;
if (name.endsWith(Convention.NAME_SUFFIX)) {
continue;
}
- switch (name) {
- case CRS_WKT:
- case SPATIAL_REF: continue; // Will be parsed
after this loop.
- case "geotransform": { // "GeoTransform" made
lower-case.
- if (parseGeoTransform(text)) {
- it.remove();
- }
- continue;
- }
- }
/*
- * In principle, we should ignore non-numeric parameters.
But in practice, some badly encoded
- * netCDF files store parameters as strings instead of
numbers. If the parameter name is known
- * to the projection method, try to parse the character
string.
+ * Some badly encoded netCDF files store parameters as
strings instead of numbers.
+ * If the parameter name is known to the projection method
and declares a numeric
+ * type, try to parse the character string.
*/
final ParameterValue<?> parameter;
try {
@@ -365,12 +389,14 @@ final class GridMapping {
}
final Class<?> type =
parameter.getDescriptor().getValueClass();
if (NumberType.isReal(type)) {
- it.remove();
parameter.setValue(Double.parseDouble(text));
} else if (NumberType.isReal(type.getComponentType())) {
- it.remove();
parameter.setValue(parseDoubles(text), null);
+ } else {
+ parameter.setValue(value);
}
+ } else { // Other types such as Bursa-Wolf parameters.
+ continue;
}
} catch (IllegalArgumentException ex) { // Includes
NumberFormatException.
warning(mapping,
@@ -383,119 +409,106 @@ final class GridMapping {
value,
ex.getLocalizedMessage());
}
+ it.remove();
}
/*
- * In principle, projection parameters do not include the semi-major
and semi-minor axis lengths.
- * But if those information are provided, then we use them for
building the geodetic reference frame.
- * Otherwise, a default reference frame will be used.
+ * The remaining of this method will create prime meridian, ellipsoid,
datum, geographic and projected CRS.
+ * The ellipsoid needs the semi-major and semi-minor values, which is
why above parameter parsing was first.
+ * The `defaultDefinitions` enumeration value is used for providing
default objects if there is no WKT.
*/
- final CommonCRS defaultDefinitions =
decoder.convention().defaultHorizontalCRS(false);
+ final CommonCRS defaultDefinitions =
convention.defaultHorizontalCRS(false);
final DatumFactory datumFactory = decoder.getDatumFactory();
- boolean hasBuiltSomeCustomObjects = false;
- /*
- * Prime meridian built from "longitude_of_prime_meridian".
- */
- final PrimeMeridian meridian;
- if (greenwichLongitude instanceof Number) {
- final double longitude = ((Number)
greenwichLongitude).doubleValue();
- final Map<String,?> properties = properties(definition, false,
Convention.PRIME_MERIDIAN_NAME, () -> {
- // Fallback if `definition` does not contain a name for the
prime meridian.
- PrimeMeridian template =
DatumOrEnsemble.getPrimeMeridian(fromWKT).orElse(null);
- if (template == null) {
- if (longitude != 0) return null;
- template = defaultDefinitions.primeMeridian();
+ final GeodeticDatum datum;
+ { // For keeping variables in local scope.
+ PrimeMeridian meridian = null;
+ if (greenwichLongitude != null) {
+ final double longitude = ((Number)
greenwichLongitude).doubleValue();
+ final Optional<PrimeMeridian> template =
DatumOrEnsemble.getPrimeMeridian(fromWKT);
+ if (template.isPresent() || longitude != 0) {
+ initializeProperties(properties, definition,
Convention.PRIME_MERIDIAN_NAME, () -> {
+ // Fallback if `definition` does not contain a name
for the prime meridian.
+ return
template.orElseGet(defaultDefinitions::primeMeridian).getName();
+ });
+ meridian = datumFactory.createPrimeMeridian(properties,
longitude, Units.DEGREE);
}
- return template.getName();
- });
- meridian = datumFactory.createPrimeMeridian(properties, longitude,
Units.DEGREE);
- hasBuiltSomeCustomObjects = true;
- } else {
- meridian = defaultDefinitions.primeMeridian();
- }
- /*
- * Ellipsoid built from "semi_major_axis" and "semi_minor_axis"
parameters. Note that it is okay
- * to use the OGC name (e.g. "semi_major") instead of the netCDF name
(e.g. ""semi_major_axis").
- * The Apache SIS implementation of parameter value group understands
the aliases. Using the OGC
- * names is safer because they should be understood by most map
projection implementations.
- */
- Ellipsoid ellipsoid;
- try {
- final ParameterValue<?> p =
parameters.parameter(Constants.SEMI_MAJOR);
- final Unit<Length> axisUnit = p.getUnit().asType(Length.class);
- final double semiMajor = p.doubleValue();
- boolean isIvfDefinitive;
- try {
- isIvfDefinitive =
parameters.parameter(Constants.IS_IVF_DEFINITIVE).booleanValue();
- } catch (ParameterNotFoundException e) {
- // Ignore - may be normal if the map projection is not an
Apache SIS implementation.
- isIvfDefinitive = false;
- }
- final double secondDefiningParameter;
- final boolean isSphere;
- if (isIvfDefinitive) {
- secondDefiningParameter =
parameters.parameter(Constants.INVERSE_FLATTENING).doubleValue();
- isSphere = (secondDefiningParameter == 0) ||
Double.isInfinite(secondDefiningParameter);
- } else {
- secondDefiningParameter =
parameters.parameter(Constants.SEMI_MINOR).doubleValue(axisUnit);
- isSphere = secondDefiningParameter == semiMajor;
}
- final Map<String,?> properties = properties(definition, false,
Convention.ELLIPSOID_NAME, () -> {
- // Fallback if `definition` does not contain a name for the
ellipsoid.
- return
DatumOrEnsemble.getEllipsoid(fromWKT).<Object>map(Ellipsoid::getName).orElseGet(()
-> {
- final Locale locale = decoder.getLocale();
- final String name =
Vocabulary.forLocale(locale).getString(isSphere ? Vocabulary.Keys.Sphere :
Vocabulary.Keys.Ellipsoid);
- final NumberFormat f =
NumberFormat.getNumberInstance(locale);
- f.setMaximumFractionDigits(5); // Centimetric
precision.
- return
f.format(axisUnit.getConverterTo(Units.KILOMETRE).convert(semiMajor),
- new StringBuffer(name).append(isSphere ? "
R=" : " a="),
- new FieldPosition(0))
- .append(" km").toString();
+ /*
+ * Ellipsoid built from "semi_major_axis" and "semi_minor_axis"
parameters. Note that it is okay
+ * to use the OGC name (e.g. "semi_major") instead of the netCDF
name (e.g. ""semi_major_axis").
+ * The Apache SIS implementation of parameter value group
understands the aliases. Using the OGC
+ * names is safer because they should be understood by most map
projection implementations.
+ */
+ Ellipsoid ellipsoid = null;
+ try {
+ final ParameterValue<?> p =
parameters.parameter(Constants.SEMI_MAJOR);
+ final Unit<Length> axisUnit = p.getUnit().asType(Length.class);
+ final double semiMajor = p.doubleValue();
+ boolean isIvfDefinitive;
+ try {
+ isIvfDefinitive =
parameters.parameter(Constants.IS_IVF_DEFINITIVE).booleanValue();
+ } catch (ParameterNotFoundException e) {
+ // Ignore - may be normal if the map projection is not an
Apache SIS implementation.
+ ignorable(e);
+ isIvfDefinitive = false;
+ }
+ final double secondDefiningParameter;
+ final boolean isSphere;
+ if (isIvfDefinitive) {
+ secondDefiningParameter =
parameters.parameter(Constants.INVERSE_FLATTENING).doubleValue();
+ isSphere = (secondDefiningParameter == 0) ||
Double.isInfinite(secondDefiningParameter);
+ } else {
+ secondDefiningParameter =
parameters.parameter(Constants.SEMI_MINOR).doubleValue(axisUnit);
+ isSphere = secondDefiningParameter == semiMajor;
+ }
+ initializeProperties(properties, definition,
Convention.ELLIPSOID_NAME, () -> {
+ // Fallback if `definition` does not contain a name for
the ellipsoid.
+ return
DatumOrEnsemble.getEllipsoid(fromWKT).<Object>map(Ellipsoid::getName).orElseGet(()
-> {
+ final Locale locale = decoder.getLocale();
+ final String name =
Vocabulary.forLocale(locale).getString(isSphere ? Vocabulary.Keys.Sphere :
Vocabulary.Keys.Ellipsoid);
+ final NumberFormat f =
NumberFormat.getNumberInstance(locale);
+ f.setMaximumFractionDigits(5); // Centimetric
precision.
+ return
f.format(axisUnit.getConverterTo(Units.KILOMETRE).convert(semiMajor),
+ new StringBuffer(name).append(isSphere
? " R=" : " a="),
+ new FieldPosition(0))
+ .append(" km").toString();
+ });
});
- });
- if (isIvfDefinitive) {
- ellipsoid = datumFactory.createFlattenedSphere(properties,
semiMajor, secondDefiningParameter, axisUnit);
- } else {
- ellipsoid = datumFactory.createEllipsoid(properties,
semiMajor, secondDefiningParameter, axisUnit);
+ if (isIvfDefinitive) {
+ ellipsoid = datumFactory.createFlattenedSphere(properties,
semiMajor, secondDefiningParameter, axisUnit);
+ } else {
+ ellipsoid = datumFactory.createEllipsoid(properties,
semiMajor, secondDefiningParameter, axisUnit);
+ }
+ } catch (IllegalStateException e) {
+ warningInMapping(mapping, e,
Resources.Keys.MissingEllipsoid_3, e.getLocalizedMessage());
}
- hasBuiltSomeCustomObjects = true;
- } catch (IllegalStateException e) {
- warningInMapping(mapping, e, Resources.Keys.MissingEllipsoid_3,
e.getLocalizedMessage());
- ellipsoid = defaultDefinitions.ellipsoid();
- }
- /*
- * Geodetic reference frame built from "towgs84" and above properties.
- * The class of the "towgs84" entry will be verified by the datum
constuctor.
- */
- final GeodeticDatum datum;
- DatumEnsemble<GeodeticDatum> ensemble = null;
- final Object bursaWolf = definition.remove(Convention.TOWGS84);
- if (hasBuiltSomeCustomObjects || bursaWolf != null) {
- Map<String, Object> properties = properties(definition, false,
Convention.GEODETIC_DATUM_NAME, () -> {
- // Fallback if `definition` does not contain a name for the
geodetic reference frame.
- return
CRS.getGeodeticReferenceFrame(fromWKT).map(GeodeticDatum::getName).orElse(null);
- });
- if (bursaWolf != null) {
- properties = new HashMap<>(properties);
+ /*
+ * Geodetic reference frame built from "towgs84" and above
properties.
+ * The class of the "towgs84" entry will be verified by the datum
constuctor.
+ */
+ final Object bursaWolf = definition.remove(Convention.TOWGS84);
+ if (bursaWolf != null || ellipsoid != null || meridian != null) {
+ if (ellipsoid == null) ellipsoid =
defaultDefinitions.ellipsoid();
+ if (meridian == null) meridian =
defaultDefinitions.primeMeridian();
+ initializeProperties(properties, definition,
Convention.GEODETIC_DATUM_NAME, () -> {
+ // Fallback if `definition` does not contain a name for
the geodetic reference frame.
+ return
CRS.getGeodeticReferenceFrame(fromWKT).map(GeodeticDatum::getName).orElse(null);
+ });
properties.put(DefaultGeodeticDatum.BURSA_WOLF_KEY, bursaWolf);
- hasBuiltSomeCustomObjects = true;
- }
- datum = datumFactory.createGeodeticDatum(properties, ellipsoid,
meridian);
- } else {
- datum = defaultDefinitions.datum(false);
- if (datum == null) {
- ensemble = defaultDefinitions.datumEnsemble();
+ datum = datumFactory.createGeodeticDatum(properties,
ellipsoid, meridian);
+ } else {
+ datum = null;
}
}
/*
- * Geographic or projected CRS built from above properties.
- * The geographic CRS will always have (latitude, longitude) axes in
that order and in degrees.
- * The swapping to (longitude, latitude) axis order will be done by
the `baseToCRS` transform.
+ * Base of the map projection (including the geographic case, handled
as Pseudo Plate Carrée).
+ * The axis order of `baseCRS` is (latitude, longitude) even if
recreated because it uses the
+ * `defaultDefinitions.geographic()` coordinate system which is
invariant by contract.
*/
final boolean wantGeographicCRS = (method instanceof
PseudoPlateCarree);
GeographicCRS baseCRS = defaultDefinitions.geographic();
final CRSFactory crsFactory = decoder.getCRSFactory();
- if (hasBuiltSomeCustomObjects) {
- final Map<String,?> properties = properties(definition,
wantGeographicCRS, Convention.GEOGRAPHIC_CRS_NAME, () -> {
+ if (datum != null) {
+ initializeProperties(properties, definition,
Convention.GEOGRAPHIC_CRS_NAME, () -> {
// Fallback if `definition` does not contain a name for the
geodetic CRS.
IdentifiedObject base = CRS.getHorizontalComponent(fromWKT);
if (base == null) {
@@ -505,10 +518,17 @@ final class GridMapping {
}
return base.getName();
});
- baseCRS = crsFactory.createGeographicCRS(properties, datum,
ensemble, baseCRS.getCoordinateSystem());
+ if (wantGeographicCRS) {
+ properties.put(GeographicCRS.REMARKS_KEY,
definition.remove(ACDD.comment));
+ }
+ baseCRS = crsFactory.createGeographicCRS(properties, datum, null,
baseCRS.getCoordinateSystem());
}
- // Only swap axis order from (latitude, longitude) to (longitude,
latitude).
- MathTransform baseToCRS = MathTransforms.swapTwoFirstAxes(2);
+ /*
+ * Geographic or projected CRS built from above objects declared in
the grid mapping variable.
+ * The geographic CRS will always have (latitude, longitude) axes in
that order and in degrees.
+ * We use that order because it will be the order that users see by
default in user interfaces.
+ * It also increases the chances that we can assign an EPSG code to
the CRS.
+ */
if (wantGeographicCRS) {
crs = baseCRS;
} else {
@@ -518,69 +538,103 @@ final class GridMapping {
* Reminder: this parameter group has been created from a subset
of `definition` at the
* beginning of this method.
*/
- Map<String,?> properties = properties(definition, false,
Convention.CONVERSION_NAME, () -> {
+ initializeProperties(properties, definition,
Convention.CONVERSION_NAME, () -> {
if (fromWKT instanceof ProjectedCRS) {
return ((ProjectedCRS)
fromWKT).getConversionFromBase().getName();
}
- return mapping.getName(); // Variable on which projection
parameters are defined as attributes.
+ // Variable on which projection parameters are defined as
attributes.
+ return (variableName != null) ? variableName :
mapping.getName();
});
final Conversion conversion =
decoder.getCoordinateOperationFactory()
.createDefiningConversion(properties, method, parameters);
- /*
- * Projected CRS. The "base to CRS" transform is the conversion
from base directly.
- */
- properties = properties(definition, true,
Convention.PROJECTED_CRS_NAME, () -> {
- return (fromWKT != null) ? fromWKT.getName() :
conversion.getName();
- });
- final ProjectedCRS projected =
crsFactory.createProjectedCRS(properties, baseCRS, conversion,
decoder.getStandardProjectedCS());
- baseToCRS = MathTransforms.concatenate(baseToCRS,
projected.getConversionFromBase().getMathTransform());
- crs = projected;
+ final IdentifiedObject template = (fromWKT != null) ? fromWKT :
conversion;
+ initializeProperties(properties, definition,
Convention.PROJECTED_CRS_NAME, template::getName);
+ properties.put(ProjectedCRS.REMARKS_KEY,
definition.remove(ACDD.comment));
+ properties.put(ProjectedCRS.IDENTIFIERS_KEY, crsIdentifier);
+ crs = crsFactory.createProjectedCRS(properties, baseCRS,
conversion, decoder.getStandardProjectedCS());
}
- /*
- * The CF-Convention said that even if a WKT definition is provided,
other attributes shall be present
- * and have precedence over the WKT definition. Consequently, the
purpose of WKT in netCDF files is not
- * obvious (except for CompoundCRS).
- */
- if (fromWKT != null) verifyCRS(fromWKT);
- setOrVerifyWKT(alreadyParsedWKT, definition, CRS_WKT);
- setOrVerifyWKT(alreadyParsedWKT, definition, SPATIAL_REF);
/*
* Report all projection parameters that have not been used. If the
map is not rendered
* at expected location, it may be because we have ignored some
important parameters.
*/
- definition.remove(CF.LONG_NAME);
if (!definition.isEmpty()) {
- warningInMapping(mapping, null,
Resources.Keys.UnknownProjectionParameters_3,
- String.join(", ", definition.keySet()));
+ String remaining = String.join(", ", definition.keySet());
+ warningInMapping(mapping, null,
Resources.Keys.UnknownProjectionParameters_3, remaining);
+ }
+ /*
+ * Verifies if the CRS parsed from WKT agrees with the CRS built from
the netCDF attributes.
+ * The comparison ignores axis order and units, because those two
characteristics may differ
+ * even if the WKT and the grid mapping attributes are otherwise
consistent.
+ */
+ if (fromWKT != null && verifyCRS(fromWKT)) {
+ // Assume that the WKT has the intended axis order and units.
+ crs = fromWKT;
+ } else if (wantGeographicCRS && gridToCRS != null) {
+ MathTransform swap =
MathTransforms.swapTwoFirstAxes(gridToCRS.getTargetDimensions());
+ gridToCRS = MathTransforms.concatenate(gridToCRS, swap);
+ }
+ /*
+ * Replace by the CRS from authority code if possible, because that
CRS has more metadata.
+ * Again (as for WKT), the axis order and units may differ. If the
CRSs are inconsistent,
+ * the authority code is ignored (grid mapping attributes or WKT have
precedence).
+ *
+ * Note: we do not use `CRS.equivalent(…)` because the WGS84 is often
not recognized as the same,
+ * because of numerous variations in the name and because of datum
versus ensemble differences.
+ */
+ if (crsIdentifier != null) try {
+ final CoordinateReferenceSystem c =
CRS.forCode(crsIdentifier.toString());
+ final MathTransform swap = CRS.findOperation(crs, c,
null).getMathTransform();
+ if (Matrices.isAffine(MathTransforms.getMatrix(swap))) { // See
above note on `CRS.equivalent(…)`.
+ /*
+ * If a "grid to CRS" transform was explicitly specified by
the GDAL's "GeoTransform" attribute,
+ * the output axes may need to be adjusted for matching the
selected CRS axis order and units.
+ * Reminder: `gridToCRS` is two-dimensional, so we need to
work only on horizontal components.
+ */
+ if (gridToCRS != null) {
+ gridToCRS = MathTransforms.concatenate(gridToCRS, swap);
+ }
+ crs = c; // Set only on success.
+ } else {
+ inconsistentCRS(null, 1); // 1 is for "EPSG code".
+ }
+ } catch (OperationNotFoundException e) {
+ inconsistentCRS(e, 1);
+ } catch (FactoryException e) {
+ ignorable(e); // May be because there is no EPSG database, which
is not a netCDF problem.
}
/*
- * Build the "grid to CRS" if present. This is not defined by
CF-convention,
- * but may be present in some non-CF conventions.
+ * If no "grid to CRS" transform was explicitly specified by the
conventions known to this class
+ * (for example, the GDAL's "GeoTransform" attribute), test if a
transform is specified by other
+ * conventions added as a plugin. Example: GCOM-C from the Japan
profile.
*/
if (gridToCRS == null) {
- gridToCRS = decoder.convention().gridToCRS(mapping, baseToCRS);
- // Map pixel corners by `convention().gridToCRS(…)` contract.
- } else {
- gridToCRS = MathTransforms.concatenate(gridToCRS, baseToCRS);
+ if (convention.getClass() != Convention.class) { // Quick check
for common case.
+ final MathTransform geographicToCRS; // From
(latitude, longitude) in degrees to the CRS.
+ geographicToCRS = CRS.findOperation(baseCRS, crs,
null).getMathTransform();
+ if ((gridToCRS = convention.gridToCRS(mapping,
geographicToCRS)) != null) {
+ anchorFromConvention = PixelInCell.CELL_CORNER; // By
`onvention.gridToCRS(…)` contract.
+ }
+ }
}
return true;
}
/**
- * Returns the {@code properties} argument value to give to the factory
methods of geodetic objects.
+ * Sets the {@code properties} argument to the values to give to the
factory methods of geodetic objects.
* The returned map contains at least an entry for {@value
IdentifiedObject#NAME_KEY} with the name
* fetched from the value of the attribute named {@code nameAttribute}.
*
+ * @param properties the map to clear and set. Will be modified
in-place.
* @param definition map containing the attribute values.
- * @param takeComment whether to consume the {@code comment} attribute.
* @param nameAttribute name of the attribute from which to get the name.
* @param nameFallback can return {@link String}, {@link Identifier} or
{@code null}.
*/
- private static Map<String, Object> properties(final Map<String, Object>
definition,
- final boolean
takeComment,
- final String
nameAttribute,
- final Supplier<?>
nameFallback)
+ private static void initializeProperties(final Map<String, Object>
properties,
+ final Map<String, Object>
definition,
+ final String
nameAttribute,
+ final Supplier<?>
nameFallback)
{
+ properties.clear();
Object name = definition.remove(nameAttribute);
if (name == null) {
name = nameFallback.get();
@@ -589,14 +643,7 @@ final class GridMapping {
name = Vocabulary.format(Vocabulary.Keys.Unnamed);
}
}
- if (takeComment) {
- Object comment = definition.remove(ACDD.comment);
- if (comment != null) {
- return Map.of(IdentifiedObject.NAME_KEY, name,
- IdentifiedObject.REMARKS_KEY, comment);
- }
- }
- return Map.of(IdentifiedObject.NAME_KEY, name);
+ properties.put(IdentifiedObject.NAME_KEY, name);
}
/**
@@ -638,17 +685,26 @@ final class GridMapping {
* Verifies that the given <abbr>CRS</abbr> is consistent with the {@link
#crs} attribute.
* If not, a warning will be logger. This method does not change the state
of this object.
*
- * @param fromWKT the object parsed from <abbr>WKT</abbr>.
+ * @param fromWKT the object parsed from <abbr>WKT</abbr>.
+ * @return whether the <abbr>CRS</abbr>s are equivalent.
*/
- private void verifyCRS(final CoordinateReferenceSystem fromWKT) {
- if (!Utilities.deepEquals(crs, fromWKT, ComparisonMode.ALLOW_VARIANT))
{
- warning(mapping, // Node
- null, // Exception
- null, // Resources
- Resources.Keys.InconsistentCRS_2,
- mapping.decoder.getFilename(),
- mapping.getName());
+ private boolean verifyCRS(final CoordinateReferenceSystem fromWKT) {
+ if (Utilities.deepEquals(crs, fromWKT, ComparisonMode.ALLOW_VARIANT)) {
+ return true;
}
+ inconsistentCRS(null, 0); // 0 is for "WKT".
+ return false;
+ }
+
+ /**
+ * Logs a warning with a message saying that a <abbr>CRS</abbr> is
inconsistent with other properties.
+ *
+ * @param cause the exception which has been thrown, or {@code null}.
+ * @param type 0 if the code was parsed from WKT, 1 if from EPSG code.
+ */
+ private void inconsistentCRS(final Exception cause, final int type) {
+ warning(mapping, cause, null, Resources.Keys.InconsistentCRS_3,
+ mapping.decoder.getFilename(), mapping.getName(), type);
}
/**
@@ -659,11 +715,10 @@ final class GridMapping {
* @return whether this method found grid geometry attributes.
*/
private boolean parseGDAL() throws ParseException {
- boolean found =
parseGeoTransform(mapping.getAttributeAsString(GEOTRANSFORM));
+ boolean found =
parseGeoTransform(mapping.getAttributeValue(GEOTRANSFORM));
final String wkt = mapping.getAttributeAsString(SPATIAL_REF);
if (wkt != null) {
- crs = createFromWKT(wkt);
- isWKT = true;
+ crs = createFromWKT(wkt); // Axis order assumed consistent
with `gridToCRS`.
found = true;
}
return found;
@@ -679,26 +734,45 @@ final class GridMapping {
* Y = c[3] + P*c[4] + L*c[5];
* }
*
- * The result is stored in the {@link #gridToCRS} field.
+ * The target <abbr>CRS</abbr> is presumed to have (longitude, latitude)
or (easting, northing) axis order.
+ * Note that this is not necessarily the axis order of the
<abbr>CRS</abbr> in the {@link #crs} field.
+ * Adjustment for change of axis order is done by the code which assigns a
value to {@link #crs}.
+ *
+ * <p>The result is stored in the {@link #gridToCRS} field.
+ * Note that the axis at dimension {@value #SOURCE_AXIS_TO_FLIP} is
flipped compared to usual netCDF direction.
+ * Setting the transform in the expected way depends on the size of the
variable associated with this transform,
+ * and therefore cannot be done here. See {@link #gridToCRS(Variable)}.</p>
*
* @return whether this method found grid geometry attributes.
*/
- @SuppressWarnings("UseSpecificCatch")
- private boolean parseGeoTransform(final String gtr) {
+ private boolean parseGeoTransform(final Object gtr) {
if (gtr != null) {
- final double[] c = parseDoubles(gtr);
- if (c.length == 6) {
- /*
- * GDAL convention maps pixel corners and see the data as if
it was an image.
- * The row which is visually on the top is handled as if its
index was zero,
- * ignoring the fact that this is usually the last row in a
netCDF variable.
- */
- gridToCRS = new AffineTransform2D(c[1], c[4], c[2], c[5],
c[0], c[3]); // X_DIMENSION, Y_DIMENSION
- anchorFromConvention = PixelInCell.CELL_CORNER;
- return true;
+ Exception warning;
+ try {
+ final double[] c;
+ if (gtr instanceof double[]) {
+ c = (double[]) gtr;
+ } else if (gtr instanceof Vector) {
+ c = ((Vector) gtr).doubleValues();
+ } else {
+ c = parseDoubles(gtr.toString());
+ }
+ if (c.length == 6) {
+ /*
+ * GDAL convention maps pixel corners and see the data as
if it was an image.
+ * The row which is visually on the top is handled as if
its index was zero,
+ * ignoring the fact that this is usually the last row in
a netCDF variable.
+ */
+ gridToCRS = new AffineTransform2D(c[1], c[4], c[2], c[5],
c[0], c[3]); // X_DIMENSION, Y_DIMENSION
+ anchorFromConvention = PixelInCell.CELL_CORNER;
+ sourceOfGridToCRS = GEOTRANSFORM;
+ return true;
+ }
+ warning = new
DataStoreContentException(mapping.errors().getString(Errors.Keys.UnexpectedArrayLength_2,
6, c.length));
+ } catch (NumberFormatException e) {
+ warning = e;
}
- var e = new
DataStoreContentException(mapping.errors().getString(Errors.Keys.UnexpectedArrayLength_2,
6, c.length));
- cannotCreateGridOrCRS(mapping, e, true);
+ cannotCreateGridOrCRS(mapping, warning, true);
}
return false;
}
@@ -720,9 +794,9 @@ final class GridMapping {
*/
private boolean parseESRI() throws ParseException, FactoryException {
String code = mapping.getAttributeAsString("ESRI_pe_string");
- isWKT = (code != null);
+ final boolean isWKT = (code != null);
if (code == null) {
- code = mapping.getAttributeAsString("EPSG_code");
+ code = mapping.getAttributeAsString(EPSG_CODE);
if (code == null) {
return false;
}
@@ -752,15 +826,34 @@ final class GridMapping {
final var parsed = (CoordinateReferenceSystem) f.parseObject(wkt);
final Warnings warnings = f.getWarnings();
if (warnings != null) {
- final var record = new LogRecord(Level.WARNING,
warnings.toString());
- record.setLoggerName(Modules.NETCDF);
- record.setSourceClassName(Variable.class.getName());
- record.setSourceMethodName("getGridGeometry");
- mapping.decoder.listeners.warning(record);
+ warning(new LogRecord(Level.WARNING, warnings.toString()));
}
return parsed;
}
+ /**
+ * Logs the given warning, pretending that the warning come form {@link
Variable#getGridGeometry()}.
+ *
+ * @param record the warning to log.
+ */
+ private void warning(final LogRecord record) {
+ record.setLoggerName(Modules.NETCDF);
+ record.setSourceClassName(Variable.class.getName());
+ record.setSourceMethodName("getGridGeometry");
+ mapping.decoder.listeners.warning(record);
+ }
+
+ /**
+ * Logs an exception that can be ignored.
+ *
+ * @param e the exception to ignore.
+ */
+ private void ignorable(final Exception e) {
+ final var record = new LogRecord(Level.FINER, e.getMessage());
+ record.setThrown(e);
+ warning(record);
+ }
+
/**
* Logs a warning with a message saying that we cannot create the grid or
the <abbr>CRS</abbr>.
*
@@ -802,16 +895,22 @@ final class GridMapping {
/**
* Returns the Coordinate Reference System inferred from grid mapping
attribute values, or {@code null} if none.
+ * This method should be invoked for information purposes only, such as
<abbr>ISO</abbr> 19115 metadata,
+ * because the axis order of this <abbr>CRS</abbr> may not agree with the
axis order netCDF variables.
*/
final CoordinateReferenceSystem crs() {
return crs;
}
/**
- * Returns the "grid to CRS", handling the reversal of <var>y</var> axis
direction.
+ * Returns the explicitly declared "grid to CRS" adapted to the specified
variable.
+ * The "grid to CRS" transform can be declared in a {@value #GEOTRANSFORM}
attribute
+ * or computed by a sub-class of {@link Convention} such as
<abbr>GCOM</abbr>.
+ * The adaptation consists in a flip of the <var>y</var> axis direction,
+ * which depends on the shape of the variable.
*
* @param variable the variable for which to obtain the transform.
- * @return the transform for the given variable.
+ * @return the transform for the given variable, or {@code null} if none.
*/
private MathTransform gridToCRS(final Variable variable) {
MathTransform flipped = gridToCRS;
@@ -832,10 +931,12 @@ final class GridMapping {
}
/**
- * Creates a new grid geometry with the extent of the given variable and a
potentially null <abbr>CRS</abbr>.
+ * Creates a new grid geometry with the grid extent of the given netCDF
variable.
+ * Some information such as the <abbr>CRS</abbr> and the "grid to CRS"
transform may be missing.
* This method should be invoked only as a fallback when no existing
{@link GridGeometry} can be used.
- * The CRS and "grid to CRS" transform are null, unless some partial
information was found for example
- * as <abbr>WKT</abbr> string.
+ *
+ * @param variable the variable for which to create a grid geometry.
+ * @return the default (often incomplete) grid geometry for the given
variable.
*/
final GridGeometry createGridCRS(final Variable variable) {
final List<Dimension> dimensions = variable.getGridDimensions();
@@ -845,17 +946,17 @@ final class GridMapping {
final int d = (srcDim - 1) - i; // Convert CRS dimension
to netCDF dimension.
upper[i] = dimensions.get(d).length();
}
- MathTransform implicitG2C = gridToCRS(variable);
- CoordinateReferenceSystem implicitCRS = crs;
- if (implicitG2C != null) {
- int tgtDim = CRS.getDimensionOrZero(implicitCRS);
+ MathTransform explicitG2C = gridToCRS(variable);
+ if (explicitG2C != null) {
+ int tgtDim = CRS.getDimensionOrZero(crs);
if (tgtDim == 0) tgtDim = srcDim;
- MathTransform step1 = changeOfDimension(srcDim,
implicitG2C.getSourceDimensions());
- MathTransform step3 =
changeOfDimension(implicitG2C.getTargetDimensions(), tgtDim);
- implicitG2C = MathTransforms.concatenate(step1, implicitG2C,
step3);
+ MathTransform step1 = changeOfDimension(srcDim,
explicitG2C.getSourceDimensions());
+ MathTransform step3 =
changeOfDimension(explicitG2C.getTargetDimensions(), tgtDim);
+ explicitG2C = MathTransforms.concatenate(step1, explicitG2C,
step3);
}
+ final var anchor = (anchorFromConvention != null) ?
anchorFromConvention : PixelInCell.CELL_CORNER;
final var extent = new GridExtent(null, null, upper, false);
- return new GridGeometry(extent, PixelInCell.CELL_CORNER, implicitG2C,
implicitCRS);
+ return new GridGeometry(extent, anchor, explicitG2C, crs);
}
/**
@@ -883,19 +984,27 @@ final class GridMapping {
*/
final GridGeometry adaptGridCRS(final Variable variable, final
GridGeometry implicit, PixelInCell anchor) {
/*
- * The CRS and grid geometry built from grid mapping attributes are
called "explicit" in this method.
- * This is by contrast with CRS derived from coordinate variables,
which is only implicit.
+ * Some conventions, such as GDAL's "GeoTransform", imply hard-coded
pixel corner/center convention.
+ * If the `gridToCRS` value came from such convention, the associated
"pixel in cell" has precedence.
*/
- CoordinateReferenceSystem explicitCRS = crs;
- MathTransform explicitG2C = gridToCRS(variable);
if (anchorFromConvention != null) {
- // GDAL "GeoTransform" uses pixel corner convention.
anchor = anchorFromConvention;
}
+ final CoordinateReferenceSystem implicitCRS =
+ implicit.isDefined(GridGeometry.CRS) ?
implicit.getCoordinateReferenceSystem() : null;
+ final MathTransform implicitG2C =
+ implicit.isDefined(GridGeometry.GRID_TO_CRS) ?
implicit.getGridToCRS(anchor) : null;
+ /*
+ * The CRS and "grid to CRS" transform built from grid mapping
attributes are called "explicit" below.
+ * This is by contrast with CRS and transform derived from coordinate
variables, which are "implicit".
+ */
+ CoordinateReferenceSystem explicitCRS = crs;
+ MathTransform explicitG2C = gridToCRS(variable);
+ if (explicitG2C == null) {
+ explicitG2C = implicitG2C;
+ }
int firstAffectedCoordinate = 0;
- boolean isSameGrid = true;
- if (implicit.isDefined(GridGeometry.CRS)) {
- final CoordinateReferenceSystem implicitCRS =
implicit.getCoordinateReferenceSystem();
+ if (implicitCRS != null) {
if (explicitCRS == null) {
explicitCRS = implicitCRS;
} else {
@@ -909,90 +1018,123 @@ final class GridMapping {
* If we cannot find which component to replace, assume that
grid mapping describes the first dimensions.
* We have no guarantees that this latter assumption is
correct, but it seems to match common practice.
*/
- Matrix swapAxisOrder = null;
+ CoordinateReferenceSystem match = explicitCRS;
final CoordinateSystem cs = implicitCRS.getCoordinateSystem();
firstAffectedCoordinate = AxisDirections.indexOfColinear(cs,
explicitCRS.getCoordinateSystem());
if (firstAffectedCoordinate < 0) {
- final CoordinateReferenceSystem beforeSwap = explicitCRS;
- explicitCRS =
AbstractCRS.castOrCopy(explicitCRS).forConvention(AxesConvention.RIGHT_HANDED);
- firstAffectedCoordinate =
AxisDirections.indexOfColinear(cs, explicitCRS.getCoordinateSystem());
+ match =
AbstractCRS.castOrCopy(explicitCRS).forConvention(AxesConvention.RIGHT_HANDED);
+ firstAffectedCoordinate =
AxisDirections.indexOfColinear(cs, match.getCoordinateSystem());
if (firstAffectedCoordinate < 0) {
firstAffectedCoordinate = 0;
- if (isWKT && crs != null) {
- explicitCRS = crs; // If specified by
WKT, use the CRS verbatim.
- }
- }
- if (explicitCRS != beforeSwap) try {
- swapAxisOrder =
CoordinateSystems.swapAndScaleAxes(beforeSwap.getCoordinateSystem(),
-
explicitCRS.getCoordinateSystem());
- } catch (IllegalArgumentException |
IncommensurableException e) {
- cannotCreateGridOrCRS(variable, e, false);
- return null;
}
}
/*
- * Replace the grid CRS (or a component of it) by the CRS
parsed from WKT or EPSG code with same (if possible)
- * axis order. If the grid CRS contains more axes (for example
elevation or time axis), we try to keep them.
+ * Replace the implicit CRS (derived from the netCDF
dimensions of the variable) by the explicit CRS
+ * specified by the grid mapping attributes (analyzed by this
class). The axis order and units may be
+ * changed from `implicitCRS` order to `explicitCRS` order.
*/
+ MathTransform swap = null;
try {
+ // Following is useful only if a "grid to CRS" exists.
+ if (explicitG2C != null && match != explicitCRS) {
+ swap = MathTransforms.linear(
+ CoordinateSystems.swapAndScaleAxes(
+ match.getCoordinateSystem(),
+ explicitCRS.getCoordinateSystem()));
+ }
explicitCRS = new CRSMerger(variable.decoder)
.replaceComponent(implicitCRS,
firstAffectedCoordinate, explicitCRS);
- } catch (FactoryException e) {
+ } catch (IllegalArgumentException | IncommensurableException |
FactoryException e) {
cannotCreateGridOrCRS(variable, e, false);
return null;
}
- isSameGrid = implicitCRS.equals(explicitCRS);
- if (isSameGrid) {
- explicitCRS = implicitCRS; // Keep existing
instance if appropriate.
- }
- /*
- * If we have run the
`AbstractCRS.castOrCopy(…).forConvention(…)` code above,
- * the axis order of the CRS may become different than the
axis order which was
- * assumed when the "grid to CRS" transform was built.
- */
- if (swapAxisOrder != null) {
- MathTransform swap = MathTransforms.linear(swapAxisOrder);
- int numTrailingCoordinates =
explicitG2C.getTargetDimensions() - firstAffectedCoordinate;
- swap = MathTransforms.passThrough(firstAffectedCoordinate,
swap, numTrailingCoordinates);
+ if (swap != null) { // Implies `explicitG2C != null`
+ swap = MathTransforms.passThrough(firstAffectedCoordinate,
swap,
+ explicitG2C.getTargetDimensions() -
(firstAffectedCoordinate + swap.getSourceDimensions()));
explicitG2C = MathTransforms.concatenate(explicitG2C,
swap);
}
+ if (implicitCRS.equals(explicitCRS)) {
+ explicitCRS = implicitCRS; // Keep existing
instance.
+ }
}
}
/*
- * Perform the same substitution as above, but in the "grid to CRS"
transform. Note that the "grid to CRS"
- * is usually not specified, so the block performing substitution will
rarely be executed. If executed,
- * then we need to perform selection in target dimensions (not source
dimensions) because the first affected
- * coordinate computed above is in CRS dimension, which is the target
of "grid to CRS" transform.
+ * Perform the same substitution as above, but in the "grid to CRS"
transform instead of in the target CRS.
+ * Note that the "grid to CRS" transform is usually not specified (it
comes from non-standard sources such
+ * as GDAL's "GeoTransform"), so the following will usually be skip
because `(explicitG2C == implicitG2C)`.
+ * If a substitution is executed, it works with target dimensions (not
source dimensions) because the value
+ * of `firstAffectedCoordinate` is in CRS dimension, which is the
target of the "grid to CRS" transform.
*/
- if (implicit.isDefined(GridGeometry.GRID_TO_CRS)) {
- final MathTransform implicitG2C = implicit.getGridToCRS(anchor);
- if (explicitG2C == null) {
- explicitG2C = implicitG2C;
- } else try {
+ if (explicitG2C != implicitG2C && implicitG2C != null) { // Implies
`explicitG2C != null` (see method start).
+ try {
final var sep = new TransformSeparator(implicitG2C,
variable.decoder.getMathTransformFactory());
final int end = firstAffectedCoordinate +
explicitG2C.getTargetDimensions();
explicitG2C = sep.replace(firstAffectedCoordinate, end,
explicitG2C);
- if (explicitG2C != implicitG2C) {
- isSameGrid = false;
- warningInMapping(variable, null,
Resources.Keys.InconsistentTransform_3, GEOTRANSFORM);
- // In current version, GDAL's GeoTransform is the only
supported attribute.
- }
} catch (FactoryException e) {
cannotCreateGridOrCRS(variable, e, true);
return null;
}
+ if (sourceOfGridToCRS != null) {
+ TransformException unexpected = null;
+ boolean consistent;
+ try {
+ consistent =
checkConsistency(MathTransforms.concatenate(implicitG2C,
explicitG2C.inverse()), implicit.getExtent());
+ } catch (TransformException e) {
+ unexpected = e;
+ consistent = false;
+ }
+ if (!consistent) {
+ warningInMapping(variable, unexpected,
Resources.Keys.InconsistentTransform_3, sourceOfGridToCRS);
+ }
+ }
}
/*
* At this point we finished to compute the grid geometry components.
* If any of them have changed, create the new grid geometry.
*/
- if (isSameGrid) {
+ if (implicitCRS == explicitCRS && implicitG2C == explicitG2C) {
return implicit;
} else {
return new GridGeometry(implicit.getExtent(), anchor, explicitG2C,
explicitCRS);
}
}
+ /**
+ * Verifies two "grid to CRS" transforms created in different ways are
approximately equivalent.
+ * The check uses a roundtrip using one transform then the inverse of the
other transform,
+ * and checking whether we get approximately the original values.
+ *
+ * @param roundtrip one transform followed by the inverse of the other
transform.
+ * @param extent range of values to test. A subsampling will be
applied.
+ * @return whether the transforms are equivalent, ignoring rounding errors.
+ */
+ private static boolean checkConsistency(final MathTransform roundtrip,
final GridExtent extent) throws TransformException {
+ if (!roundtrip.isIdentity()) {
+ final int dimension = roundtrip.getSourceDimensions();
+ if (roundtrip.getTargetDimensions() != dimension) {
+ return false;
+ }
+ final int length = dimension * NUM_SAMPLING_POINTS;
+ final double[] source = new double[length];
+ final double[] target = new double[length];
+ for (int c=0; c<dimension; c++) {
+ final double scale = extent.getSize(c, true) /
(NUM_SAMPLING_POINTS - 1);
+ final double offset = extent.getLow(c);
+ for (int i=0; i<NUM_SAMPLING_POINTS; i++) {
+ source[i * dimension + c] = Math.fma(i, scale, offset);
+ }
+ }
+ roundtrip.transform(source, 0, target, 0, NUM_SAMPLING_POINTS);
+ for (int i=0; i<length; i++) {
+ // Condition uses `!` for catching NaN values.
+ if (!(Math.abs(target[i] - source[i]) <=
Numerics.COMPARISON_THRESHOLD)) {
+ return false;
+ }
+ }
+ }
+ return true;
+ }
+
/**
* Returns a string representation for debugging purposes.
*
@@ -1001,8 +1143,8 @@ final class GridMapping {
@Override
public String toString() {
return Strings.toString(getClass(),
- null, mapping.getName(),
- "crs", IdentifiedObjects.getName(crs, null),
- "isWKT", isWKT);
+ null, mapping.getName(),
+ "gridToCRS", sourceOfGridToCRS,
+ "crs", IdentifiedObjects.getName(crs, null));
}
}
diff --git
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources.java
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources.java
index d4556f0320..3a809455ad 100644
---
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources.java
+++
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources.java
@@ -157,10 +157,10 @@ public class Resources extends IndexedResourceBundle {
public static final short IllegalValueRange_4 = 16;
/**
- * The CRS declared by WKT is inconsistent with the attributes of
“{1}” in the “{0}” netCDF
- * file.
+ * The CRS declared by {2,choice,0#WKT|1#EPSG code} is inconsistent
with the attributes of
+ * “{1}” in the “{0}” netCDF file.
*/
- public static final short InconsistentCRS_2 = 29;
+ public static final short InconsistentCRS_3 = 29;
/**
* The “{2}” attribute does not match the transform inferred from the
axes of “{1}” in the
diff --git
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources.properties
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources.properties
index dc03103727..ac66848910 100644
---
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources.properties
+++
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources.properties
@@ -35,7 +35,7 @@ DuplicatedAxis_2 = Duplicated axis
\u201c{1}\u201d in a grid of
DuplicatedAxisType_4 = Axes \u201c{2}\u201d and \u201c{3}\u201d
have the same type \u201c{1}\u201d in netCDF file \u201c{0}\u201d.
IllegalAttributeValue_3 = Illegal value \u201c{2}\u201d for
attribute \u201c{1}\u201d in netCDF file \u201c{0}\u201d.
IllegalValueRange_4 = Illegal value range {2,number} \u2026
{3,number} for variable \u201c{1}\u201d in netCDF file \u201c{0}\u201d.
-InconsistentCRS_2 = The CRS declared by WKT is inconsistent
with the attributes of \u201c{1}\u201d in the \u201c{0}\u201d netCDF file.
+InconsistentCRS_3 = The CRS declared by {2,choice,0#WKT|1#EPSG
code} is inconsistent with the attributes of \u201c{1}\u201d in the
\u201c{0}\u201d netCDF file.
InconsistentTransform_3 = The \u201c{2}\u201d attribute does not
match the transform inferred from the axes of \u201c{1}\u201d in the
\u201c{0}\u201d netCDF file.
GridLongitudeSpanTooWide_2 = The grid spans {0}\u00b0 of longitude,
which may be too wide for the \u201c{1}\u201d domain.
MismatchedAttributeLength_5 = Attributes \u201c{1}\u201d and
\u201c{2}\u201d on variable \u201c{0}\u201d have different lengths: {3} and {4}
respectively.
diff --git
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources_fr.properties
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources_fr.properties
index 24f66085b3..06a31a66bf 100644
---
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources_fr.properties
+++
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/internal/Resources_fr.properties
@@ -40,7 +40,7 @@ DuplicatedAxis_2 = Axe
\u00ab\u202f{1}\u202f\u00bb dupliqu\u00e
DuplicatedAxisType_4 = Les axes \u00ab\u202f{2}\u202f\u00bb et
\u00ab\u202f{3}\u202f\u00bb d\u00e9clarent le m\u00eame type
\u00ab\u202f{1}\u202f\u00bb dans le fichier netCDF \u00ab\u202f{0}\u202f\u00bb.
IllegalAttributeValue_3 = La valeur \u00ab\u202f{2}\u202f\u00bb est
ill\u00e9gale pour l\u2019attribut \u00ab\u202f{1}\u202f\u00bb dans le fichier
netCDF \u00ab\u202f{0}\u202f\u00bb.
IllegalValueRange_4 = Plage de valeurs {2,number} \u2026
{3,number} ill\u00e9gale pour la variable \u00ab\u202f{1}\u202f\u00bb dans le
fichier netCDF \u00ab\u202f{0}\u202f\u00bb.
-InconsistentCRS_2 = Le syst\u00e8me de r\u00e9f\u00e9rence
d\u00e9clar\u00e9 par WKT est incoh\u00e9rent avec les attributs de
\u00ab\u202f{1}\u202f\u00bb dans le fichier netCDF \u00ab\u202f{0}\u202f\u00bb.
+InconsistentCRS_3 = Le syst\u00e8me de r\u00e9f\u00e9rence
d\u00e9clar\u00e9 par {2,choice,0#WKT|1#code EPSG} est incoh\u00e9rent avec les
attributs de \u00ab\u202f{1}\u202f\u00bb dans le fichier netCDF
\u00ab\u202f{0}\u202f\u00bb.
InconsistentTransform_3 = L\u2019attribut
\u00ab\u202f{2}\u202f\u00bb ne correspond pas \u00e0 la transformation
d\u00e9riv\u00e9e des axes de la variable \u00ab\u202f{1}\u202f\u00bb du
fichier netCDF \u00ab\u202f{0}\u202f\u00bb.
GridLongitudeSpanTooWide_2 = La grille s\u2019\u00e9tend sur {0}\u00b0
de longitude, ce qui peut \u00eatre trop pour le domaine de
\u00ab\u202f{1}\u202f\u00bb.
MismatchedAttributeLength_5 = Les attributs \u201c{1}\u201d et
\u201c{2}\u201d de la variable \u201c{0}\u201d ont des longueurs
diff\u00e9rentes\u00a0: {3} et {4} respectivement.
diff --git
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/package-info.java
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/package-info.java
index f7dff8776b..2c8a6f6021 100644
---
a/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/package-info.java
+++
b/endorsed/src/org.apache.sis.storage.netcdf/main/org/apache/sis/storage/netcdf/package-info.java
@@ -40,7 +40,7 @@
* Care must be taken for avoiding confusion when using SIS and UCAR libraries
together.
*
* @author Martin Desruisseaux (IRD, Geomatys)
- * @version 1.7
+ * @version 1.8
* @since 0.3
*/
package org.apache.sis.storage.netcdf;
diff --git
a/endorsed/src/org.apache.sis.storage.netcdf/test/org/apache/sis/storage/netcdf/NetcdfStoreTest.java
b/endorsed/src/org.apache.sis.storage.netcdf/test/org/apache/sis/storage/netcdf/NetcdfStoreTest.java
index 29f56e7cb8..e8c9b7060e 100644
---
a/endorsed/src/org.apache.sis.storage.netcdf/test/org/apache/sis/storage/netcdf/NetcdfStoreTest.java
+++
b/endorsed/src/org.apache.sis.storage.netcdf/test/org/apache/sis/storage/netcdf/NetcdfStoreTest.java
@@ -125,7 +125,9 @@ public final class NetcdfStoreTest extends
DataStoreTestCase {
}
/**
- * Tests the reading of a file a <abbr>CRS</abbr> declared using
<abbr>GDAL</abbr>-specific properties.
+ * Tests the reading of a file a <abbr>CRS</abbr> declared using
<abbr>GDAL</abbr>-specific attributes.
+ * In this test, the "GeoTransform" attribute is inconsistent with the
transform implied by netCDF dimensions
+ * (the coordinate values do not match). Therefore, a warning is expected.
*
* @throws DataStoreException if an error occurred while reading the
netCDF file.
*/
diff --git
a/endorsed/src/org.apache.sis.storage/main/org/apache/sis/storage/esri/WritableStore.java
b/endorsed/src/org.apache.sis.storage/main/org/apache/sis/storage/esri/WritableStore.java
index 66ca6828ea..b74438fc13 100644
---
a/endorsed/src/org.apache.sis.storage/main/org/apache/sis/storage/esri/WritableStore.java
+++
b/endorsed/src/org.apache.sis.storage/main/org/apache/sis/storage/esri/WritableStore.java
@@ -91,7 +91,7 @@ final class WritableStore extends AsciiGridStore implements
WritableGridCoverage
*/
private static double distanceFromIntegers(final MathTransform gridToCRS)
throws TransformException {
final Matrix m = MathTransforms.getMatrix(gridToCRS.inverse());
- if (m != null && Matrices.isAffine(m)) {
+ if (Matrices.isAffine(m)) {
final int last = m.getNumCol() - 1;
double sum = 0;
for (int j=0; j<last; j++) {
diff --git
a/endorsed/src/org.apache.sis.util/main/org/apache/sis/math/Vector.java
b/endorsed/src/org.apache.sis.util/main/org/apache/sis/math/Vector.java
index 9a850a09e5..b19efa9c66 100644
--- a/endorsed/src/org.apache.sis.util/main/org/apache/sis/math/Vector.java
+++ b/endorsed/src/org.apache.sis.util/main/org/apache/sis/math/Vector.java
@@ -1528,6 +1528,8 @@ search: for (;;) {
* to {@link #size()} exclusive. Subclasses may override with more
efficient implementation.</p>
*
* @return a copy of all floating point values in this vector.
+ * @throws NullPointerException if the values are stored as {@code
Object}s and one value is {@code null}.
+ * @throws NumberFormatException if the values are stored as {@code
String}s and one value cannot be parsed.
*
* @see #doubleValue(int)
*/
@@ -1547,6 +1549,8 @@ search: for (;;) {
* to {@link #size()} exclusive. Subclasses may override with more
efficient implementation.</p>
*
* @return a copy of all floating point values in this vector.
+ * @throws NullPointerException if the values are stored as {@code
Object}s and one value is {@code null}.
+ * @throws NumberFormatException if the values are stored as {@code
String}s and one value cannot be parsed.
*
* @see #floatValue(int)
* @see #isSinglePrecision()
