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: Fix the case of straight line at equator and add more tests.
Date Wed, 15 May 2019 16:39:35 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 321d86b  Fix the case of straight line at equator and add more tests.
321d86b is described below

commit 321d86b58946b3a172e8cd65d75eae2e106e661d
Author: Martin Desruisseaux <martin.desruisseaux@geomatys.com>
AuthorDate: Wed May 15 15:27:02 2019 +0200

    Fix the case of straight line at equator and add more tests.
---
 .../referencing/GeodeticObjectBuilder.java         |  57 ++++++
 .../internal/referencing/j2d/ShapeUtilities.java   |  17 +-
 .../apache/sis/referencing/GeodeticCalculator.java |  19 ++
 .../sis/referencing/GeodeticCalculatorTest.java    | 214 ++++++++++++++++-----
 4 files changed, 251 insertions(+), 56 deletions(-)

diff --git a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/GeodeticObjectBuilder.java
b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/GeodeticObjectBuilder.java
index de8fbeb..3bdf76a 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/GeodeticObjectBuilder.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/GeodeticObjectBuilder.java
@@ -21,6 +21,7 @@ import java.util.Date;
 import java.util.Collections;
 import javax.measure.Unit;
 import javax.measure.quantity.Time;
+import javax.measure.quantity.Length;
 import org.opengis.util.FactoryException;
 import org.opengis.parameter.ParameterValueGroup;
 import org.opengis.parameter.ParameterNotFoundException;
@@ -34,6 +35,9 @@ import org.opengis.referencing.cs.AxisDirection;
 import org.opengis.referencing.cs.CSFactory;
 import org.opengis.referencing.cs.CartesianCS;
 import org.opengis.referencing.cs.TimeCS;
+import org.opengis.referencing.datum.Ellipsoid;
+import org.opengis.referencing.datum.DatumFactory;
+import org.opengis.referencing.datum.GeodeticDatum;
 import org.opengis.referencing.datum.TemporalDatum;
 import org.opengis.referencing.operation.OperationMethod;
 import org.opengis.referencing.operation.Conversion;
@@ -62,6 +66,11 @@ import org.apache.sis.referencing.CommonCRS;
  */
 public class GeodeticObjectBuilder extends Builder<GeodeticObjectBuilder> {
     /**
+     * The geodetic datum, or {@code null} if none.
+     */
+    private GeodeticDatum datum;
+
+    /**
      * The name of the conversion to use for creating a {@code ProjectedCRS} or {@code DerivedCRS}.
      * This name is for information purpose; its value does not impact the numerical results
of coordinate operations.
      *
@@ -126,6 +135,26 @@ public class GeodeticObjectBuilder extends Builder<GeodeticObjectBuilder>
{
     }
 
     /**
+     * Creates a geodetic datum with an ellipsoid of the given shape.
+     *
+     * @param  name               ellipsoid and datum name.
+     * @param  semiMajorAxis      equatorial radius in supplied linear units.
+     * @param  inverseFlattening  eccentricity of ellipsoid. An infinite value creates a
sphere.
+     * @param  units              linear units of major axis.
+     * @return {@code this}, for method call chaining.
+     * @throws FactoryException if the datum can not be created.
+     */
+    public GeodeticObjectBuilder setFlattenedSphere(final String name, final double semiMajorAxis,
+            final double inverseFlattening, final Unit<Length> units) throws FactoryException
+    {
+        final DatumFactory factory = factories.getDatumFactory();
+        final Ellipsoid ellipsoid = factory.createFlattenedSphere(
+                Collections.singletonMap(Ellipsoid.NAME_KEY, name), semiMajorAxis, inverseFlattening,
units);
+        datum = factory.createGeodeticDatum(name(ellipsoid), ellipsoid, CommonCRS.WGS84.primeMeridian());
+        return this;
+    }
+
+    /**
      * Sets the conversion method to use for creating a {@code ProjectedCRS} or {@code DerivedCRS}.
      * The method is typically a map projection method. Examples:
      *
@@ -328,6 +357,34 @@ public class GeodeticObjectBuilder extends Builder<GeodeticObjectBuilder>
{
     }
 
     /**
+     * Creates a projected CRS with base CRS on the specified datum and with default axes.
+     * The base CRS uses the ellipsoid specified by {@link #setFlattenedSphere(String, double,
double, Unit)}.
+     *
+     * @return the projected CRS.
+     * @throws FactoryException if an error occurred while building the projected CRS.
+     */
+    public ProjectedCRS createProjectedCRS() throws FactoryException {
+        GeographicCRS crs = CommonCRS.WGS84.geographic();
+        if (datum != null) {
+            crs = factories.getCRSFactory().createGeographicCRS(name(datum), datum, crs.getCoordinateSystem());
+        }
+        return createProjectedCRS(crs, ReferencingUtilities.standardProjectedCS(factories.getCSAuthorityFactory()));
+    }
+
+    /**
+     * Creates a geographic CRS.
+     *
+     * @param  normalized  whether axes should be in (longitude, latitude) order instead
than (latitude, longitude).
+     * @return the geographic coordinate reference system.
+     * @throws FactoryException if an error occurred while building the geographic CRS.
+     */
+    public GeographicCRS createGeographicCRS(final boolean normalized) throws FactoryException
{
+        final GeographicCRS crs = normalized ? CommonCRS.WGS84.geographic() : CommonCRS.defaultGeographic();
+        if (datum != null) properties.putIfAbsent(GeographicCRS.NAME_KEY, datum.getName());
+        return factories.getCRSFactory().createGeographicCRS(properties, datum, crs.getCoordinateSystem());
+    }
+
+    /**
      * Creates a temporal CRS from the given origin and temporal unit. For this method, the
CRS name is optional:
      * if no {@code addName(…)} method has been invoked, then a default name will be used.
      *
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/ShapeUtilities.java
b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/ShapeUtilities.java
index 3a240c7..2226f4e 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/ShapeUtilities.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/ShapeUtilities.java
@@ -405,18 +405,23 @@ public final class ShapeUtilities extends Static {
         ym -= y1;
         final double Δx = x2 - x1;
         final double Δy = y2 - y1;
-        if ((isInfinite(α1) || abs(Δx*α1 - Δy) <= εy) &&
-            (isInfinite(α2) || abs(Δx*α2 - Δy) <= εy))
+        final boolean isVertical = abs(Δx) <= εx;
+        if (isVertical || ((isInfinite(α1) || abs(Δx*α1 - Δy) <= εy)
+                       &&  (isInfinite(α2) || abs(Δx*α2 - Δy) <= εy)))
         {
             /*
              * Following tests are partially redundant with above tests, but may detect a
larger
-             * error than above tests did. They are also necessary is a derivative was infinite.
+             * error than above tests did. They are also necessary if a derivative was infinite,
+             * or if `isVertical` was true.
              */
-            if ((α1 == 0 || abs(Δy/α1 - Δx) <= εx) &&
-                (α2 == 0 || abs(Δy/α2 - Δx) <= εx))
+            final boolean isHorizontal = abs(Δy) <= εy;
+            if (isHorizontal || ((α1 == 0 || abs(Δy/α1 - Δx) <= εx)
+                             &&  (α2 == 0 || abs(Δy/α2 - Δx) <= εx)))
             {
                 final double slope = Δy / Δx;
-                if (abs(xm*slope - ym) <= εy && abs(ym/slope - xm) <= εx)
{
+                if ((isVertical   || abs(xm*slope - ym) <= εy) &&
+                    (isHorizontal || abs(ym/slope - xm) <= εx))
+                {
                     return new Line2D.Double(x1, y1, x2, y2);
                 }
             }
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodeticCalculator.java
b/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodeticCalculator.java
index 27b4687..cbaff00 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodeticCalculator.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodeticCalculator.java
@@ -223,6 +223,7 @@ public class GeodeticCalculator {
      * @throws IllegalArgumentException if the latitude is out of bounds.
      *
      * @see #setEndPoint(double, double)
+     * @see #moveToEndPoint()
      */
     public void setStartPoint(final double latitude, final double longitude) {
         ArgumentChecks.ensureBetween("latitude", Latitude.MIN_VALUE, Latitude.MAX_VALUE,
latitude);
@@ -525,6 +526,24 @@ public class GeodeticCalculator {
     }
 
     /**
+     * Sets the start point and starting azimuth to the current end point and ending azimuth
values.
+     * The {@linkplain #getEndingAzimuth() ending azimuths}, the {@linkplain #getGeodesicDistance()
+     * geodesic distance} and the {@linkplain #getEndPoint() end point} are discarded by
this method call;
+     * some of them will need to be specified again.
+     *
+     * @see #setStartPoint(double, double)
+     */
+    public void moveToEndPoint() {
+        if ((validity & END_POINT) == 0) {
+            computeEndPoint();
+        }
+        φ1 = φ2;
+        λ1 = λ2;
+        α1 = α2;
+        validity = START_POINT | STARTING_AZIMUTH;
+    }
+
+    /**
      * Creates an approximation of the geodesic track from start point to end point as a
Java2D object.
      * The coordinates are expressed in the coordinate reference system specified at creation
time.
      * The approximation uses linear, quadratic or cubic Bézier curves.
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/referencing/GeodeticCalculatorTest.java
b/core/sis-referencing/src/test/java/org/apache/sis/referencing/GeodeticCalculatorTest.java
index 1a5c9a2..1d2daff 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/referencing/GeodeticCalculatorTest.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/referencing/GeodeticCalculatorTest.java
@@ -18,6 +18,7 @@ package org.apache.sis.referencing;
 
 import java.awt.Shape;
 import java.awt.geom.Point2D;
+import java.awt.geom.PathIterator;
 import java.util.Random;
 import org.opengis.geometry.DirectPosition;
 import org.opengis.referencing.operation.TransformException;
@@ -31,7 +32,8 @@ import org.apache.sis.test.TestCase;
 import org.junit.Test;
 
 import static java.lang.StrictMath.*;
-import static org.junit.Assert.*;
+import static org.opengis.test.Assert.*;
+import static org.apache.sis.internal.metadata.ReferencingServices.NAUTICAL_MILE;
 
 
 /**
@@ -43,11 +45,40 @@ import static org.junit.Assert.*;
  */
 public final strictfp class GeodeticCalculatorTest extends TestCase {
     /**
-     * Returns the calculator to use for testing purpose.
-     * This classes uses calculator for a sphere.
+     * Verifies that the given point is equals to the given latitude and longitude.
+     *
+     * @param φ  the expected latitude value, in degrees.
+     * @param λ  the expected longitude value, in degrees.
+     * @param p  the actual position to verify.
+     * @param ε  the tolerance threshold.
      */
-    private static GeodeticCalculator create() {
-        return new GeodeticCalculator(CommonCRS.SPHERE.geographic());
+    private static void assertPositionEquals(final double φ, final double λ, final DirectPosition
p, final double ε) {
+        assertEquals("φ", φ, p.getOrdinate(0), ε);
+        assertEquals("λ", λ, p.getOrdinate(1), ε);
+    }
+
+    /**
+     * Asserts that a Java2D point is equal to the expected value. Used for verifying geodesic
paths.
+     *
+     * @param x  the expected <var>x</var> coordinates.
+     * @param y  the expected <var>y</var> coordinates.
+     * @param p  the actual position to verify.
+     * @param ε  the tolerance threshold.
+     */
+    private static void assertPointEquals(final double x, final double y, final Point2D p,
final double ε) {
+        assertEquals("x", x, p.getX(), ε);
+        assertEquals("y", y, p.getY(), ε);
+    }
+
+    /**
+     * Returns the calculator to use for testing purpose. This classes uses a calculator
for a sphere.
+     *
+     * @param  normalized  whether to force (longitude, latitude) axis order.
+     */
+    private static GeodeticCalculator create(final boolean normalized) {
+        return new GeodeticCalculator(normalized
+                ? CommonCRS.SPHERE.normalizedGeographic()
+                : CommonCRS.SPHERE.geographic());
     }
 
     /**
@@ -57,46 +88,46 @@ public final strictfp class GeodeticCalculatorTest extends TestCase {
      */
     @Test
     public void testCardinalAzimuths() {
-        final GeodeticCalculator c = create();
-        final double EPS = 0.2;
+        final GeodeticCalculator c = create(false);
+        final double tolerance = 0.2;
         c.setStartPoint(20, 12);
-        c.setEndPoint(20, 13);  assertEquals("East",   90, c.getStartingAzimuth(), EPS);
-        c.setEndPoint(21, 12);  assertEquals("North",   0, c.getStartingAzimuth(), EPS);
-        c.setEndPoint(20, 11);  assertEquals("West",  -90, c.getStartingAzimuth(), EPS);
-        c.setEndPoint(19, 12);  assertEquals("South", 180, c.getStartingAzimuth(), EPS);
+        c.setEndPoint(20, 13);  assertEquals("East",   90, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint(21, 12);  assertEquals("North",   0, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint(20, 11);  assertEquals("West",  -90, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint(19, 12);  assertEquals("South", 180, c.getStartingAzimuth(), tolerance);
     }
 
     /**
      * Tests azimuths at poles.
      */
     @Test
-    public void testPoles() {
-        final GeodeticCalculator c = create();
-        final double EPS = 0.2;
+    public void testAzimuthAtPoles() {
+        final GeodeticCalculator c = create(false);
+        final double tolerance = 0.2;
         c.setStartPoint( 90,  30);
-        c.setEndPoint  ( 20,  20);  assertEquals(-170, c.getStartingAzimuth(), EPS);
-        c.setEndPoint  ( 20,  40);  assertEquals( 170, c.getStartingAzimuth(), EPS);
-        c.setEndPoint  ( 20,  30);  assertEquals( 180, c.getStartingAzimuth(), EPS);
-        c.setEndPoint  (-20,  30);  assertEquals( 180, c.getStartingAzimuth(), EPS);
-        c.setEndPoint  (-90,  30);  assertEquals( 180, c.getStartingAzimuth(), EPS);
+        c.setEndPoint  ( 20,  20);  assertEquals(-170, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint  ( 20,  40);  assertEquals( 170, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint  ( 20,  30);  assertEquals( 180, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint  (-20,  30);  assertEquals( 180, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint  (-90,  30);  assertEquals( 180, c.getStartingAzimuth(), tolerance);
 
         c.setStartPoint( 90,   0);
-        c.setEndPoint  ( 20,  20);  assertEquals( 160, c.getStartingAzimuth(), EPS);
-        c.setEndPoint  ( 20, -20);  assertEquals(-160, c.getStartingAzimuth(), EPS);
-        c.setEndPoint  ( 20,   0);  assertEquals( 180, c.getStartingAzimuth(), EPS);
-        c.setEndPoint  (-90,   0);  assertEquals( 180, c.getStartingAzimuth(), EPS);
+        c.setEndPoint  ( 20,  20);  assertEquals( 160, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint  ( 20, -20);  assertEquals(-160, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint  ( 20,   0);  assertEquals( 180, c.getStartingAzimuth(), tolerance);
+        c.setEndPoint  (-90,   0);  assertEquals( 180, c.getStartingAzimuth(), tolerance);
     }
 
     /**
      * Tests geodesic distances on the equator.
      */
     @Test
-    public void testEquator() {
+    public void testDistanceAtEquator() {
         final Random random = TestUtilities.createRandomNumberGenerator();
-        final GeodeticCalculator c = create();
+        final GeodeticCalculator c = create(false);
         final double r = c.ellipsoid.getSemiMajorAxis() * (PI / 180);
         c.setStartPoint(0, 0);
-        for (double i=0; i<100; i++) {
+        for (int i=0; i<100; i++) {
             final double x = 180 * random.nextDouble();
             c.setEndPoint(0, x);
             assertEquals(x * r, c.getGeodesicDistance(), Formulas.LINEAR_TOLERANCE);
@@ -108,10 +139,12 @@ public final strictfp class GeodeticCalculatorTest extends TestCase
{
      * This computes the great circle route from Valparaíso (33°N 71.6W) to Shanghai (31.4°N
121.8°E).
      *
      * @throws TransformException if an error occurred while transforming coordinates.
+     *
+     * @see <a href="https://en.wikipedia.org/wiki/Great-circle_navigation#Example">Great-circle
navigation on Wikipedia</a>
      */
     @Test
     public void testWikipediaExample() throws TransformException {
-        final GeodeticCalculator c = create();
+        final GeodeticCalculator c = create(false);
         c.setStartPoint(-33.0, -71.6);          // Valparaíso
         c.setEndPoint  ( 31.4, 121.8);          // Shanghai
         /*
@@ -140,19 +173,6 @@ public final strictfp class GeodeticCalculatorTest extends TestCase {
     }
 
     /**
-     * Verifies that the given point is equals to the given latitude and longitude.
-     *
-     * @param φ  the expected latitude value, in degrees.
-     * @param λ  the expected longitude value, in degrees.
-     * @param p  the position to verify.
-     * @param ε  the tolerance.
-     */
-    private static void assertPositionEquals(final double φ, final double λ, final DirectPosition
p, final double ε) {
-        assertEquals("φ", φ, p.getOrdinate(0), ε);
-        assertEquals("λ", λ, p.getOrdinate(1), ε);
-    }
-
-    /**
      * Tests geodetic calculator involving a coordinate operation.
      * This test uses a simple CRS with only the axis order interchanged.
      * The coordinates are the same than {@link #testWikipediaExample()}.
@@ -162,7 +182,7 @@ public final strictfp class GeodeticCalculatorTest extends TestCase {
     @Test
     @DependsOnMethod("testWikipediaExample")
     public void testUsingTransform() throws TransformException {
-        final GeodeticCalculator c = new GeodeticCalculator(CommonCRS.SPHERE.normalizedGeographic());
+        final GeodeticCalculator c = create(true);
         final double φ = -33.0;
         final double λ = -71.6;
         c.setStartPoint(new DirectPosition2D(λ, φ));
@@ -176,26 +196,120 @@ public final strictfp class GeodeticCalculatorTest extends TestCase
{
     /**
      * Tests {@link GeodeticCalculator#toGeodesicPath2D(double)}. This method uses a CRS
      * that swap axis order as a way to verify that user-specified CRS is taken in account.
+     * The tested coordinates are from Wikipedia example.
      *
      * @throws TransformException if an error occurred while transforming coordinates.
      */
     @Test
     @DependsOnMethod("testUsingTransform")
-    public void testToGeodesicPath2D() throws TransformException {
-        final GeodeticCalculator c = new GeodeticCalculator(CommonCRS.SPHERE.normalizedGeographic());
+    public void testGeodesicPath2D() throws TransformException {
+        final GeodeticCalculator c = create(true);
+        final double tolerance = 0.05;
         c.setStartPoint(-33.0, -71.6);                  // Valparaíso
         c.setEndPoint  ( 31.4, 121.8);                  // Shanghai
         final Shape path = c.toGeodesicPath2D(1000);
-        assertPointEquals( -71.6, -33.0, ShapeUtilities.pointOnBezier(path, 0));
-        assertPointEquals(-238.2,  31.4, ShapeUtilities.pointOnBezier(path, 1));        
   // λ₂ = 121.8° - 360°
-        assertPointEquals(-159.2,  -6.8, ShapeUtilities.pointOnBezier(path, 0.5));
+        assertPointEquals( -71.6, -33.0, ShapeUtilities.pointOnBezier(path, 0),   tolerance);
+        assertPointEquals(-238.2,  31.4, ShapeUtilities.pointOnBezier(path, 1),   tolerance);
      // λ₂ = 121.8° - 360°
+        assertPointEquals(-159.2,  -6.8, ShapeUtilities.pointOnBezier(path, 0.5), tolerance);
+    }
+
+    /**
+     * Verifies that all <var>y</var> coordinates are zero for a geodesic path
on equator.
+     *
+     * @throws TransformException if an error occurred while transforming coordinates.
+     */
+    @Test
+    public void testGeodesicPathOnEquator() throws TransformException {
+        final GeodeticCalculator c = create(false);
+        final double tolerance = 1E-12;
+        c.setStartPoint(0, 20);
+        c.setEndPoint  (0, 12);
+        assertEquals(-90, c.getStartingAzimuth(), tolerance);
+        assertEquals(-90, c.getEndingAzimuth(),   tolerance);
+        final Shape geodeticCurve = c.toGeodesicPath2D(1);
+        final double[] coords = new double[2];
+        for (final PathIterator it = geodeticCurve.getPathIterator(null, 1); !it.isDone();
it.next()) {
+            it.currentSegment(coords);
+            assertEquals ("φ",  0, coords[0], tolerance);
+            assertBetween("λ", 12, 20, coords[1]);
+        }
     }
 
     /**
-     * Asserts that a Java2D point is equal to the expected value.
+     * Tests path on the parallel at 45°N. This tests Data for this test have been taken
from
+     * <a href="http://perso.univ-lemans.fr/~hainry/articles/loxonavi.html">Orthodromie
et loxodromie</a> page.
+     *
+     * @throws TransformException if an error occurred while transforming coordinates.
+     */
+    @Test
+    public void testOnParallel45() throws TransformException {
+        /*
+         * Following numbers assume a radius R = 60 * 180/π nautical miles ≈ 6366 km.
+         * (compare to 6371 km for authalic sphere, or (6356 / 6378) km for WGS84 semi-minor/major
axis length).
+         *
+         * Column 1: Longitude difference in degrees.
+         * Column 2: Geodesic distance in kilometers.
+         * Column 3: Rhumb line distance in kilometers.
+         */
+        final double[] data = {
+              0.00,      0,      0,
+             11.25,    883,    884,
+             22.50,   1762,   1768,
+             33.75,   2632,   2652,
+             45.00,   3489,   3536,
+             56.25,   4327,   4419,
+             67.50,   5140,   5303,
+             78.75,   5923,   6187,
+             90.00,   6667,   7071,
+            101.25,   7363,   7955,
+            112.50,   8002,   8839,
+            123.75,   8573,   9723,
+            135.00,   9064,  10607,
+            146.25,   9463,  11490,
+            157.50,   9758,  12374,
+            168.75,   9939,  13258,
+            180.00,  10000,  14142
+        };
+        final GeodeticCalculator c = create(false);
+        final double toTestValue = (60 * NAUTICAL_MILE * 180/PI) / 6371007 / 1000;
+        for (int i=0; i<data.length; i+=3) {
+            c.setStartPoint(45, 0);
+            c.setEndPoint(45, data[i]);
+            double geodesic  = c.getGeodesicDistance();
+//          double rhumbLine = c.getRhumbLineDistance();
+            geodesic  *= toTestValue;
+//          rhumbLine *= toTestValue;
+            final double tolerance = data[i+1] / 1000;           // In kilometres.
+            assertEquals("Geodesic distance",   data[i+1], geodesic,  tolerance);
+//          assertEquals("Rhumb line distance", data[i+2], rhumbLine, tolerance);
+            assertEquals("Distance measured along geodesic path", geodesic, length(c) * toTestValue,
tolerance * 20);
+        }
+    }
+
+    /**
+     * Measures an estimation of the length of the path returned by {@link GeodeticCalculator#toGeodesicPath2D(double)}.
+     * This method iterates over line segments and use the given calculator for computing
the geodesic distance of each
+     * segment. The state of the given calculator is modified by this method.
      */
-    private static void assertPointEquals(final double x, final double y, final Point2D actual)
{
-        assertEquals("x", x, actual.getX(), 0.05);
-        assertEquals("y", y, actual.getY(), 0.05);
+    private static double length(final GeodeticCalculator c) throws TransformException {
+        final PathIterator iterator = c.toGeodesicPath2D(10).getPathIterator(null, 100);
+        final double[] buffer = new double[2];
+        double length=0;
+        while (!iterator.isDone()) {
+            switch (iterator.currentSegment(buffer)) {
+                default: fail("Unexpected path"); break;
+                case PathIterator.SEG_MOVETO: {
+                    c.setStartPoint(buffer[0], buffer[1]);
+                    break;
+                }
+                case PathIterator.SEG_LINETO: {
+                    c.setEndPoint(buffer[0], buffer[1]);
+                    length += c.getGeodesicDistance();
+                    c.moveToEndPoint();
+                }
+            }
+            iterator.next();
+        }
+        return length;
     }
 }


Mime
View raw message