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: Add spherical case of Cassini-Soldner projection. Formulas are taken from Wikipedia. This commit completes https://issues.apache.org/jira/browse/SIS-218
Date Wed, 29 Apr 2020 15:00:27 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 887b72a  Add spherical case of Cassini-Soldner projection. Formulas are taken from
Wikipedia. This commit completes https://issues.apache.org/jira/browse/SIS-218
887b72a is described below

commit 887b72a3d0c2bcaa48e0b0d9626d66fcce4b64a4
Author: Martin Desruisseaux <martin.desruisseaux@geomatys.com>
AuthorDate: Wed Apr 29 16:59:10 2020 +0200

    Add spherical case of Cassini-Soldner projection. Formulas are taken from Wikipedia.
    This commit completes https://issues.apache.org/jira/browse/SIS-218
---
 .../operation/projection/CassiniSoldner.java       | 103 +++++++++++++++++++++
 .../operation/projection/CassiniSoldnerTest.java   |   3 +-
 2 files changed, 105 insertions(+), 1 deletion(-)

diff --git a/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/projection/CassiniSoldner.java
b/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/projection/CassiniSoldner.java
index d506692..e857a29 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/projection/CassiniSoldner.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/referencing/operation/projection/CassiniSoldner.java
@@ -17,7 +17,10 @@
 package org.apache.sis.referencing.operation.projection;
 
 import java.util.EnumMap;
+import org.opengis.util.FactoryException;
 import org.opengis.parameter.ParameterDescriptor;
+import org.opengis.referencing.operation.MathTransform;
+import org.opengis.referencing.operation.MathTransformFactory;
 import org.opengis.referencing.operation.Matrix;
 import org.opengis.referencing.operation.OperationMethod;
 import org.apache.sis.referencing.operation.matrix.Matrix2;
@@ -93,6 +96,33 @@ public class CassiniSoldner extends MeridianArcBased {
     }
 
     /**
+     * Creates a new projection initialized to the same parameters than the given one.
+     */
+    CassiniSoldner(final CassiniSoldner other) {
+        super(other);
+    }
+
+    /**
+     * Returns the sequence of <cite>normalization</cite> → {@code this} →
<cite>denormalization</cite> transforms
+     * as a whole. The transform returned by this method expects (<var>longitude</var>,
<var>latitude</var>)
+     * coordinates in <em>degrees</em> and returns (<var>x</var>,<var>y</var>)
coordinates in <em>metres</em>.
+     * The non-linear part of the returned transform will be {@code this} transform, except
if the ellipsoid
+     * is spherical. In the later case, {@code this} transform will be replaced by a simplified
implementation.
+     *
+     * @param  factory  the factory to use for creating the transform.
+     * @return the map projection from (λ,φ) to (<var>x</var>,<var>y</var>)
coordinates.
+     * @throws FactoryException if an error occurred while creating a transform.
+     */
+    @Override
+    public MathTransform createMapProjection(final MathTransformFactory factory) throws FactoryException
{
+        CassiniSoldner kernel = this;
+        if (eccentricity == 0) {
+            kernel = new Spherical(this);
+        }
+        return context.completeTransform(factory, kernel);
+    }
+
+    /**
      * Converts the specified (λ,φ) coordinate (units in radians) and stores the result
in {@code dstPts}.
      * In addition, opportunistically computes the projection derivative if {@code derivate}
is {@code true}.
      *
@@ -181,4 +211,77 @@ public class CassiniSoldner extends MeridianArcBased {
         dstPts[dstOff  ]   = (D - T1*(D2*D)/3 + (1 + 3*T1)*T1*(D4*D)/15) / cosφ1;
         dstPts[dstOff+1]   = φ1 - (tanφ1 / ρ1_ν1) * (D2/2 - (1 + 3*T1)*D4 / 24);
     }
+
+
+    /**
+     * Provides the transform equations for the spherical case of the Cassini-Soldner projection.
+     * Formulas are available on <a href="https://en.wikipedia.org/wiki/Cassini_projection">Wikipedia</a>.
+     *
+     * @author  Martin Desruisseaux (Geomatys)
+     * @author  Rémi Maréchal (Geomatys)
+     * @version 1.1
+     * @since   1.1
+     * @module
+     */
+    static final class Spherical extends CassiniSoldner {
+        /**
+         * For cross-version compatibility.
+         */
+        private static final long serialVersionUID = -8131887916538379858L;
+
+        /**
+         * Constructs a new map projection from the parameters of the given projection.
+         *
+         * @param  other  the other projection (usually ellipsoidal) from which to copy the
parameters.
+         */
+        protected Spherical(final CassiniSoldner other) {
+            super(other);
+        }
+
+        /**
+         * {@inheritDoc}
+         */
+        @Override
+        public Matrix transform(final double[] srcPts, final int srcOff,
+                                final double[] dstPts, final int dstOff,
+                                final boolean derivate)
+        {
+            final double λ    = srcPts[srcOff  ];
+            final double φ    = srcPts[srcOff+1];
+            final double sinλ = sin(λ);
+            final double cosλ = cos(λ);
+            final double sinφ = sin(φ);
+            final double cosφ = cos(φ);
+            final double tanφ = sinφ / cosφ;
+            if (dstPts != null) {
+                dstPts[dstOff  ] = asin (cosφ * sinλ);
+                dstPts[dstOff+1] = atan2(tanφ,  cosλ);
+            }
+            if (!derivate) {
+                return null;
+            }
+            final double cosφ2 = cosφ * cosφ;
+            final double sinλ2 = sinλ * sinλ;
+            final double dxden = sqrt(1 - sinλ2*cosφ2);
+            final double dyden = tanφ*tanφ + cosλ*cosλ;
+            return new Matrix2(
+                    cosλ*cosφ  / dxden,
+                   -sinλ*sinφ  / dxden,
+                    sinλ*tanφ  / dyden,
+                    cosλ/cosφ2 / dyden);
+        }
+
+        /**
+         * {@inheritDoc}
+         */
+        @Override
+        protected void inverseTransform(final double[] srcPts, final int srcOff,
+                                        final double[] dstPts, final int dstOff)
+        {
+            final double x = srcPts[srcOff  ];
+            final double y = srcPts[srcOff+1];
+            dstPts[dstOff  ] = atan2(tan(x),  cos(y));
+            dstPts[dstOff+1] = asin (sin(y) * cos(x));
+        }
+    }
 }
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/referencing/operation/projection/CassiniSoldnerTest.java
b/core/sis-referencing/src/test/java/org/apache/sis/referencing/operation/projection/CassiniSoldnerTest.java
index bf4dd67..c3f1b28 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/referencing/operation/projection/CassiniSoldnerTest.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/referencing/operation/projection/CassiniSoldnerTest.java
@@ -59,11 +59,12 @@ public final strictfp class CassiniSoldnerTest extends MapProjectionTestCase
{
         if (ellipse) {
             pg.parameter("semi-major").setValue(WGS84_A);
             pg.parameter("semi-minor").setValue(WGS84_B);
+            return new CassiniSoldner(provider, pg);
         } else {
             pg.parameter("semi-major").setValue(RADIUS);
             pg.parameter("semi-minor").setValue(RADIUS);
+            return new CassiniSoldner.Spherical(new CassiniSoldner(provider, pg));
         }
-        return new CassiniSoldner(provider, pg);
     }
 
     /**


Mime
View raw message