Oblique Mercator Projection Improvements (#72)

* no_uoff, gamma, and oblique mercator projection updates and support * test fix for newly passing epsg codes * minor documentation change * Upd CHANGELOG.md Co-authored-by: Grigory Pomadchin <[email protected]>
locationtech · Jun 17, 2021 · 21c21b1 · 21c21b1
1 parent f5a026e
commit 21c21b1
Show file tree

Hide file tree

Showing 8 changed files with 80 additions and 39 deletions.
diff --git a/CHANGELOG.md b/CHANGELOG.md
@@ -6,6 +6,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [Unreleased]
 
+### Fixed
+- Problem with omerc projection e.g. EPSG: 3375 [#21](https://github.com/locationtech/proj4j/issues/21)
+
 ## [1.1.2] - 2021-04-12
 
 ### Fixed

diff --git a/src/main/java/org/locationtech/proj4j/parser/Proj4Keyword.java b/src/main/java/org/locationtech/proj4j/parser/Proj4Keyword.java
@@ -67,7 +67,7 @@ public class Proj4Keyword {
     public static final String nadgrids = "nadgrids";
     public static final String no_defs = "no_defs";
     public static final String wktext = "wktext";
-    public static final String no_uoff = "no_uoff"; // TODO: Implement no_uoff parameter
+    public static final String no_uoff = "no_uoff";
 
 
     private static Set<String> supportedParams = null;
@@ -111,12 +111,12 @@ public static synchronized Set supportedParameters() {
             supportedParams.add(axis);
 
             supportedParams.add(gamma);       // Just for Oblique Mercator projection
+            supportedParams.add(no_uoff);     // Just for Oblique Mercator projection
             supportedParams.add(zone);        // Just for Transverse Mercator projection
 
             supportedParams.add(title);       // no-op
             supportedParams.add(no_defs);     // no-op
             supportedParams.add(wktext);      // no-op
-            supportedParams.add(no_uoff);     // no-op
         }
         return supportedParams;
     }

diff --git a/src/main/java/org/locationtech/proj4j/parser/Proj4Parser.java b/src/main/java/org/locationtech/proj4j/parser/Proj4Parser.java
@@ -122,9 +122,11 @@ private Projection parseProjection(Map params, Ellipsoid ellipsoid) {
             s = (String) params.get(Proj4Keyword.k);
         if (s != null)
             projection.setScaleFactor(Double.parseDouble(s));
+        if (params.containsKey(Proj4Keyword.no_uoff))
+            projection.setNoUoff(true);
         s = (String) params.get(Proj4Keyword.gamma);
         if (s != null)
-            projection.setGamma(Double.parseDouble(s) * ProjectionMath.DTR);
+            projection.setGammaDegrees(Double.parseDouble(s));
 
         s = (String) params.get(Proj4Keyword.units);
         if (s != null) {

diff --git a/src/main/java/org/locationtech/proj4j/proj/ObliqueMercatorProjection.java b/src/main/java/org/locationtech/proj4j/proj/ObliqueMercatorProjection.java
@@ -34,7 +34,7 @@ public class ObliqueMercatorProjection extends CylindricalProjection {
 	private final static double TOL	= 1.0e-7;
 
 	private double lamc, lam1, phi1, lam2, phi2, Gamma, al, bl, el, singam, cosgam, sinrot, cosrot, u_0;
-	private boolean ellips, rot;
+	private boolean ellips, rot, no_uoff;
 
 	public ObliqueMercatorProjection() {
 		ellipsoid = Ellipsoid.WGS84;
@@ -64,14 +64,16 @@ public ObliqueMercatorProjection(Ellipsoid ellipsoid, double lon_0, double lat_0
 
 	public void initialize() {
 		super.initialize();
-		double con, com, cosphi0, d, f, h, l, sinphi0, p, j;
+		double con, com, cosphi0, d, f, h, l, sinphi0, p, j, gamma0;
 
 		//FIXME-setup rot, alpha, longc,lon/lat1/2
 		rot = true;
     lamc = lonc;
 
     // true if alpha provided
     int azi = Double.isNaN(alpha) ? 0 : 1;
+    // true if gamma provided
+    int gzi = Double.isNaN(Gamma) ? 0 : 1;
 		if (azi != 0) { // alpha specified
 			if (Math.abs(alpha) <= TOL ||
 				Math.abs(Math.abs(projectionLatitude) - ProjectionMath.HALFPI) <= TOL ||
@@ -116,10 +118,18 @@ public void initialize() {
 			al = scaleFactor;
 			el = d = f = 1.;
 		}
-		if (azi != 0) {
-			Gamma = Math.asin(Math.sin(alpha) / d);
+		if (azi != 0 || gzi != 0) {
+			if (azi != 0) {
+				gamma0 = Math.asin(Math.sin(alpha) / d);
+				if(gzi == 0) {
+					Gamma = alpha;
+				}
+			}else {
+				gamma0 = Gamma;
+				alpha = Math.asin(d * Math.sin(gamma0));
+			}
 			projectionLongitude = lamc - Math.asin((.5 * (f - 1. / f)) *
-			   Math.tan(Gamma)) / bl;
+					   Math.tan(gamma0)) / bl;
 		} else {
 			if (!spherical) {
 				h = Math.pow(ProjectionMath.tsfn(phi1, Math.sin(phi1), e), bl);
@@ -138,18 +148,17 @@ else if (con > Math.PI)
 				lam2 += ProjectionMath.TWOPI;
 			projectionLongitude = ProjectionMath.normalizeLongitude(.5 * (lam1 + lam2) - Math.atan(
 			   j * Math.tan(.5 * bl * (lam1 - lam2)) / p) / bl);
-			Gamma = Math.atan(2. * Math.sin(bl * ProjectionMath.normalizeLongitude(lam1 - projectionLongitude)) /
+			gamma0 = Math.atan(2. * Math.sin(bl * ProjectionMath.normalizeLongitude(lam1 - projectionLongitude)) /
 			   (f - 1. / f));
-			alpha = Math.asin(d * Math.sin(Gamma));
+			Gamma = Math.asin(d * Math.sin(gamma0));
+			alpha = Gamma;
 		}
-		singam = Math.sin(Gamma);
-		cosgam = Math.cos(Gamma);
-//		f = MapMath.param(params, "brot_conv").i ? Gamma : alpha;
-		f = alpha;//FIXME
-		sinrot = Math.sin(f);
-		cosrot = Math.cos(f);
-//		u_0 = MapMath.param(params, "bno_uoff").i ? 0. :
-		u_0 = false ? 0. ://FIXME
+		singam = Math.sin(gamma0);
+		cosgam = Math.cos(gamma0);
+		sinrot = Math.sin(Gamma);
+		cosrot = Math.cos(Gamma);
+
+		u_0 = no_uoff ? 0. :
 			Math.abs(al * Math.atan(Math.sqrt(d * d - 1.) / cosrot) / bl);
 		if (projectionLatitude < 0.)
 			u_0 = - u_0;
@@ -158,6 +167,10 @@ else if (con > Math.PI)
     @Override public void setGamma(double gamma) {
         this.Gamma = gamma;
     }
+
+    @Override public void setNoUoff(boolean no_uoff) {
+    	this.no_uoff = no_uoff;
+    }
 
 	public ProjCoordinate project(double lam, double phi, ProjCoordinate xy) {
 		double con, q, s, ul, us, vl, vs;

diff --git a/src/main/java/org/locationtech/proj4j/proj/Projection.java b/src/main/java/org/locationtech/proj4j/proj/Projection.java
@@ -639,6 +639,13 @@ public double getProjectionLatitude2Degrees() {
         return projectionLatitude2*RTD;
     }
 
+    /**
+     * Sets the alpha value.
+     */
+    public void setAlpha( double alpha ) {
+        this.alpha = alpha;
+    }
+
     /**
      * Sets the alpha value.
      */
@@ -655,12 +662,19 @@ public double getAlpha()
     {
         return alpha;
     }
-
+
+    /**
+     * Sets the lonc value.
+     */
+    public void setLonC( double lonc ) {
+        this.lonc = lonc;
+    }
+
     /**
      * Sets the lonc value.
      */
     public void setLonCDegrees( double lonc ) {
-        this.lonc = DTR * lonc;
+    	this.lonc = DTR * lonc;
     }
 
     /**
@@ -819,6 +833,14 @@ public static double normalizeLongitudeRadians( double angle ) {
     public void setGamma(double gamma) {
         // no-op, overridden for Oblique Mercator
     }
+
+    public void setGammaDegrees(double gamma) {
+    	setGamma(DTR * gamma);
+    }
+
+    public void setNoUoff(boolean no_uoff) {
+        // no-op, overridden for Oblique Mercator
+    }
 
     /** Is this "projection" longlat? Overridden in LongLatProjection. */
     public Boolean isGeographic() {

diff --git a/src/test/java/org/locationtech/proj4j/FeatureTest.java b/src/test/java/org/locationtech/proj4j/FeatureTest.java
@@ -44,10 +44,10 @@ public void NOTSUPPORTED_testPrimeMeridian() {
         checkTransformFromGeo("EPSG:27563", 3.005, 43.89, 653704.865208, 176887.660037);
     }
 
-//    @Test
-    public void NOTSUPPORTED_testGamma() {
+    @Test
+    public void testGamma() {
         // from Proj4.JS
-        checkTransformFromGeo("EPSG:2057", -53.0, 5.0, -1.160832226E7, 1.828261223E7, 0.1);
+    	checkTransformFromWGS84("EPSG:2057", -53.0, 5.0, -1.160832226E7, 1.828261223E7, 0.1);
     }
 
     @Test

diff --git a/src/test/java/org/locationtech/proj4j/Proj4JSTest.java b/src/test/java/org/locationtech/proj4j/Proj4JSTest.java
@@ -44,6 +44,7 @@ public void testGood() {
         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);
+        checkTransformFromGeo("EPSG:3375", 101.70979078430528, 3.06268465621428, 412597.532715, 338944.957259, 0.1);
     }
 
 //    @Test