You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@sis.apache.org by de...@apache.org on 2020/10/05 11:04:29 UTC
[sis] 04/04: Apply on transformations of `Rectangle2D` objects the
same wraparound checks than we did in previous commit for `Envelope`
objects.
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
commit 579bc7e59e3e6236ffc0112b6474e6a41d79de32
Author: Martin Desruisseaux <ma...@geomatys.com>
AuthorDate: Mon Oct 5 12:41:07 2020 +0200
Apply on transformations of `Rectangle2D` objects the same wraparound checks than we did in previous commit for `Envelope` objects.
---
.../java/org/apache/sis/geometry/Shapes2D.java | 291 +++++++++++----------
.../apache/sis/geometry/WraparoundInEnvelope.java | 15 +-
2 files changed, 161 insertions(+), 145 deletions(-)
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/geometry/Shapes2D.java b/core/sis-referencing/src/main/java/org/apache/sis/geometry/Shapes2D.java
index c037999..7b565a7 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/geometry/Shapes2D.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/geometry/Shapes2D.java
@@ -49,7 +49,7 @@ import org.apache.sis.util.Static;
*
* @author Martin Desruisseaux (IRD, Geomatys)
* @author Johann Sorel (Geomatys)
- * @version 1.0
+ * @version 1.1
* @since 0.8
* @module
*/
@@ -200,159 +200,162 @@ public final class Shapes2D extends Static {
double ymin = Double.POSITIVE_INFINITY;
double xmax = Double.NEGATIVE_INFINITY;
double ymax = Double.NEGATIVE_INFINITY;
- /*
- * Notation (as if we were applying a map projection, but this is not necessarily the case):
- * - (λ,φ) are coordinate values before projection.
- * - (x,y) are coordinate values after projection.
- * - D[00|01|10|11] are the ∂x/∂λ, ∂x/∂φ, ∂y/∂λ and ∂y/∂φ derivatives respectively.
- * - Variables with indice 0 are for the very first point in iteration order.
- * - Variables with indice 1 are for the values of the previous iteration.
- * - Variables with indice 2 are for the current values in the iteration.
- * - P1-P2 form a line segment to be checked for curvature.
- */
- double x0=0, y0=0, λ0=0, φ0=0;
- double x1=0, y1=0, λ1=0, φ1=0;
- Matrix D0=null, D1=null, D2=null;
- // x2 and y2 defined inside the loop.
- boolean isDerivativeSupported = true;
- final CurveExtremum extremum = new CurveExtremum();
- for (int i=0; i<=8; i++) {
+ final WraparoundInEnvelope.Controller wc = new WraparoundInEnvelope.Controller(transform);
+ do {
/*
- * Iteration order (center must be last):
- *
- * (6)────(5)────(4)
- * | |
- * (7) (8) (3)
- * | |
- * (0)────(1)────(2)
+ * Notation (as if we were applying a map projection, but this is not necessarily the case):
+ * - (λ,φ) are coordinate values before projection.
+ * - (x,y) are coordinate values after projection.
+ * - D[00|01|10|11] are the ∂x/∂λ, ∂x/∂φ, ∂y/∂λ and ∂y/∂φ derivatives respectively.
+ * - Variables with indice 0 are for the very first point in iteration order.
+ * - Variables with indice 1 are for the values of the previous iteration.
+ * - Variables with indice 2 are for the current values in the iteration.
+ * - P1-P2 form a line segment to be checked for curvature.
*/
- double λ2, φ2;
- switch (i) {
- case 0: case 6: case 7: λ2 = envelope.getMinX(); break;
- case 1: case 5: case 8: λ2 = envelope.getCenterX(); break;
- case 2: case 3: case 4: λ2 = envelope.getMaxX(); break;
- default: throw new AssertionError(i);
- }
- switch (i) {
- case 0: case 1: case 2: φ2 = envelope.getMinY(); break;
- case 3: case 7: case 8: φ2 = envelope.getCenterY(); break;
- case 4: case 5: case 6: φ2 = envelope.getMaxY(); break;
- default: throw new AssertionError(i);
- }
- point[0] = λ2;
- point[1] = φ2;
- try {
- D1 = D2;
- D2 = Envelopes.derivativeAndTransform(transform, point, point, 0, isDerivativeSupported && i != 8);
- } catch (TransformException e) {
- if (!isDerivativeSupported) {
- throw e; // Derivative were already disabled, so something went wrong.
+ double x0=0, y0=0, λ0=0, φ0=0;
+ double x1=0, y1=0, λ1=0, φ1=0;
+ Matrix D0=null, D1=null, D2=null;
+ // x2 and y2 defined inside the loop.
+ boolean isDerivativeSupported = true;
+ final CurveExtremum extremum = new CurveExtremum();
+ for (int i=0; i<=8; i++) {
+ /*
+ * Iteration order (center must be last):
+ *
+ * (6)────(5)────(4)
+ * | |
+ * (7) (8) (3)
+ * | |
+ * (0)────(1)────(2)
+ */
+ double λ2, φ2;
+ switch (i) {
+ case 0: case 6: case 7: λ2 = envelope.getMinX(); break;
+ case 1: case 5: case 8: λ2 = envelope.getCenterX(); break;
+ case 2: case 3: case 4: λ2 = envelope.getMaxX(); break;
+ default: throw new AssertionError(i);
+ }
+ switch (i) {
+ case 0: case 1: case 2: φ2 = envelope.getMinY(); break;
+ case 3: case 7: case 8: φ2 = envelope.getCenterY(); break;
+ case 4: case 5: case 6: φ2 = envelope.getMaxY(); break;
+ default: throw new AssertionError(i);
}
- isDerivativeSupported = false; D2 = null;
point[0] = λ2;
point[1] = φ2;
- transform.transform(point, 0, point, 0, 1);
- Envelopes.recoverableException(Shapes2D.class, e); // Log only if the above call was successful.
- }
- double x2 = point[0];
- double y2 = point[1];
- if (x2 < xmin) xmin = x2;
- if (x2 > xmax) xmax = x2;
- if (y2 < ymin) ymin = y2;
- if (y2 > ymax) ymax = y2;
- switch (i) {
- case 0: { // Remember the first point.
- λ0=λ2; x0=x2;
- φ0=φ2; y0=y2;
- D0=D2;
- break;
- }
- case 8: { // Close the iteration with the first point.
- λ2=λ0; x2=x0; // Discard P2 because it is the rectangle center.
- φ2=φ0; y2=y0;
- D2=D0;
- break;
+ try {
+ D1 = D2;
+ D2 = Envelopes.derivativeAndTransform(wc.transform, point, point, 0, isDerivativeSupported && i != 8);
+ } catch (TransformException e) {
+ if (!isDerivativeSupported) {
+ throw e; // Derivative were already disabled, so something went wrong.
+ }
+ isDerivativeSupported = false; D2 = null;
+ point[0] = λ2;
+ point[1] = φ2;
+ wc.transform.transform(point, 0, point, 0, 1);
+ Envelopes.recoverableException(Shapes2D.class, e); // Log only if the above call was successful.
}
- }
- /*
- * At this point, we expanded the rectangle using the projected points. Now try
- * to use the information provided by derivatives at those points, if available.
- * For the following block, notation is:
- *
- * - s are coordinate values in the source space (λ or φ)
- * - t are coordinate values in the target space (x or y)
- *
- * They are not necessarily in the same dimension. For example would could have
- * s=λ while t=y. This is typically the case when inspecting the top or bottom
- * line segment of the rectangle.
- *
- * The same technic is also applied in the transform(MathTransform, Envelope) method.
- * The general method is more "elegant", at the cost of more storage requirement.
- */
- if (D1 != null && D2 != null) {
- final int srcDim;
- final double s1, s2; // Coordinate values in source space (before projection)
+ double x2 = point[0];
+ double y2 = point[1];
+ if (x2 < xmin) xmin = x2;
+ if (x2 > xmax) xmax = x2;
+ if (y2 < ymin) ymin = y2;
+ if (y2 > ymax) ymax = y2;
switch (i) {
- case 1: case 2: case 5: case 6: {assert φ2==φ1; srcDim=0; s1=λ1; s2=λ2; break;} // Horizontal segment
- case 3: case 4: case 7: case 8: {assert λ2==λ1; srcDim=1; s1=φ1; s2=φ2; break;} // Vertical segment
- default: throw new AssertionError(i);
+ case 0: { // Remember the first point.
+ λ0=λ2; x0=x2;
+ φ0=φ2; y0=y2;
+ D0=D2;
+ break;
+ }
+ case 8: { // Close the iteration with the first point.
+ λ2=λ0; x2=x0; // Discard P2 because it is the rectangle center.
+ φ2=φ0; y2=y0;
+ D2=D0;
+ break;
+ }
}
- final double min, max;
- if (s1 < s2) {min=s1; max=s2;}
- else {min=s2; max=s1;}
- int tgtDim = 0;
- do { // Executed exactly twice, for dimensions 0 and 1 in the projected space.
- extremum.resolve(s1, (tgtDim == 0) ? x1 : y1, D1.getElement(tgtDim, srcDim),
- s2, (tgtDim == 0) ? x2 : y2, D2.getElement(tgtDim, srcDim));
- /*
- * At this point we found the extremum of the projected line segment
- * using a cubic curve t = A + Bs + Cs² + Ds³ approximation. Before
- * to add those extremum into the projected bounding box, we need to
- * ensure that the source coordinate is inside the the original
- * (unprojected) bounding box.
- */
- boolean isP2 = false;
- do { // Executed exactly twice, one for each point.
- final double se = isP2 ? extremum.ex2 : extremum.ex1;
- if (se > min && se < max) {
- final double te = isP2 ? extremum.ey2 : extremum.ey1;
- if ((tgtDim == 0) ? (te < xmin || te > xmax) : (te < ymin || te > ymax)) {
- /*
- * At this point, we have determined that adding the extremum point
- * to the rectangle would have expanded it. However we will not add
- * that point directly, because maybe its position is not quite right
- * (since we used a cubic curve approximation). Instead, we project
- * the point on the rectangle border which is located vis-à-vis the
- * extremum. Our tests show that the correction can be as much as 50
- * metres.
- */
- final double oldX = point[0];
- final double oldY = point[1];
- if (srcDim == 0) {
- point[0] = se;
- point[1] = φ1; // == φ2 since we have an horizontal segment.
- } else {
- point[0] = λ1; // == λ2 since we have a vertical segment.
- point[1] = se;
+ /*
+ * At this point, we expanded the rectangle using the projected points. Now try
+ * to use the information provided by derivatives at those points, if available.
+ * For the following block, notation is:
+ *
+ * - s are coordinate values in the source space (λ or φ)
+ * - t are coordinate values in the target space (x or y)
+ *
+ * They are not necessarily in the same dimension. For example would could have
+ * s=λ while t=y. This is typically the case when inspecting the top or bottom
+ * line segment of the rectangle.
+ *
+ * The same technic is also applied in the transform(MathTransform, Envelope) method.
+ * The general method is more "elegant", at the cost of more storage requirement.
+ */
+ if (D1 != null && D2 != null) {
+ final int srcDim;
+ final double s1, s2; // Coordinate values in source space (before projection)
+ switch (i) {
+ case 1: case 2: case 5: case 6: {assert φ2==φ1; srcDim=0; s1=λ1; s2=λ2; break;} // Horizontal segment
+ case 3: case 4: case 7: case 8: {assert λ2==λ1; srcDim=1; s1=φ1; s2=φ2; break;} // Vertical segment
+ default: throw new AssertionError(i);
+ }
+ final double min, max;
+ if (s1 < s2) {min=s1; max=s2;}
+ else {min=s2; max=s1;}
+ int tgtDim = 0;
+ do { // Executed exactly twice, for dimensions 0 and 1 in the projected space.
+ extremum.resolve(s1, (tgtDim == 0) ? x1 : y1, D1.getElement(tgtDim, srcDim),
+ s2, (tgtDim == 0) ? x2 : y2, D2.getElement(tgtDim, srcDim));
+ /*
+ * At this point we found the extremum of the projected line segment
+ * using a cubic curve t = A + Bs + Cs² + Ds³ approximation. Before
+ * to add those extremum into the projected bounding box, we need to
+ * ensure that the source coordinate is inside the the original
+ * (unprojected) bounding box.
+ */
+ boolean isP2 = false;
+ do { // Executed exactly twice, one for each point.
+ final double se = isP2 ? extremum.ex2 : extremum.ex1;
+ if (se > min && se < max) {
+ final double te = isP2 ? extremum.ey2 : extremum.ey1;
+ if ((tgtDim == 0) ? (te < xmin || te > xmax) : (te < ymin || te > ymax)) {
+ /*
+ * At this point, we have determined that adding the extremum point
+ * to the rectangle would have expanded it. However we will not add
+ * that point directly, because maybe its position is not quite right
+ * (since we used a cubic curve approximation). Instead, we project
+ * the point on the rectangle border which is located vis-à-vis the
+ * extremum. Our tests show that the correction can be as much as 50
+ * metres.
+ */
+ final double oldX = point[0];
+ final double oldY = point[1];
+ if (srcDim == 0) {
+ point[0] = se;
+ point[1] = φ1; // == φ2 since we have an horizontal segment.
+ } else {
+ point[0] = λ1; // == λ2 since we have a vertical segment.
+ point[1] = se;
+ }
+ wc.transform.transform(point, 0, point, 0, 1);
+ final double x = point[0];
+ final double y = point[1];
+ if (x < xmin) xmin = x;
+ if (x > xmax) xmax = x;
+ if (y < ymin) ymin = y;
+ if (y > ymax) ymax = y;
+ point[0] = oldX;
+ point[1] = oldY;
}
- transform.transform(point, 0, point, 0, 1);
- final double x = point[0];
- final double y = point[1];
- if (x < xmin) xmin = x;
- if (x > xmax) xmax = x;
- if (y < ymin) ymin = y;
- if (y > ymax) ymax = y;
- point[0] = oldX;
- point[1] = oldY;
}
- }
- } while ((isP2 = !isP2) == true);
- } while (++tgtDim == 1);
+ } while ((isP2 = !isP2) == true);
+ } while (++tgtDim == 1);
+ }
+ λ1=λ2; x1=x2;
+ φ1=φ2; y1=y2;
+ D1=D2;
}
- λ1=λ2; x1=x2;
- φ1=φ2; y1=y2;
- D1=D2;
- }
+ } while (wc.translate());
if (destination != null) {
destination.setRect(xmin, ymin, xmax - xmin, ymax - ymin);
} else {
@@ -455,7 +458,7 @@ public final class Shapes2D extends Static {
* method for comments about the algorithm. The code below is the same algorithm adapted for
* the 2D case and the related objects (Point2D, Rectangle2D, etc.).
*
- * The 'border' variable in the loop below is used in order to compress 2 dimensions
+ * The `border` variable in the loop below is used in order to compress 2 dimensions
* and 2 extremums in a single loop, in this order: (xmin, xmax, ymin, ymax).
*/
TransformException warning = null;
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/geometry/WraparoundInEnvelope.java b/core/sis-referencing/src/main/java/org/apache/sis/geometry/WraparoundInEnvelope.java
index cbced19..4c59995 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/geometry/WraparoundInEnvelope.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/geometry/WraparoundInEnvelope.java
@@ -37,7 +37,20 @@ import org.apache.sis.util.ArraysExt;
* from source and target envelopes, which do not need to be the same.</p>
*
* <p>The final result is that envelopes transformed using {@code WraparoundInEnvelope} may be larger
- * than envelopes transformed using {@link WraparoundTransform} but should never be smaller.</p>
+ * than envelopes transformed using {@link WraparoundTransform} but should never be smaller.
+ * For example when transforming the following envelope with wraparound on the dashed line:</p>
+ *
+ * {@preformat text
+ * ┌─┆───────────────┆───┐ ┆ Envelope to transform.
+ * │ ┆ ┆ │ ┆
+ * └─┆───────────────┆───┘ ┆
+ * ┆ ┆ ┌─────────┐ ┆ Result we would got without `WraparoundInEnvelope`.
+ * ┆ ┆ │ │ ┆
+ * ┆ ┆ └─────────┘ ┆
+ * ┆ ┌─┆───┐ ┌─┆───┐ Better result (union to be done by caller).
+ * ┆ │ ┆ │ │ ┆ │
+ * ┆ └─┆───┘ └─┆───┘
+ * }
*
* <h2>Mutability</h2>
* <b>This class is mutable.</b> This class records the translations that {@link #shift(double)} wanted to apply