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: Apply Clenshaw summation on Bennett (1996) equation 2.
Date Thu, 08 Aug 2019 06:59:59 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 dc53abf  Apply Clenshaw summation on Bennett (1996) equation 2.
dc53abf is described below

commit dc53abfa02d287783d3a323d6ccdac6b8e101855
Author: Martin Desruisseaux <martin.desruisseaux@geomatys.com>
AuthorDate: Thu Aug 8 07:09:00 2019 +0200

    Apply Clenshaw summation on Bennett (1996) equation 2.
---
 .../sis/referencing/GeodesicsOnEllipsoid.java      | 25 ++++++++++++----------
 .../apache/sis/referencing/ClenshawSummation.java  | 13 +++++++++++
 2 files changed, 27 insertions(+), 11 deletions(-)

diff --git a/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodesicsOnEllipsoid.java
b/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodesicsOnEllipsoid.java
index 9d4c05c..024ee73 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodesicsOnEllipsoid.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodesicsOnEllipsoid.java
@@ -210,7 +210,7 @@ class GeodesicsOnEllipsoid extends GeodeticCalculator {
     private double A1, A2, A3, C31, C32, C33, C34, C35;
 
     /**
-     * Coefficients for Rhumb line calculation from Bennett (1996) equation 2.
+     * Coefficients for Rhumb line calculation from Bennett (1996) equation 2, modified with
Clenshaw summation.
      */
     private final double R0, R2, R4, R6;
 
@@ -231,13 +231,13 @@ class GeodesicsOnEllipsoid extends GeodeticCalculator {
         axisRatio = b / a;
         semiMinor = b;
         /*
-         * For rhumb-line calculation.
+         * For rhumb-line calculation from Bennett (1996) equation 2, modified with Clenshaw
summation.
          */
         double fe;
-        R0 = 1  +  (eccentricitySquared)*(-1./4   + eccentricitySquared*(-3./64 + eccentricitySquared*(-5./256)));
-        R2 = (fe  = eccentricitySquared)*( 3./8   + eccentricitySquared*( 3./32 + eccentricitySquared*(45./1024)));
-        R4 = (fe *= eccentricitySquared)*(15./256 + eccentricitySquared*(45./1024));
-        R6 = (fe *  eccentricitySquared)*(35./3072);
+        R0 = 1  -  (eccentricitySquared)*(  1./4   + eccentricitySquared*( 3./64 + eccentricitySquared*(
5./256)));
+        R2 = (fe  = eccentricitySquared)*( -3./8   - eccentricitySquared*( 3./32 + eccentricitySquared*(25./768)));
+        R4 = (fe *= eccentricitySquared)*( 15./128 + eccentricitySquared*(45./512));
+        R6 = (fe *  eccentricitySquared)*(-35./768);
     }
 
     /**
@@ -977,8 +977,8 @@ class GeodesicsOnEllipsoid extends GeodeticCalculator {
             final double sinφ = sin(φm);
             S = Δλ * cos(φm) / (sin(azimuth) * sqrt(1 - eccentricitySquared*(sinφ*sinφ)));
     // Bennett equation 4.
         } else {
-            final double m1 = m(φ1);
-            final double m2 = m(φ2);
+            final double m1 = m(φ1, sinφ1);
+            final double m2 = m(φ2, sinφ2);
             S = (m2 - m1) / cos(azimuth);
             if (STORE_LOCAL_VARIABLES) {
                 store("m₁", m1);
@@ -995,9 +995,12 @@ class GeodesicsOnEllipsoid extends GeodeticCalculator {
     }
 
     /**
-     * Computes Bennett (1996) equation 2.
+     * Computes Bennett (1996) equation 2 modified with Clenshaw summation.
      */
-    private double m(final double φ) {
-        return R0*φ - R2*sin(2*φ) + R4*sin(4*φ) - R6*sin(6*φ);
+    private double m(final double φ, final double sinφ) {
+        final double cosφ = cos(φ);
+        final double sinθ = 2*sinφ*cosφ;                        // sin(2φ)
+        final double cosθ = (cosφ + sinφ) * (cosφ - sinφ);      // cos(2φ)  =  cos²φ
- sin²φ
+        return R0*φ + sinθ*(R2 + cosθ*(R4 + cosθ*R6));
     }
 }
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/referencing/ClenshawSummation.java
b/core/sis-referencing/src/test/java/org/apache/sis/referencing/ClenshawSummation.java
index b073559..fa1447a 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/referencing/ClenshawSummation.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/referencing/ClenshawSummation.java
@@ -431,6 +431,19 @@ public final class ClenshawSummation {
     }
 
     /**
+     * Coefficients for Rhumb line calculation from Bennett (1996) equation 2.
+     *
+     * @return Bennett (1996) equation 2.
+     */
+    public static ClenshawSummation Bennett2() {
+        return new ClenshawSummation(
+            new Coefficient(t(-3, 8), t(-3,  32), t(-45, 1024)),
+            new Coefficient(null,     t(15, 256), t( 45, 1024)),
+            new Coefficient(null,     null,       t(-35, 3072))
+        );
+    }
+
+    /**
      * Computes the Clenshaw summation of a series and display the code to standard output.
      * This method has been used for generating the Java code implemented in methods such
as
      * {@link GeodesicsOnEllipsoid#computeSeriesExpansionCoefficients()} and may be executed


Mime
View raw message