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: Change DatumShiftGrid.replaceOutsideGridCoordinates(…) API for processing all dimensions at once. This is needed when a shift of 360° of longitude (for example) corresponds in a change of grid indices in all dimensions, because the grid is inclined.
Date Thu, 20 Feb 2020 14:50:11 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 e11726f  Change DatumShiftGrid.replaceOutsideGridCoordinates(…) API for processing
all dimensions at once. This is needed when a shift of 360° of longitude (for example) corresponds
in a change of grid indices in all dimensions, because the grid is inclined.
e11726f is described below

commit e11726f8b4148c7a387f013364216afad8792058
Author: Martin Desruisseaux <martin.desruisseaux@geomatys.com>
AuthorDate: Thu Feb 20 15:47:00 2020 +0100

    Change DatumShiftGrid.replaceOutsideGridCoordinates(…) API for processing all dimensions
at once. This is needed when a shift of 360° of longitude (for example) corresponds in a
change of grid indices in all dimensions, because the grid is inclined.
---
 .../provider/DatumShiftGridCompressed.java         |  47 ++++----
 .../referencing/provider/DatumShiftGridFile.java   |  56 ++++-----
 .../sis/referencing/datum/DatumShiftGrid.java      |  96 +++++++++-------
 .../operation/builder/ResidualGrid.java            | 127 +++++++++++++--------
 4 files changed, 183 insertions(+), 143 deletions(-)

diff --git a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/provider/DatumShiftGridCompressed.java
b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/provider/DatumShiftGridCompressed.java
index 139760e..0f6ec07 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/provider/DatumShiftGridCompressed.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/provider/DatumShiftGridCompressed.java
@@ -42,13 +42,7 @@ final class DatumShiftGridCompressed<C extends Quantity<C>, T
extends Quantity<T
     /**
      * Serial number for inter-operability with different versions.
      */
-    private static final long serialVersionUID = 4847888093457104917L;
-
-    /**
-     * Maximal grid index along the <var>y</var> axis, inclusive.
-     * This is the number of grid cells minus 2.
-     */
-    private final int ymax;
+    private static final long serialVersionUID = 1889111858140209014L;
 
     /**
      * An "average" value for the offset in each dimension.
@@ -75,7 +69,6 @@ final class DatumShiftGridCompressed<C extends Quantity<C>, T extends
Quantity<T
             final short[][] data, final double scale)
     {
         super(grid);
-        this.ymax     = getGridSize()[1] - 2;
         this.averages = averages;
         this.data     = data;
         this.scale    = scale;
@@ -181,7 +174,7 @@ final class DatumShiftGridCompressed<C extends Quantity<C>, T
extends Quantity<T
      */
     @Override
     public double getCellValue(final int dim, final int gridX, final int gridY) {
-        return data[dim][gridX + gridY*nx] * scale + averages[dim];
+        return data[dim][gridX + gridY*scanlineStride] * scale + averages[dim];
     }
 
     /**
@@ -190,24 +183,26 @@ final class DatumShiftGridCompressed<C extends Quantity<C>,
T extends Quantity<T
      */
     @Override
     public void interpolateInCell(double gridX, double gridY, double[] vector) {
-        boolean skipX = false;
-        boolean skipY = false;                          // Whether to skip derivative calculation
for X or Y.
-        if (gridX < 0) {gridX = 0; skipX = true;}
-        if (gridY < 0) {gridY = 0; skipY = true;}
-        int ix = (int) gridX;  gridX -= ix;
-        int iy = (int) gridY;  gridY -= iy;
-        if (ix > nx - 2) {
-            skipX |= (ix != nx-1 || gridX != 0);        // Keep value 'false' if gridX ==
gridSize[0] - 1.
-            ix     = nx - 2;
-            gridX  = 1;
-        }
-        if (iy > ymax) {                                // Subtraction of 2 already done
by the constructor.
-            skipY |= (iy != ymax+1 || gridY != 0);      // Keep value 'false' if gridY ==
gridSize[1] - 1.
-            iy     = ymax;
-            gridY  = 1;
+        final int xmax = getGridSize(0) - 2;
+        final int ymax = getGridSize(1) - 2;
+        int ix = (int) gridX;                               // Really want rounding toward
zero (not floor).
+        int iy = (int) gridY;
+        if (ix < 0 || ix > xmax || iy < 0 || iy > ymax) {
+            final double[] gridCoordinates = {gridX, gridY};
+            replaceOutsideGridCoordinates(gridCoordinates);
+            gridX = gridCoordinates[0];
+            gridY = gridCoordinates[1];
+            ix = Math.max(0, Math.min(xmax, (int) gridX));
+            iy = Math.max(0, Math.min(ymax, (int) gridY));
         }
-        final int p00 = nx*iy + ix;
-        final int p10 = nx + p00;
+        gridX -= ix;                                        // If was negative, will continue
to be negative.
+        gridY -= iy;
+        boolean skipX = (gridX < 0); if (skipX) gridX = 0;
+        boolean skipY = (gridY < 0); if (skipY) gridY = 0;
+        if (gridX > 1)  {gridX = 1; skipX = true;}
+        if (gridY > 1)  {gridY = 1; skipY = true;}
+        final int p00 = scanlineStride * iy + ix;
+        final int p10 = scanlineStride + p00;
         final int n   = data.length;
         boolean derivative = (vector.length >= n + INTERPOLATED_DIMENSIONS * INTERPOLATED_DIMENSIONS);
         for (int dim = 0; dim < n; dim++) {
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/provider/DatumShiftGridFile.java
b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/provider/DatumShiftGridFile.java
index 6de48b8..be13d1a 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/provider/DatumShiftGridFile.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/provider/DatumShiftGridFile.java
@@ -77,7 +77,7 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T extends
Quantity<T>>
     /**
      * Serial number for inter-operability with different versions.
      */
-    private static final long serialVersionUID = -4471670781277328193L;
+    private static final long serialVersionUID = -5801692909082130314L;
 
     /**
      * Cache of grids loaded so far. Those grids will be stored by soft references until
the amount of
@@ -111,10 +111,11 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T
extends Quantity<T>>
     private final Path[] files;
 
     /**
-     * Number of grid cells along the <var>x</var> axis.
-     * This is <code>{@linkplain #getGridSize()}[0]</code> as a field for performance
reasons.
+     * Number of cells between the start of adjacent rows in the grid. This is usually {@code
getGridSize(0)},
+     * stored as a field for performance reasons. Value could be greater than {@code getGridSize(0)}
if there
+     * is some elements to ignore at the end of each row.
      */
-    protected final int nx;
+    protected final int scanlineStride;
 
     /**
      * Number of cells that the grid would have if it was spanning 360° of longitude, or
0 if no wraparound
@@ -122,7 +123,7 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T extends
Quantity<T>>
      * an integer number of cells in 360°. This value is used for longitude values that
are on the other side
      * of the ±180° meridian compared to the region where the grid is defined.
      *
-     * @see #replaceOutsideGridCoordinate(int, double)
+     * @see #replaceOutsideGridCoordinates(double[])
      */
     private final double periodX;
 
@@ -182,9 +183,9 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T extends
Quantity<T>>
     {
         super(coordinateUnit, new AffineTransform2D(Δx, 0, 0, Δy, x0, y0).inverse(),
               new int[] {nx, ny}, isCellValueRatio, translationUnit);
-        this.descriptor = descriptor;
-        this.files      = files;
-        this.nx         = nx;
+        this.descriptor     = descriptor;
+        this.files          = files;
+        this.scanlineStride = nx;
         if (Units.isAngular(coordinateUnit)) {
             periodX = Math.rint((Longitude.MAX_VALUE - Longitude.MIN_VALUE) / Math.abs(Δx));
         } else {
@@ -205,12 +206,12 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T
extends Quantity<T>>
      */
     protected DatumShiftGridFile(final DatumShiftGridFile<C,T> other) {
         super(other);
-        descriptor = other.descriptor;
-        files      = other.files;
-        nx         = other.nx;
-        accuracy   = other.accuracy;
-        subgrids   = other.subgrids;
-        periodX    = other.periodX;
+        descriptor     = other.descriptor;
+        files          = other.files;
+        scanlineStride = other.scanlineStride;
+        accuracy       = other.accuracy;
+        subgrids       = other.subgrids;
+        periodX        = other.periodX;
     }
 
     /**
@@ -230,11 +231,11 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T
extends Quantity<T>>
     {
         super(other.getCoordinateUnit(), gridToCRS.inverse(), new int[] {nx, ny},
               other.isCellValueRatio(), other.getTranslationUnit());
-        descriptor = other.descriptor;
-        files      = other.files;
-        this.nx    = nx;
-        periodX    = (other.periodX == 0) ? 0 : Math.rint((Longitude.MAX_VALUE - Longitude.MIN_VALUE)
-                                              / AffineTransforms2D.getScaleX0(gridToCRS));
+        scanlineStride = nx;
+        descriptor     = other.descriptor;
+        files          = other.files;
+        periodX        = (other.periodX == 0) ? 0 : Math.rint((Longitude.MAX_VALUE - Longitude.MIN_VALUE)
+                                                  / AffineTransforms2D.getScaleX0(gridToCRS));
         // Accuracy to be set by caller. Initial value needs to be zero.
     }
 
@@ -395,11 +396,10 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T
extends Quantity<T>>
      * that coordinate inside the range.
      */
     @Override
-    protected double replaceOutsideGridCoordinate(final int dimension, final double gridCoordinate)
{
-        if (dimension == 0 && periodX != 0) {
-            return Math.IEEEremainder(gridCoordinate, periodX);
+    protected void replaceOutsideGridCoordinates(final double[] gridCoordinates) {
+        if (periodX != 0) {
+            gridCoordinates[0] = Math.IEEEremainder(gridCoordinates[0], periodX);
         }
-        return super.replaceOutsideGridCoordinate(dimension, gridCoordinate);
     }
 
     /**
@@ -495,7 +495,7 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T extends
Quantity<T>>
         /**
          * The translation values. {@code offsets.length} is the number of dimensions, and
{@code offsets[dim].length}
          * shall be the same for all {@code dim} value. Component {@code dim} of the translation
vector at coordinate
-         * {@code gridX}, {@code gridY} is {@code offsets[dim][gridX + gridY*nx]}.
+         * {@code gridX}, {@code gridY} is {@code offsets[dim][gridX + gridY*scanlineStride]}.
          */
         final float[][] offsets;
 
@@ -564,13 +564,13 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T
extends Quantity<T>>
          * from an ASCII file, or any other medium that format numbers in base 10.
          *
          * @param  dim    the dimension for which to get an average value.
-         * @param  gridX  the grid index along the <var>x</var> axis, from 0
inclusive to {@link #nx} exclusive.
-         * @param  gridY  the grid index along the <var>y</var> axis, from 0
inclusive to {@code  ny} exclusive.
+         * @param  gridX  the grid index along the <var>x</var> axis, from 0
inclusive to {@code nx} exclusive.
+         * @param  gridY  the grid index along the <var>y</var> axis, from 0
inclusive to {@code ny} exclusive.
          * @return the offset at the given dimension in the grid cell at the given index.
          */
         @Override
         public final double getCellValue(final int dim, final int gridX, final int gridY)
{
-            return DecimalFunctions.floatToDouble(offsets[dim][gridX + gridY*nx]);
+            return DecimalFunctions.floatToDouble(offsets[dim][gridX + gridY*scanlineStride]);
         }
 
         /**
@@ -673,7 +673,7 @@ abstract class DatumShiftGridFile<C extends Quantity<C>, T extends
Quantity<T>>
          */
         @Override
         public final double getCellValue(final int dim, final int gridX, final int gridY)
{
-            return offsets[dim][gridX + gridY*nx];
+            return offsets[dim][gridX + gridY*scanlineStride];
         }
     }
 }
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/referencing/datum/DatumShiftGrid.java
b/core/sis-referencing/src/main/java/org/apache/sis/referencing/datum/DatumShiftGrid.java
index 87f074f..ab08860 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/referencing/datum/DatumShiftGrid.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/referencing/datum/DatumShiftGrid.java
@@ -124,7 +124,7 @@ import org.apache.sis.measure.Units;
  * in which case a longitude of 170° needs to be replaced by −190° before it can be processed
by the grid.
  * The default {@code DatumShiftGrid} class does not apply longitude wraparound automatically
  * (it does not even know which axis, if any, is longitude),
- * but subclasses can add this support by overriding the {@link #replaceOutsideGridCoordinate(int,
double)} method.
+ * but subclasses can add this support by overriding the {@link #replaceOutsideGridCoordinates(double[])}
method.
  *
  * <h2>Sub-grids</h2>
  * Some datum shift grid files provide a grid valid on a wide region, refined with denser
sub-grids in smaller regions.
@@ -283,6 +283,20 @@ public abstract class DatumShiftGrid<C extends Quantity<C>,
T extends Quantity<T
     }
 
     /**
+     * Returns the number of cells in the specified dimension.
+     * Invoking this method is equivalent to {@code getGridSize()[dimension]}.
+     *
+     * @param  dimension  the dimension for which to get the grid size.
+     * @return the number of grid cells in the specified dimension.
+     * @throws IndexOutOfBoundsException if the given index is out of bounds.
+     *
+     * @since 1.1
+     */
+    public int getGridSize(final int dimension) {
+        return gridSize[dimension];
+    }
+
+    /**
      * Returns the domain of validity of input coordinates that can be specified to the
      * {@link #interpolateAt interpolateAt(…)} method. Coordinates outside that domain
      * will still be accepted, but results may be extrapolations far from reality.
@@ -506,13 +520,13 @@ public abstract class DatumShiftGrid<C extends Quantity<C>,
T extends Quantity<T
         final int xmax = gridSize[0] - 2;
         final int ymax = gridSize[1] - 2;
         int ix = (int) gridX;                               // Really want rounding toward
zero (not floor).
-        if (ix < 0 || ix > xmax) {
-            gridX = replaceOutsideGridCoordinate(0, gridX);
-            ix = Math.max(0, Math.min(xmax, (int) gridX));
-        }
         int iy = (int) gridY;
-        if (iy < 0 || iy > ymax) {
-            gridY = replaceOutsideGridCoordinate(1, gridY);
+        if (ix < 0 || ix > xmax || iy < 0 || iy > ymax) {
+            final double[] gridCoordinates = {gridX, gridY};
+            replaceOutsideGridCoordinates(gridCoordinates);
+            gridX = gridCoordinates[0];
+            gridY = gridCoordinates[1];
+            ix = Math.max(0, Math.min(xmax, (int) gridX));
             iy = Math.max(0, Math.min(ymax, (int) gridY));
         }
         gridX -= ix;                                        // If was negative, will continue
to be negative.
@@ -608,13 +622,13 @@ public abstract class DatumShiftGrid<C extends Quantity<C>,
T extends Quantity<T
         final int xmax = gridSize[0] - 2;
         final int ymax = gridSize[1] - 2;
         int ix = (int) gridX;
-        if (ix < 0 || ix > xmax) {
-            gridX = replaceOutsideGridCoordinate(0, gridX);
-            ix = Math.max(0, Math.min(xmax, (int) gridX));
-        }
         int iy = (int) gridY;
-        if (iy < 0 || iy > ymax) {
-            gridY = replaceOutsideGridCoordinate(1, gridY);
+        if (ix < 0 || ix > xmax || iy < 0 || iy > ymax) {
+            final double[] gridCoordinates = {gridX, gridY};
+            replaceOutsideGridCoordinates(gridCoordinates);
+            gridX = gridCoordinates[0];
+            gridY = gridCoordinates[1];
+            ix = Math.max(0, Math.min(xmax, (int) gridX));
             iy = Math.max(0, Math.min(ymax, (int) gridY));
         }
         gridX -= ix;
@@ -745,17 +759,16 @@ public abstract class DatumShiftGrid<C extends Quantity<C>,
T extends Quantity<T
      * @since 1.0
      */
     public boolean isCellInGrid(double gridX, double gridY) {
-        double max = gridSize[1] - 1;
-        if (!(gridY >= 0 && gridY <= max)) {
-            gridY = replaceOutsideGridCoordinate(1, gridY);
-            if (!(gridY >= 0 && gridY <= max)) return false;
-        }
-        max = gridSize[0] - 1;
-        if (!(gridX >= 0 && gridX <= max)) {
-            gridX = replaceOutsideGridCoordinate(0, gridX);
-            if (!(gridX >= 0 && gridX <= max)) return false;
+        final double xmax = gridSize[0] - 1;
+        final double ymax = gridSize[1] - 1;
+        if (gridX >= 0 && gridX <= xmax && gridY >= 0 &&
gridY <= ymax) {
+            return true;
         }
-        return true;
+        final double[] gridCoordinates = {gridX, gridY};
+        replaceOutsideGridCoordinates(gridCoordinates);
+        gridX = gridCoordinates[0];
+        gridY = gridCoordinates[1];
+        return (gridX >= 0 && gridX <= xmax && gridY >= 0 &&
gridY <= ymax);
     }
 
     /**
@@ -765,35 +778,40 @@ public abstract class DatumShiftGrid<C extends Quantity<C>,
T extends Quantity<T
      * longitude range, so a longitude of 170° needs to be converted to a longitude of −190°
before it can be
      * processed by that grid.
      *
-     * <p>The default implementation returns {@code gridCoordinate} unchanged. Subclasses
need to override this
-     * method if they want to handle longitude cycles. Note that the coordinate value is
a grid index, not a
-     * longitude value. So the period to add or remove is the number of cells that the grid
would have if it
-     * was spanning 360° of longitude.</p>
+     * <p>The default implementation does nothing. Subclasses need to override this
method if they want to handle
+     * longitude wraparounds. Note that the coordinate values are grid indices, not longitude
or latitude values.
+     * So the period to add or remove is the number of cells that the grid would have if
it was spanning 360° of
+     * longitude.</p>
      *
      * <div class="note"><b>Example:</b>
-     * this method may be implemented as below:
+     * if longitude values are mapped to {@code gridX} coordinates (in dimension 0), and
if a shift of 360° in
+     * longitude values is equivalent to a shift of {@code periodX} cells in the grid, then
this method can be
+     * implemented as below:
      *
      * {@preformat java
+     *     private final double periodX = ...;      // Number of grid cells in 360° of longitude.
+     *
      *     &#64;Override
-     *     protected double replaceOutsideGridCoordinate(int dimension, double gridCoordinate)
{
-     *         if (dimension == 0) {
-     *             return Math.IEEEremainder(gridCoordinate, periodX);
-     *         } else {
-     *             return super.replaceOutsideGridCoordinate(dimension, gridCoordinate);
-     *         }
+     *     protected void replaceOutsideGridCoordinates(double[] gridCoordinates) {
+     *         gridCoordinates[0] = Math.IEEEremainder(gridCoordinates[0], periodX);
      *     }
      * }</div>
      *
-     * @param  dimension       0 if the given coordinate is {@code gridX}, or 1 if it is
{@code gridY}.
-     * @param  gridCoordinate  the grid coordinate which is outside the range of valid values.
-     * @return the grid coordinate to use. May still be outside the range of valid grid values.
+     * This method receives all grid coordinates in the {@code gridCoordinates} argument
and can modify any
+     * of them, possibly many at once. The reason is because a shift of 360° of longitude
(for example) may
+     * correspond to an offset in both {@code gridX} and {@code gridY} indices if the grid
is inclined.
+     * The {@code gridCoordinates} array length is the number of grid dimensions,
+     * typically {@value #INTERPOLATED_DIMENSIONS}.
+     *
+     * @param  gridCoordinates  on input, the cell indices of the point which is outside
the grid.
+     *         On output, the cell indices of an equivalent point inside the grid if possible.
+     *         Coordinate values are modified in-place.
      *
      * @see Math#IEEEremainder(double, double)
      *
      * @since 1.1
      */
-    protected double replaceOutsideGridCoordinate(final int dimension, final double gridCoordinate)
{
-        return gridCoordinate;
+    protected void replaceOutsideGridCoordinates(double[] gridCoordinates) {
     }
 
     /**
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/builder/ResidualGrid.java
b/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/builder/ResidualGrid.java
index 85ffe2a..2b31c52 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/builder/ResidualGrid.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/builder/ResidualGrid.java
@@ -51,7 +51,7 @@ final class ResidualGrid extends DatumShiftGrid<Dimensionless,Dimensionless>
{
     /**
      * For cross-version compatibility.
      */
-    private static final long serialVersionUID = 5207799661806374259L;
+    private static final long serialVersionUID = -3668228260650927123L;
 
     /**
      * Number of source dimensions of the residual grid.
@@ -104,10 +104,11 @@ final class ResidualGrid extends DatumShiftGrid<Dimensionless,Dimensionless>
{
     }
 
     /**
-     * Number of grid cells along the <var>x</var> axis.
-     * This is <code>{@linkplain #getGridSize()}[0]</code> as a field for performance
reasons.
+     * Number of cells between the start of adjacent rows in the grid. This is usually {@code
getGridSize(0)},
+     * stored as a field for performance reasons. Value could be greater than {@code getGridSize(0)}
if there
+     * is some elements to ignore at the end of each row.
      */
-    private final int nx;
+    private final int scanlineStride;
 
     /**
      * The residual data, as translations to apply on the result of affine transform.
@@ -134,16 +135,27 @@ final class ResidualGrid extends DatumShiftGrid<Dimensionless,Dimensionless>
{
     private final double accuracy;
 
     /**
-     * If grid coordinates in some target dimensions are cyclic, their periods in number
of cells.
-     * Otherwise {@code null}. If non-null, non-cyclic dimensions shall have a period of
zero.
+     * If grid coordinates in some target dimensions are cyclic, the period in number of
cells.
+     * For each scalar value in the {@link LocalizationGridBuilder#periods} array (in units
of
+     * target CRS), the corresponding period in number of cells is a vector. For example
a 360°
+     * shift in longitude does not necessarily correspond to an horizontal or vertical offset
+     * in grid indices; it may be a combination of both if the grid is inclined.
      *
-     * <p>Note that the unit of measurement is different than {@link LocalizationGridBuilder#periods}.
-     * This fields contain values in number of cells, which must have been converted from
values in
-     * degrees (for example). Those numbers may be non-integers.</p>
+     * <p>We should have as many vectors as non-zero values in {@link LocalizationGridBuilder#periods}.
+     * Each {@code periodVector} (in cell units) should be computed from a {@code periods}
vector with
+     * exactly one non-zero value (in CRS units) for allowing shifts in different CRS dimensions
to be
+     * applied independently. Consequently this field should actually be of type {@code double[][]}.
+     * But current version uses only one vector for avoiding the complexity of searching
how to combine
+     * multiple vectors. It is okay for the usual case where only one CRS axis has wraparound
range,
+     * but may need to be revisited in the future.</p>
+     *
+     * <p>This array is {@code null} if no period has been specified, or if a period
has been specified
+     * but we can not convert it from CRS units to a constant number of cells.</p>
      *
      * @see LocalizationGridBuilder#periods
+     * @see #replaceOutsideGridCoordinates(double[])
      */
-    private final double[] cellPeriods;
+    private final double[] periodVector;
 
     /**
      * Creates a new residual grid.
@@ -159,43 +171,34 @@ final class ResidualGrid extends DatumShiftGrid<Dimensionless,Dimensionless>
{
             final double[] periods) throws TransformException
     {
         super(Units.UNITY, sourceToGrid, new int[] {nx, ny}, true, Units.UNITY);
-        this.gridToTarget = gridToTarget;
-        this.offsets      = residuals;
-        this.accuracy     = precision;
-        this.nx           = nx;
-        double[] cellPeriods = null;
+        this.gridToTarget   = gridToTarget;
+        this.offsets        = residuals;
+        this.accuracy       = precision;
+        this.scanlineStride = nx;
+        double[] periodVector = null;
         if (periods != null && gridToTarget.isAffine()) {
             /*
              * We require the transform to be affine because it makes the Jacobian independent
of
-             * coordinate values. It allows us to replace a period in target units by periods
in
-             * grid units without having to take the coordinate values in account.
+             * coordinate values. It allows us to replace a period in target CRS units by
periods
+             * in grid units without having to take the coordinate values in account.
              */
             final MatrixSIS m = MatrixSIS.castOrCopy(gridToTarget.inverse().derivative(null));
-            cellPeriods = m.multiply(periods);
+            periodVector = m.multiply(periods);
             /*
-             * Ideally we would be done. But if we have more than one wraparound dimension
-             * (i.e. if the `periods` array contains more than 1 non-zero value) and if the
-             * `gridToTarget` transform combines the coordinates in many dimensions, we may
-             * have the contributions of many dimensions mixed in a way that we will not
be
-             * able to use. This loop verify that each `cellPeriods` value depends on only
-             * one wraparound dimension of the target CRS. If this is not the case, we have
-             * ambiguity; value will be reset to 0 for disabling wraparound on that dimension.
+             * Above code is not really right. We should verify that `periods` contains exactly
+             * one non-zero element, and if not the case execute above code in a loop with
one
+             * non-zero element by iteration, creating one independent vector each time.
+             * We don't do that for now because `replaceOutsideGridCoordinates(…)` is not
yet
+             * capable to search the best combination of many vectors.
+             *
+             * With current implementation, if the `periods` array contains more than 1 non-zero
value,
+             * then `replaceOutsideGridCoordinates(…)` always shift all wraparound dimensions
together.
+             * For example if the first dimension has a period of 360° and the second dimension
has a
+             * period of 12 months, then `replaceOutsideGridCoordinates(…)` will only shift
by 360°
+             * AND 12 months together, never 360° only or 12 months only.
              */
-            for (int j=0; j<cellPeriods.length; j++) {
-                int contributions = 0;
-                final double cp = Math.abs(cellPeriods[j]);
-                if (cp > 0 && cp >= getGridSize()[j]) {
-                    final double tolerance = Math.ulp(cp);
-                    for (int i=0; i<periods.length; i++) {
-                        if (!(Math.abs(periods[i] * m.getElement(j,i)) <= tolerance))
{     // ! for counting NaNs.
-                            if (++contributions >= 2) break;                         
      // No need to continue.
-                        }
-                    }
-                }
-                cellPeriods[j] = (contributions == 1) ? cp : 0;
-            }
         }
-        this.cellPeriods = cellPeriods;
+        this.periodVector = periodVector;
     }
 
     /**
@@ -240,7 +243,7 @@ final class ResidualGrid extends DatumShiftGrid<Dimensionless,Dimensionless>
{
      */
     @Override
     public double getCellValue(int dim, int gridX, int gridY) {
-        return offsets[(gridX + gridY*nx) * SOURCE_DIMENSION + dim];
+        return offsets[(gridX + gridY*scanlineStride) * SOURCE_DIMENSION + dim];
     }
 
     /**
@@ -249,14 +252,39 @@ final class ResidualGrid extends DatumShiftGrid<Dimensionless,Dimensionless>
{
      * that coordinate inside the range.
      */
     @Override
-    protected double replaceOutsideGridCoordinate(final int dimension, final double gridCoordinate)
{
-        if (cellPeriods != null) {
-            final double p = cellPeriods[dimension];
-            if (p != 0) {
-                return Math.IEEEremainder(gridCoordinate, p);
+    protected void replaceOutsideGridCoordinates(final double[] gridCoordinates) {
+        if (periodVector != null) {
+            /*
+             * We will try to shift the point by an integral multiple of `periodVector`.
+             * Test along each dimension which range of multiples could bring the point
+             * inside the grid, and take the range intersection in all dimensions.
+             */
+            double min = Double.NEGATIVE_INFINITY;
+            double max = Double.POSITIVE_INFINITY;
+            for (int i=0; i<gridCoordinates.length; i++) {
+                final double period = periodVector[i];
+                double toLower = gridCoordinates[i];                // If we subtract a shift
≤ toLower, then coordinate ≥ lower bound (which is 0).
+                double toUpper = toLower - (getGridSize(i) - 1);    // If we subtract a shift
≥ toUpper, then coordinate ≤ upper bound (which is gridSize−1).
+                toLower = Math.floor(toLower / period);             // Convert to integral
number of periods, reverse may be ≤ original value.
+                toUpper = Math.ceil (toUpper / period);             // Convert to integral
number of periods, reverse may be ≥ original value.
+                if (toLower < max) max = toLower;                   // Must shift by no
more than this value otherwise coordinate < lower bound.
+                if (toUpper > min) min = toUpper;                   // Must shift by no
less than this value otherwise coordinate > upper bound.
+            }
+            if (min <= max) {
+                /*
+                 * If at least one multiple exists, take the multiple closest to zero
+                 * (i.e. apply the smallest displacement). Note that the range should
+                 * not include zero, otherwise this point would be inside the grid and
+                 * this method should not have been invoked.
+                 */
+                final double n = (min >= 0) ? min : max;
+                if (Double.isFinite(n)) {
+                    for (int i=0; i<gridCoordinates.length; i++) {
+                        gridCoordinates[i] -= periodVector[i] * n;
+                    }
+                }
             }
         }
-        return super.replaceOutsideGridCoordinate(dimension, gridCoordinate);
     }
 
     /**
@@ -288,8 +316,8 @@ final class ResidualGrid extends DatumShiftGrid<Dimensionless,Dimensionless>
{
         @SuppressWarnings({"CloneInNonCloneableClass", "CloneDoesntCallSuperClone"})
         @Override public Matrix  clone()                            {return this;}
         @Override public boolean isIdentity()                       {return false;}
-        @Override public int     getNumCol()                        {return nx;}
-        @Override public int     getNumRow()                        {return getGridSize()[1];}
+        @Override public int     getNumCol()                        {return getGridSize(0);}
+        @Override public int     getNumRow()                        {return getGridSize(1);}
         @Override public Number  apply     (int[] p)                {return getElement(p[1],
p[0]);}
         @Override public void    setElement(int y, int x, double v) {throw new UnsupportedOperationException();}
 
@@ -305,10 +333,9 @@ final class ResidualGrid extends DatumShiftGrid<Dimensionless,Dimensionless>
{
          * in the table formatted for {@link ParameterDescriptorGroup} string representation.
          */
         @Override public String toString() {
-            final int[] size = getGridSize();
             return new StringBuilder(80).append('[')
                     .append(getElement(0, 0)).append(", …, ")
-                    .append(getElement(size[1] - 1, size[0] - 1))
+                    .append(getElement(getGridSize(1) - 1, getGridSize(0) - 1))
                     .append(']').toString();
         }
 


Mime
View raw message