sis-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From desruisse...@apache.org
Subject [sis] branch geoapi-4.0 updated: NetCDF grid now uses the new Vector methods for calculating the grid geometry.
Date Wed, 21 Nov 2018 23:46:55 GMT
This is an automated email from the ASF dual-hosted git repository.

desruisseaux 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 2e1fb8e  NetCDF grid now uses the new Vector methods for calculating the grid geometry.
2e1fb8e is described below

commit 2e1fb8ecb2afec1e7d9de939a76506279b1bb4b7
Author: Martin Desruisseaux <martin.desruisseaux@geomatys.com>
AuthorDate: Thu Nov 22 00:45:53 2018 +0100

    NetCDF grid now uses the new Vector methods for calculating the grid geometry.
---
 .../java/org/apache/sis/internal/netcdf/Axis.java  |  84 +++++++++++++
 .../java/org/apache/sis/internal/netcdf/Grid.java  | 131 +++++++++++++--------
 .../org/apache/sis/internal/netcdf/Variable.java   |   8 +-
 .../sis/internal/netcdf/ucar/VariableWrapper.java  |   4 +-
 4 files changed, 171 insertions(+), 56 deletions(-)

diff --git a/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Axis.java b/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Axis.java
index 15010fc..a1e533d 100644
--- a/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Axis.java
+++ b/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Axis.java
@@ -27,13 +27,17 @@ import org.opengis.util.FactoryException;
 import org.opengis.referencing.cs.CSFactory;
 import org.opengis.referencing.cs.AxisDirection;
 import org.opengis.referencing.cs.CoordinateSystemAxis;
+import org.opengis.referencing.operation.Matrix;
+import org.opengis.referencing.operation.MathTransform;
 import org.apache.sis.internal.metadata.AxisDirections;
+import org.apache.sis.referencing.operation.transform.MathTransforms;
 import org.apache.sis.referencing.NamedIdentifier;
 import org.apache.sis.metadata.iso.citation.Citations;
 import org.apache.sis.storage.DataStoreException;
 import org.apache.sis.util.iso.Types;
 import org.apache.sis.util.ArraysExt;
 import org.apache.sis.measure.Units;
+import org.apache.sis.math.Vector;
 import ucar.nc2.constants.CDM;
 import ucar.nc2.constants.CF;
 
@@ -301,4 +305,84 @@ public final class Axis extends NamedElement {
         }
         return factory.createCoordinateSystemAxis(properties, abbr, direction, unit);
     }
+
+    /**
+     * Sets the scale and offset coefficients in the given "grid to CRS" transform if possible.
+     * Source and target dimensions used by this method are in "natural" order (reverse of
netCDF order).
+     * Setting the coefficient is possible only if values in this variable are regular,
+     * i.e. the difference between two consecutive values is constant.
+     *
+     * <p>If this method returns {@code true}, then the {@code nonLinears} list is
left unchanged.
+     * If this method returns {@code false}, then a non-linear transform or {@code null}
has been
+     * added to the {@code nonLinears} list.</p>
+     *
+     * @param  gridToCRS   the matrix in which to set scale and offset coefficient.
+     * @param  srcEnd      number of source dimensions (grid dimensions) - 1. Identifies
the last column in the matrix.
+     * @param  tgtDim      the target dimension, which is a dimension of the CRS. Identifies
the matrix row of scale factor.
+     * @param  nonLinears  where to add a non-linear transform if we can not compute a linear
one. {@code null} may be added.
+     * @return whether this method successfully set the scale and offset coefficients.
+     * @throws IOException if an error occurred while reading the data.
+     * @throws DataStoreException if a logical error occurred.
+     */
+    final boolean trySetTransform(final Matrix gridToCRS, final int srcEnd, final int tgtDim,
+            final List<MathTransform> nonLinears) throws IOException, DataStoreException
+    {
+        /*
+         * Normal case where the axis has only one dimension.
+         */
+        if (sourceDimensions.length == 1) {
+            final int srcDim = srcEnd - sourceDimensions[0];
+            if (coordinates.trySetTransform(gridToCRS, srcDim, tgtDim, null)) {
+                return true;
+            } else {
+                nonLinears.add(MathTransforms.interpolate(null, coordinates.read().doubleValues()));
+                return false;
+            }
+        }
+        /*
+         * In netCDF files, axes are sometime associated to two-dimensional localization
grids.
+         * If this is the case, then the following block checks if we can reduce those grids
to
+         * one-dimensional vector. For example the following localisation grids:
+         *
+         *    10 10 10 10                  10 12 15 20
+         *    12 12 12 12        or        10 12 15 20
+         *    15 15 15 15                  10 12 15 20
+         *    20 20 20 20                  10 12 15 20
+         *
+         * can be reduced to a one-dimensional {10 12 15 20} vector (orientation matter however).
+         *
+         * Note: following block is currently restricted to the two-dimensional case, but
it could
+         * be generalized to n-dimensional case if we resolve the default case in the switch
statement.
+         */
+        if (sourceDimensions.length == 2) {
+            Vector data = coordinates.read();
+            final int[] repetitions = data.repetitions();           // Detects repetitions
as illustrated above.
+            for (int i=0; i<sourceDimensions.length; i++) {
+                final int srcDim = srcEnd - sourceDimensions[i];    // "Natural" order is
reverse of netCDF order.
+                final int length = sourceSizes[i];
+                int step = 1;
+                for (int j=0; j<sourceDimensions.length; j++) {
+                    int previous = srcEnd - sourceDimensions[j];
+                    if (previous < srcDim) step *= sourceSizes[j];
+                }
+                final boolean condition;
+                switch (srcDim) {
+                    case 0:  condition = repetitions.length > 1 && (repetitions[1]
% length) == 0; break;
+                    case 1:  condition = repetitions.length > 0 && (repetitions[0]
% step)   == 0; break;
+                    default: throw new AssertionError();        // I don't know yet how to
generalize to n dimensions.
+                }
+                if (condition) {                                // Repetition length shall
be grid size (or a multiple).
+                    data = data.subSampling(0, step, length);
+                    if (coordinates.trySetTransform(gridToCRS, srcDim, tgtDim, data)) {
+                        return true;
+                    } else {
+                        nonLinears.add(MathTransforms.interpolate(null, data.doubleValues()));
+                        return false;
+                    }
+                }
+            }
+        }
+        nonLinears.add(null);
+        return false;
+    }
 }
diff --git a/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Grid.java b/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Grid.java
index af17ba7..e8a4c50 100644
--- a/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Grid.java
+++ b/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Grid.java
@@ -25,13 +25,13 @@ import org.opengis.referencing.datum.PixelInCell;
 import org.opengis.referencing.operation.Matrix;
 import org.opengis.referencing.operation.MathTransform;
 import org.opengis.referencing.operation.TransformException;
+import org.opengis.referencing.operation.MathTransformFactory;
 import org.opengis.referencing.cs.CoordinateSystem;
 import org.opengis.referencing.crs.SingleCRS;
 import org.opengis.referencing.crs.CoordinateReferenceSystem;
 import org.opengis.metadata.spatial.DimensionNameType;
 import org.apache.sis.internal.metadata.AxisDirections;
 import org.apache.sis.referencing.operation.matrix.Matrices;
-import org.apache.sis.referencing.operation.transform.MathTransforms;
 import org.apache.sis.coverage.grid.GridExtent;
 import org.apache.sis.coverage.grid.GridGeometry;
 import org.apache.sis.storage.DataStoreException;
@@ -196,6 +196,43 @@ public abstract class Grid extends NamedElement {
     }
 
     /**
+     * Builds the grid extent if the shape is available. The shape may not be available
+     * if a dimension has unlimited length. The dimension names are informative only.
+     *
+     * @param  axes  value of {@link #getAxes()}. Should be the same as {@link #axes}.
+     */
+    @SuppressWarnings("fallthrough")
+    private GridExtent getExtent(final Axis[] axes) {
+        final long[] high = getShape();
+        if (high == null) {
+            return null;
+        }
+        final DimensionNameType[] names = new DimensionNameType[high.length];
+        switch (names.length) {
+            default: names[1] = DimensionNameType.ROW;      // Fall through
+            case 1:  names[0] = DimensionNameType.COLUMN;   // Fall through
+            case 0:  break;
+        }
+        for (final Axis axis : axes) {
+            if (axis.sourceDimensions.length == 1) {
+                final DimensionNameType name;
+                if (AxisDirections.isVertical(axis.direction)) {
+                    name = DimensionNameType.VERTICAL;
+                } else if (AxisDirections.isTemporal(axis.direction)) {
+                    name = DimensionNameType.TIME;
+                } else {
+                    continue;
+                }
+                final int dim = axis.sourceDimensions[0];
+                if (dim >= 0 && dim < names.length) {
+                    names[names.length - 1 - dim] = name;
+                }
+            }
+        }
+        return new GridExtent(names, new long[high.length], high, false);
+    }
+
+    /**
      * Returns an object containing the grid size, the CRS and the conversion from grid indices
to CRS coordinates.
      * This is the public object exposed to users.
      *
@@ -204,82 +241,72 @@ public abstract class Grid extends NamedElement {
      * @throws  IOException if an I/O operation was necessary but failed.
      * @throws  DataStoreException if the CRS can not be constructed.
      */
-    @SuppressWarnings("fallthrough")
     public final GridGeometry getGridGeometry(final Decoder decoder) throws IOException,
DataStoreException {
         if (!isGeometryDetermined) try {
             isGeometryDetermined = true;        // Set now for avoiding new attempts if creation
fail.
             final Axis[] axes = getAxes();      // In netCDF order (reverse of "natural"
order).
             /*
-             * Build the grid extent if the shape is available. The shape may not be available
-             * if a dimension has unlimited length. The dimension names are informative only.
-             */
-            final GridExtent extent;
-            final long[] high = getShape();
-            if (high != null) {
-                final DimensionNameType[] names = new DimensionNameType[high.length];
-                switch (names.length) {
-                    default: names[1] = DimensionNameType.ROW;      // Fall through
-                    case 1:  names[0] = DimensionNameType.COLUMN;   // Fall through
-                    case 0:  break;
-                }
-                for (final Axis axis : axes) {
-                    if (axis.sourceDimensions.length == 1) {
-                        final DimensionNameType name;
-                        if (AxisDirections.isVertical(axis.direction)) {
-                            name = DimensionNameType.VERTICAL;
-                        } else if (AxisDirections.isTemporal(axis.direction)) {
-                            name = DimensionNameType.TIME;
-                        } else {
-                            continue;
-                        }
-                        final int dim = axis.sourceDimensions[0];
-                        if (dim >= 0 && dim < names.length) {
-                            names[names.length - 1 - dim] = name;
-                        }
-                    }
-                }
-                extent = new GridExtent(names, new long[high.length], high, false);
-            } else {
-                extent = null;
-            }
-            /*
              * Creates the "grid to CRS" transform. The number of columns is the number of
dimensions in the grid
              * (the source) +1, and the number of rows is the number of dimensions in the
CRS (the target) +1.
              * The order of dimensions in the transform is the reverse of the netCDF axis
order.
              */
             int srcEnd = getSourceDimensions();
             int tgtEnd = getTargetDimensions();
+            final int[] deferred = new int[axes.length];
+            final List<MathTransform> nonLinears = new ArrayList<>();
             final Matrix affine = Matrices.createZero(tgtEnd + 1, srcEnd + 1);
             affine.setElement(tgtEnd--, srcEnd--, 1);
             for (int i=axes.length; --i >= 0;) {
-                final Axis axis = axes[i];
-                final int tgtDim = tgtEnd - i;
-                if (axis.sourceDimensions.length == 1) {
-                    final int srcDim = srcEnd - axis.sourceDimensions[0];
-                    if (axis.coordinates.trySetTransform(affine, srcDim, tgtDim)) {
-                        continue;
-                    }
+                if (!axes[i].trySetTransform(affine, srcEnd, tgtEnd - i, nonLinears)) {
+                    deferred[nonLinears.size() - 1] = i;
                 }
-                /*
-                 * If we have not been able to set coefficients in the matrix (because the
transform is non-linear),
-                 * set a single scale factor to 1 in the matrix row; we will concatenate
later with a non-linear transform.
-                 * The coefficient that we set to 1 is the one for the source dimension which
is not already taken by another row.
-                 */
-findFree:       for (int srcDim : axis.sourceDimensions) {
+            }
+            /*
+             * If we have not been able to set coefficients in the matrix (because the transform
is non-linear),
+             * set a single scale factor to 1 in the matrix row; we will concatenate later
with a non-linear transform.
+             * The coefficient that we set to 1 is the one for the source dimension which
is not already taken by another row.
+             */
+            final MathTransformFactory factory = decoder.getMathTransformFactory();
+            for (int i=0; i<nonLinears.size(); i++) {
+                final int d = deferred[i];
+findFree:       for (int srcDim : axes[d].sourceDimensions) {
                     srcDim = srcEnd - srcDim;
                     for (int j=affine.getNumRow(); --j>=0;) {
                         if (affine.getElement(j, srcDim) != 0) {
                             continue findFree;
                         }
                     }
+                    final int tgtDim = tgtEnd - d;
                     affine.setElement(tgtDim, srcDim, 1);
-                    // TODO: prepare non-linear transform here for later concatenation.
+                    /*
+                     * Prepare non-linear transform here for later concatenation.
+                     * TODO: what to do if the transform is null?
+                     */
+                    MathTransform tr = nonLinears.get(i);
+                    if (tr != null) {
+                        tr = factory.createPassThroughTransform(srcDim, tr, (srcEnd + 1)
- (srcDim + tr.getSourceDimensions()));
+                        nonLinears.set(i, tr);
+                    }
                     break;
                 }
             }
-            MathTransform gridToCRS = MathTransforms.linear(affine);
-            geometry = new GridGeometry(extent, PixelInCell.CELL_CENTER, gridToCRS, getCoordinateReferenceSystem(decoder));
-        } catch (TransformException ex) {
+            /*
+             * Final transform, as the concatenation of the non-linear transforms followed
by the affine transform.
+             * We concatenate the affine transform last because it may change axis order.
+             */
+            MathTransform gridToCRS = null;
+            nonLinears.add(factory.createAffineTransform(affine));
+            for (final MathTransform tr : nonLinears) {
+                if (tr != null) {
+                    if (gridToCRS == null) {
+                        gridToCRS = tr;
+                    } else {
+                        gridToCRS = factory.createConcatenatedTransform(gridToCRS, tr);
+                    }
+                }
+            }
+            geometry = new GridGeometry(getExtent(axes), PixelInCell.CELL_CENTER, gridToCRS,
getCoordinateReferenceSystem(decoder));
+        } catch (FactoryException | TransformException ex) {
             canNotCreate(decoder, "getGridGeometry", Resources.Keys.CanNotCreateGridGeometry_3,
ex);
         }
         return geometry;
diff --git a/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Variable.java
b/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Variable.java
index e7cc8dc..412d345 100644
--- a/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Variable.java
+++ b/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/Variable.java
@@ -363,6 +363,7 @@ public abstract class Variable extends NamedElement {
 
     /**
      * Sets the scale and offset coefficients in the given "grid to CRS" transform if possible.
+     * Source and target dimensions given to this method are in "natural" order (reverse
of netCDF order).
      * This method is invoked only for variables that represent a coordinate system axis.
      * Setting the coefficient is possible only if values in this variable are regular,
      * i.e. the difference between two consecutive values is constant.
@@ -370,14 +371,17 @@ public abstract class Variable extends NamedElement {
      * @param  gridToCRS  the matrix in which to set scale and offset coefficient.
      * @param  srcDim     the source dimension, which is a dimension of the grid. Identifies
the matrix column of scale factor.
      * @param  tgtDim     the target dimension, which is a dimension of the CRS.  Identifies
the matrix row of scale factor.
+     * @param  values     the vector to use for computing scale and offset, or {@code null}
if it has not been read yet.
      * @return whether this method successfully set the scale and offset coefficients.
      * @throws IOException if an error occurred while reading the data.
      * @throws DataStoreException if a logical error occurred.
      */
-    protected boolean trySetTransform(final Matrix gridToCRS, final int srcDim, final int
tgtDim)
+    protected boolean trySetTransform(final Matrix gridToCRS, final int srcDim, final int
tgtDim, Vector values)
             throws IOException, DataStoreException
     {
-        final Vector values = read();
+        if (values == null) {
+            values = read();
+        }
         final int n = values.size() - 1;
         if (n >= 0) {
             final double first = values.doubleValue(0);
diff --git a/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/ucar/VariableWrapper.java
b/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/ucar/VariableWrapper.java
index 10bde94..6bc71dc 100644
--- a/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/ucar/VariableWrapper.java
+++ b/storage/sis-netcdf/src/main/java/org/apache/sis/internal/netcdf/ucar/VariableWrapper.java
@@ -357,7 +357,7 @@ final class VariableWrapper extends Variable {
      * This method is invoked only for variables that represent a coordinate system axis.
      */
     @Override
-    protected boolean trySetTransform(final Matrix gridToCRS, final int srcDim, final int
tgtDim)
+    protected boolean trySetTransform(final Matrix gridToCRS, final int srcDim, final int
tgtDim, final Vector values)
             throws IOException, DataStoreException
     {
         if (variable instanceof CoordinateAxis1D) {
@@ -368,7 +368,7 @@ final class VariableWrapper extends Variable {
                 return true;
             }
         }
-        return super.trySetTransform(gridToCRS, srcDim, tgtDim);
+        return super.trySetTransform(gridToCRS, srcDim, tgtDim, values);
     }
 
     /**


Mime
View raw message