sis-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From desruisse...@apache.org
Subject [sis] 03/03: Fix derivative calculation and add tests.
Date Wed, 24 Apr 2019 13:57:19 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

commit d2c062570f551d9d89956068f5d6ac2012fd2dbb
Author: Martin Desruisseaux <martin.desruisseaux@geomatys.com>
AuthorDate: Wed Apr 24 15:49:12 2019 +0200

    Fix derivative calculation and add tests.
---
 .../projection/MeridionalDistanceBased.java        | 15 ++--
 .../projection/MeridionalDistanceTest.java         | 81 ++++++++++++++++++++--
 2 files changed, 84 insertions(+), 12 deletions(-)

diff --git a/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/projection/MeridionalDistanceBased.java
b/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/projection/MeridionalDistanceBased.java
index 192a097..137b031 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/projection/MeridionalDistanceBased.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/projection/MeridionalDistanceBased.java
@@ -65,7 +65,7 @@ abstract class MeridionalDistanceBased extends NormalizedProjection {
      * Coefficients for the formula implemented by the {@link #meridianArc(double, double,
double)} method.
      * Values are computed by {@link #computeCoefficients()} at construction time or after
deserialization.
      * The values depends on the form of equation implemented by {@code meridianArc(…)}.
We do not use the
-     * form published commonly found in publication since a few algebraic operations allow
to replace the
+     * form commonly found in publications since a few algebraic operations allow to replace
the
      * sin(2φ), sin(4φ), sin(6φ) and sin(8φ) terms by only sin²(φ), which is faster
to calculate.
      * See {@link #computeCoefficients()} comments in method body for more information.
      */
@@ -158,6 +158,7 @@ abstract class MeridionalDistanceBased extends NormalizedProjection {
         c2 = -2625./0x8000  * e10  +   175./0x6000  * e8  +  5120./0x60000 * e6  +  -15./0x20
* e4;
         c3 =   735./0x800   * e10  +  2240./0x60000 * e8  +   -35./0x60    * e6;
         c4 = -2205./0x1000  * e10  +  -315./0x400   * e8;
+     // c6 =   693./0x20000 * e10  omitted for now (not yet used).
     }
 
     /**
@@ -183,11 +184,11 @@ abstract class MeridionalDistanceBased extends NormalizedProjection
{
      *
      * @return the derivative at the specified latitude.
      */
-    final double dM_dφ(final double sinφ2, final double cosφ2) {
-        return c0 +
-               c1 * (sinφ2 -   cosφ2) + sinφ2*(
-               c2 * (sinφ2 - 3*cosφ2) + sinφ2*(
-               c3 * (sinφ2 - 5*cosφ2) + sinφ2*
-               c4 * (    1 - 7*cosφ2)));
+    final double dM_dφ(final double sinφ2) {
+        return ((((7 - 8*sinφ2)*c4 - 6*c3) * sinφ2
+                            + 5*c3 - 4*c2) * sinφ2
+                            + 3*c2 - 2*c1) * sinφ2
+                            +   c1
+                            +   c0;
     }
 }
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/referencing/operation/projection/MeridionalDistanceTest.java
b/core/sis-referencing/src/test/java/org/apache/sis/referencing/operation/projection/MeridionalDistanceTest.java
index 2db8320..c71cf39 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/referencing/operation/projection/MeridionalDistanceTest.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/referencing/operation/projection/MeridionalDistanceTest.java
@@ -17,8 +17,11 @@
 package org.apache.sis.referencing.operation.projection;
 
 import java.util.Random;
+import org.opengis.referencing.operation.TransformException;
+import org.apache.sis.referencing.operation.transform.AbstractMathTransform1D;
 import org.apache.sis.referencing.operation.DefaultOperationMethod;
 import org.apache.sis.test.TestUtilities;
+import org.apache.sis.test.DependsOnMethod;
 import org.apache.sis.test.DependsOn;
 import org.junit.Test;
 
@@ -37,12 +40,22 @@ import static org.junit.Assert.assertEquals;
 @DependsOn(NormalizedProjectionTest.class)
 public final strictfp class MeridionalDistanceTest extends MapProjectionTestCase {
     /**
-     * Creates a new test case.
+     * Threshold for comparison of floating point values.
      */
-    public MeridionalDistanceTest() {
+    private static final double STRICT = 0;
+
+    /**
+     * Creates the projection to be tested.
+     *
+     * @param  ellipsoidal   {@code false} for a sphere, or {@code true} for WGS84 ellipsoid.
+     * @return a test instance of the projection.
+     */
+    private MeridionalDistanceBased create(final boolean ellipsoidal) {
         final DefaultOperationMethod provider = new org.apache.sis.internal.referencing.provider.Sinusoidal();
-        transform = new Sinusoidal(provider, parameters(provider, true));
+        final Sinusoidal projection = new Sinusoidal(provider, parameters(provider, ellipsoidal));
         tolerance = NormalizedProjection.ANGULAR_TOLERANCE;     // = linear tolerance on
a sphere of radius 1.
+        transform = projection;
+        return projection;
     }
 
     /**
@@ -112,7 +125,7 @@ public final strictfp class MeridionalDistanceTest extends MapProjectionTestCase
      */
     @Test
     public void compareWithReference() {
-        final MeridionalDistanceBased projection = (MeridionalDistanceBased) transform;
+        final MeridionalDistanceBased projection = create(true);
         final Random random = TestUtilities.createRandomNumberGenerator();
         for (int i=0; i<100; i++) {
             final double φ = random.nextDouble() * PI - PI/2;
@@ -120,7 +133,65 @@ public final strictfp class MeridionalDistanceTest extends MapProjectionTestCase
             final double accurate  = referenceMoreAccurate(φ);
             final double actual    = projection.meridianArc(φ, sin(φ), cos(φ));
             assertEquals("Accurate formula disagrees with reference.", reference, accurate,
2E-10);
-            assertEquals("Implementation disagrees with reference.",   accurate,  actual,
  1E-13);
+            assertEquals("Implementation disagrees with our formula.", accurate,  actual,
  1E-13);
         }
     }
+
+    /**
+     * Compares {@link MeridionalDistanceBased#meridianArc(double, double, double)} with
spherical formula.
+     * In the spherical case, {@code meridianArc(φ)} should be equal to φ.
+     */
+    @Test
+    public void compareWithSphere() {
+        final MeridionalDistanceBased projection = create(false);
+        assertEquals("Expected spherical projection.", 0, projection.eccentricity, STRICT);
+        final Random random = TestUtilities.createRandomNumberGenerator();
+        for (int i=0; i<20; i++) {
+            final double φ = random.nextDouble() * PI - PI/2;
+            assertEquals("When excentricity=0, meridianArc(φ, sinφ, cosφ) simplify to
φ.",
+                    φ, projection.meridianArc(φ, sin(φ), cos(φ)), 1E-15);
+        }
+    }
+
+    /**
+     * Tests the {@link MeridionalDistanceBased#dM_dφ(double)} method on a sphere.
+     *
+     * @throws TransformException should never happen.
+     */
+    @Test
+    public void testDerivativeOnSphere() throws TransformException {
+        testDerivative(false);
+    }
+
+    /**
+     * Tests the {@link MeridionalDistanceBased#dM_dφ(double)} method on an ellipsoid.
+     *
+     * @throws TransformException should never happen.
+     */
+    @Test
+    @DependsOnMethod("testDerivativeOnSphere")
+    public void testDerivativeOnEllipsoid() throws TransformException {
+        testDerivative(true);
+    }
+
+    /**
+     * Implementation of tests for derivative method.
+     */
+    private void testDerivative(final boolean ellipsoidal) throws TransformException {
+        final MeridionalDistanceBased projection = create(ellipsoidal);
+        isInverseTransformSupported = false;
+        derivativeDeltas = new double[] {toRadians(0.01)};
+        transform = new AbstractMathTransform1D() {
+            @Override public double transform (final double φ) {
+                return projection.meridianArc(φ, sin(φ), cos(φ));
+            }
+            @Override public double derivative(final double φ) {
+                final double sinφ = sin(φ);
+                return projection.dM_dφ(sinφ*sinφ);
+            }
+        };
+        verifyInDomain(new double[] {toRadians(-89)},
+                       new double[] {toRadians(+89)},
+                       new int[] {ellipsoidal ? 100 : 20}, null);
+    }
 }


Mime
View raw message