Fix NZ Map projection and add a test for it. (#64)

* Update Proj4JSTest.java Add tests for NZ Map Grid (27200) * Update ProjectionMath.java Fix to zpoly1d translating pointer math into java * Update NewZealandMapGridProjection.java Scale constants to match original C and fix various errors from translating from C to java Signed-off-by: Phil Scadden <[email protected]>
locationtech · Jul 28, 2020 · 7100c00 · 7100c00
1 parent 40c349b
commit 7100c00
Show file tree

Hide file tree

Showing 3 changed files with 9 additions and 7 deletions.
diff --git a/src/main/java/org/locationtech/proj4j/proj/NewZealandMapGridProjection.java b/src/main/java/org/locationtech/proj4j/proj/NewZealandMapGridProjection.java
@@ -49,8 +49,8 @@ public class NewZealandMapGridProjection extends Projection {
 
     private final static double tpsi[] = { .6399175073, -.1358797613, .063294409, -.02526853, .0117879, -.0055161, .0026906, -.001333, .00067, -.00034 };
 
-    private final static double SECS_TO_RAD = ProjectionMath.DTR / 3600d;
-    private final static double RAD_TO_SECS = ProjectionMath.RTD * 3600d;
+    private final static double SECS_TO_RAD = 1E5 * ProjectionMath.DTR / 3600d;
+    private final static double RAD_TO_SECS = 1E-5 *ProjectionMath.RTD * 3600d;
 
 	public NewZealandMapGridProjection() {
 		initialize();
@@ -65,7 +65,7 @@ public ProjCoordinate project(double lplam, double lpphi, ProjCoordinate out) {
             p.r = tpsi[i] + lpphi * p.r;
         p.r *= lpphi;
         p.i = lplam;
-        zpoly1(p, bf);
+        p = zpoly1(p, bf);
         out.x = p.i;
         out.y = p.r;
         return out;
@@ -92,9 +92,9 @@ protected ProjCoordinate projectInverse(double x, double y, ProjCoordinate dst)
             dst.x = p.i;
             dst.y = tphi[tphi.length - 1];
             for (i = tphi.length - 1; i > 0; i--) {
-                dst.y = tphi[i] + p.r * dst.y;
+                dst.y = tphi[i-1] + p.r * dst.y;
             }
-            dst.y = projectionLongitude + p.r * dst.x * SECS_TO_RAD;
+            dst.y = projectionLatitude + p.r * dst.y * SECS_TO_RAD;
         } else
             dst.y = dst.x = Double.NaN;
         return dst;

diff --git a/src/main/java/org/locationtech/proj4j/util/ProjectionMath.java b/src/main/java/org/locationtech/proj4j/util/ProjectionMath.java
@@ -511,7 +511,7 @@ public static Complex zpoly1d(Complex z, Complex[] C, Complex der) {
                 b.r = a.r + z.r * (t = b.r) - z.i * b.i;
                 b.i = a.i + z.r * b.i + z.i * t;
             }
-            Complex c = C[i];
+            Complex c = C[i-1];
             a.r = c.r + z.r * (t = a.r) - z.i * a.i;
             a.i = c.i + z.r * a.i + z.i * t;
         }
@@ -525,7 +525,7 @@ public static Complex zpoly1d(Complex z, Complex[] C, Complex der) {
         return a;
     }
 
-    /**
+     /**
      * Tests whether the datum parameter-based transform
      * is the identity transform
      * (in which case datum transformation can be short-circuited,

diff --git a/src/test/java/org/locationtech/proj4j/Proj4JSTest.java b/src/test/java/org/locationtech/proj4j/Proj4JSTest.java
@@ -42,6 +42,8 @@ public void testGood() {
         checkTransformFromGeo("EPSG:32615", -93.0, 42.0, 500000.0, 4649776.22482, 0.1);
         checkTransformFromGeo("EPSG:3411", -32.0, 48.0, 1070076.44, -4635010.27, 2);
         checkTransformFromGeo("EPSG:3573", 9.84375, 61.875, 2923052.02009, 1054885.46559, 0.1);
+        checkTransformToGeo("EPSG:27200",2464770.343667, 6056137.861919,172.465,-40.7,0.1);
+        checkTransformFromGeo("EPSG:27200",172.465,-40.7, 2464780.81,6056330.22,0.1);
     }
 
 //    @Test