Export ESRI:54099 WGS_1984_Spilhaus_Ocean_Map_in_Square as PROJ string,

and fix ESRI:54098 mapping to use +proj=spilhaus with +k_0=sqrt(2) Fixes #1851
OSGeo · Feb 18, 2025 · eb0aebc · eb0aebc
1 parent 05f6183
commit eb0aebc
Show file tree

Hide file tree

Showing 2 changed files with 18 additions and 16 deletions.
diff --git a/src/iso19111/operation/conversion.cpp b/src/iso19111/operation/conversion.cpp
@@ -4105,10 +4105,11 @@ void Conversion::_exportToPROJString(
         // Look for ESRI method and parameter names (to be opposed
         // to the OGC WKT2 names we use elsewhere, because there's no mapping
         // of those parameters to OGC WKT2)
-        // We also reject non-default values for a number of parameters,
-        // because they are not implemented on PROJ side. The subset we
-        // support can handle ESRI:54098 WGS_1984_Adams_Square_II, but not
+        // We at least support ESRI:54098 WGS_1984_Adams_Square_II and
         // ESRI:54099 WGS_1984_Spilhaus_Ocean_Map_in_Square
+        // More generally, we think our implementation of +proj=spilhaus
+        // matches ESRI Adams_Square_II with just a sqrt(2) factor difference
+        // for the scale factor
         const double falseEasting = parameterValueNumeric(
             "False_Easting", common::UnitOfMeasure::METRE);
         const double falseNorthing = parameterValueNumeric(
@@ -4126,14 +4127,13 @@ void Conversion::_exportToPROJString(
             "Latitude_Of_Center", common::UnitOfMeasure::DEGREE);
         const double XYPlaneRotation = parameterValueNumeric(
             "XY_Plane_Rotation", common::UnitOfMeasure::DEGREE);
-        if (scaleFactor != 1.0 || azimuth != 0.0 || latitudeOfCenter != 0.0 ||
-            XYPlaneRotation != 0.0) {
-            throw io::FormattingException("Unsupported value for one or "
-                                          "several parameters of "
-                                          "Adams_Square_II");
-        }
-        formatter->addStep("adams_ws2");
-        formatter->addParam("lon_0", longitudeOfCenter);
+
+        formatter->addStep("spilhaus");
+        formatter->addParam("lon_c", longitudeOfCenter);
+        formatter->addParam("lat_c", latitudeOfCenter);
+        formatter->addParam("k_0", M_SQRT2 * scaleFactor);
+        formatter->addParam("azi", azimuth);
+        formatter->addParam("rot", XYPlaneRotation);
         formatter->addParam("x_0", falseEasting);
         formatter->addParam("y_0", falseNorthing);
         bConversionDone = true;

diff --git a/test/unit/test_operation.cpp b/test/unit/test_operation.cpp
@@ -4360,19 +4360,21 @@ TEST(operation, adams_ws2_export) {
     auto crs = AuthorityFactory::create(dbContext, "ESRI")
                    ->createProjectedCRS("54098");
     EXPECT_EQ(crs->exportToPROJString(PROJStringFormatter::create().get()),
-              "+proj=adams_ws2 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m "
-              "+no_defs +type=crs");
+              "+proj=spilhaus +lon_c=0 +lat_c=0 +k_0=1.4142135623731 +azi=0 "
+              "+rot=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +type=crs");
 }
 
 // ---------------------------------------------------------------------------
 
-TEST(operation, adams_ws2_export_failure) {
+TEST(operation, spilhaus_esri_export) {
     auto dbContext = DatabaseContext::create();
     // ESRI:54099 WGS_1984_Spilhaus_Ocean_Map_in_Square
     auto crs = AuthorityFactory::create(dbContext, "ESRI")
                    ->createProjectedCRS("54099");
-    EXPECT_THROW(crs->exportToPROJString(PROJStringFormatter::create().get()),
-                 FormattingException);
+    EXPECT_EQ(crs->exportToPROJString(PROJStringFormatter::create().get()),
+              "+proj=spilhaus +lon_c=66.94970198 +lat_c=-49.56371678 "
+              "+k_0=1.4142135623731 +azi=40.17823482 +rot=45 "
+              "+x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +type=crs");
 }
 
 // ---------------------------------------------------------------------------