Merge pull request #2389 from rouault/wkt_parse_implict_compound_crs_esri

WKT1_ESRI: fix import and export of CompoundCRS

3 files changed, 161 insertions, 10 deletions

diff --git a/test/cli/testprojinfo b/test/cli/testprojinfo
index 19d96735..043e5724 100755
--- a/test/cli/testprojinfo
+++ b/test/cli/testprojinfo
@@ -175,6 +175,14 @@ echo 'Testing -k operation EPSG:8457 -o PROJ -q' >> ${OUT}
 $EXE -k operation EPSG:8457 -o PROJ -q >>${OUT} 2>&1
 echo "" >>${OUT}
 
+echo 'Testing D_WGS_1984' >> ${OUT}
+$EXE D_WGS_1984 >>${OUT} 2>&1
+echo "" >>${OUT}
+
+echo 'Testing -k datum D_WGS_1984' >> ${OUT}
+$EXE -k datum  D_WGS_1984 >>${OUT} 2>&1
+echo "" >>${OUT}
+
 # Testing --searchpaths
 if ! $EXE --searchpaths > testprojinfo_out_searchpaths.txt; then
     echo "--searchpaths failed"
diff --git a/test/cli/testprojinfo_out.dist b/test/cli/testprojinfo_out.dist
index b4d19002..3025295b 100644
--- a/test/cli/testprojinfo_out.dist
+++ b/test/cli/testprojinfo_out.dist
@@ -1191,6 +1191,20 @@ Testing -k operation EPSG:8457 -o PROJ -q
   +step +proj=unitconvert +xy_in=rad +xy_out=deg
   +step +proj=axisswap +order=2,1
 
+Testing D_WGS_1984
+WKT2:2019 string:
+DATUM["World Geodetic System 1984",
+    ELLIPSOID["WGS 84",6378137,298.257223563,
+        LENGTHUNIT["metre",1]],
+    ID["EPSG",6326]]
+
+Testing -k datum D_WGS_1984
+WKT2:2019 string:
+DATUM["World Geodetic System 1984",
+    ELLIPSOID["WGS 84",6378137,298.257223563,
+        LENGTHUNIT["metre",1]],
+    ID["EPSG",6326]]
+
 Testing -s NZGD2000 -t ITRF96 -o PROJ -q
 +proj=pipeline
   +step +proj=unitconvert +xy_in=deg +xy_out=rad
diff --git a/test/unit/test_io.cpp b/test/unit/test_io.cpp
index 055b1e1d..555d1159 100644
--- a/test/unit/test_io.cpp
+++ b/test/unit/test_io.cpp
@@ -2597,11 +2597,12 @@ TEST(wkt_parse,
     EXPECT_EQ(baseCRS->nameStr(), "NAD83");
     EXPECT_EQ(baseCRS->coordinateSystem()->axisList().size(), 3U);
 
-    EXPECT_EQ(
-        crs->exportToWKT(
-            WKTFormatter::create(WKTFormatter::Convention::WKT1_GDAL, dbContext)
-                .get()),
-        wkt);
+    EXPECT_EQ(replaceAll(crs->exportToWKT(
+                             WKTFormatter::create(
+                                 WKTFormatter::Convention::WKT1_GDAL, dbContext)
+                                 .get()),
+                         "ellipsoidal height", "Up"),
+              wkt);
 }
 
 // ---------------------------------------------------------------------------
@@ -2645,11 +2646,12 @@ TEST(wkt_parse,
     EXPECT_EQ(crs->nameStr(), "WGS 84 / UTM zone 31N");
     EXPECT_EQ(crs->coordinateSystem()->axisList().size(), 3U);
 
-    EXPECT_EQ(
-        crs->exportToWKT(
-            WKTFormatter::create(WKTFormatter::Convention::WKT1_GDAL, dbContext)
-                .get()),
-        wkt);
+    EXPECT_EQ(replaceAll(crs->exportToWKT(
+                             WKTFormatter::create(
+                                 WKTFormatter::Convention::WKT1_GDAL, dbContext)
+                                 .get()),
+                         "ellipsoidal height", "Up"),
+              wkt);
 }
 
 // ---------------------------------------------------------------------------
@@ -2683,10 +2685,137 @@ TEST(wkt_parse,
     EXPECT_NEAR(crs->coordinateSystem()->axisList()[2]->unit().conversionToSI(),
                 0.304800609601219, 1e-15);
 
+    EXPECT_EQ(replaceAll(crs->exportToWKT(
+                             WKTFormatter::create(
+                                 WKTFormatter::Convention::WKT1_GDAL, dbContext)
+                                 .get()),
+                         "ellipsoidal height", "Up"),
+              wkt);
+}
+
+// ---------------------------------------------------------------------------
+
+TEST(wkt_parse, implicit_compound_CRS_ESRI) {
+    // See https://lists.osgeo.org/pipermail/gdal-dev/2020-October/052843.html
+    // and https://pro.arcgis.com/en/pro-app/arcpy/classes/spatialreference.htm
+    const char *wkt =
+        "PROJCS[\"NAD_1983_2011_StatePlane_Colorado_Central_FIPS_0502_Ft_US\","
+        "GEOGCS[\"GCS_NAD_1983_2011\",DATUM[\"D_NAD_1983_2011\","
+        "SPHEROID[\"GRS_1980\",6378137.0,298.257222101]],"
+        "PRIMEM[\"Greenwich\",0.0],"
+        "UNIT[\"Degree\",0.0174532925199433]],"
+        "PROJECTION[\"Lambert_Conformal_Conic\"],"
+        "PARAMETER[\"False_Easting\",3000000.00031608],"
+        "PARAMETER[\"False_Northing\",999999.999996],"
+        "PARAMETER[\"Central_Meridian\",-105.5],"
+        "PARAMETER[\"Standard_Parallel_1\",38.45],"
+        "PARAMETER[\"Standard_Parallel_2\",39.75],"
+        "PARAMETER[\"Latitude_Of_Origin\",37.8333333333333],"
+        "UNIT[\"US survey foot\",0.304800609601219]],"
+        "VERTCS[\"CGVD2013_height\","
+        "VDATUM[\"Canadian_Geodetic_Vertical_Datum_of_2013\"],"
+        "PARAMETER[\"Vertical_Shift\",0.0],"
+        "PARAMETER[\"Direction\",1.0],"
+        "UNIT[\"Meter\",1.0]]";
+    auto dbContext = DatabaseContext::create();
+    auto obj = WKTParser().attachDatabaseContext(dbContext).createFromWKT(wkt);
+    auto crs = nn_dynamic_pointer_cast<CompoundCRS>(obj);
+    ASSERT_TRUE(crs != nullptr);
+    EXPECT_EQ(crs->nameStr(), "NAD83(2011) / Colorado Central (ftUS) + "
+                              "CGVD2013(CGG2013) height");
+
+    EXPECT_EQ(
+        crs->exportToWKT(
+            WKTFormatter::create(WKTFormatter::Convention::WKT1_ESRI, dbContext)
+                .get()),
+        wkt);
+}
+
+// ---------------------------------------------------------------------------
+
+TEST(wkt_parse, VERTCS_with_ellipsoidal_height_ESRI) {
+    const char *wkt = "VERTCS[\"WGS_1984\",DATUM[\"D_WGS_1984\","
+                      "SPHEROID[\"WGS_1984\",6378137.0,298.257223563]],"
+                      "PARAMETER[\"Vertical_Shift\",0.0],"
+                      "PARAMETER[\"Direction\",1.0],UNIT[\"Meter\",1.0]]";
+    auto dbContext = DatabaseContext::create();
+    auto obj = WKTParser().attachDatabaseContext(dbContext).createFromWKT(wkt);
+    auto crs = nn_dynamic_pointer_cast<VerticalCRS>(obj);
+    ASSERT_TRUE(crs != nullptr);
+
+    EXPECT_EQ(
+        crs->exportToWKT(
+            WKTFormatter::create(WKTFormatter::Convention::WKT1_ESRI, dbContext)
+                .get()),
+        wkt);
+
+    const char *expected_wkt1 =
+        "VERT_CS[\"WGS_1984\",\n"
+        "    VERT_DATUM[\"World Geodetic System 1984\",2002],\n"
+        "    UNIT[\"metre\",1,\n"
+        "        AUTHORITY[\"EPSG\",\"9001\"]],\n"
+        "    AXIS[\"ellipsoidal height\",UP]]";
     EXPECT_EQ(
         crs->exportToWKT(
             WKTFormatter::create(WKTFormatter::Convention::WKT1_GDAL, dbContext)
                 .get()),
+        expected_wkt1);
+}
+
+// ---------------------------------------------------------------------------
+
+TEST(wkt_parse, implicit_compound_CRS_geographic_with_ellipsoidal_height_ESRI) {
+    const char *wkt =
+        "GEOGCS[\"GCS_WGS_1984\",DATUM[\"D_WGS_1984\","
+        "SPHEROID[\"WGS_1984\",6378137.0,298.257223563]],"
+        "PRIMEM[\"Greenwich\",0.0],UNIT[\"Degree\",0.0174532925199433]],"
+        "VERTCS[\"WGS_1984\",DATUM[\"D_WGS_1984\","
+        "SPHEROID[\"WGS_1984\",6378137.0,298.257223563]],"
+        "PARAMETER[\"Vertical_Shift\",0.0],"
+        "PARAMETER[\"Direction\",1.0],UNIT[\"Meter\",1.0]]";
+    auto dbContext = DatabaseContext::create();
+    auto obj = WKTParser().attachDatabaseContext(dbContext).createFromWKT(wkt);
+    auto crs = nn_dynamic_pointer_cast<GeographicCRS>(obj);
+    ASSERT_TRUE(crs != nullptr);
+    EXPECT_EQ(crs->coordinateSystem()->axisList().size(), 3U);
+
+    EXPECT_EQ(
+        crs->exportToWKT(
+            WKTFormatter::create(WKTFormatter::Convention::WKT1_ESRI, dbContext)
+                .get()),
+        wkt);
+}
+
+// ---------------------------------------------------------------------------
+
+TEST(wkt_parse, implicit_compound_CRS_projected_with_ellipsoidal_height_ESRI) {
+    const char *wkt =
+        "PROJCS[\"WGS_1984_UTM_Zone_31N\",GEOGCS[\"GCS_WGS_1984\","
+        "DATUM[\"D_WGS_1984\","
+        "SPHEROID[\"WGS_1984\",6378137.0,298.257223563]],"
+        "PRIMEM[\"Greenwich\",0.0],UNIT[\"Degree\",0.0174532925199433]],"
+        "PROJECTION[\"Transverse_Mercator\"],"
+        "PARAMETER[\"False_Easting\",500000.0],"
+        "PARAMETER[\"False_Northing\",0.0],"
+        "PARAMETER[\"Central_Meridian\",3.0],"
+        "PARAMETER[\"Scale_Factor\",0.9996],"
+        "PARAMETER[\"Latitude_Of_Origin\",0.0],"
+        "UNIT[\"Meter\",1.0]],"
+        "VERTCS[\"WGS_1984\","
+        "DATUM[\"D_WGS_1984\",SPHEROID[\"WGS_1984\",6378137.0,298.257223563]],"
+        "PARAMETER[\"Vertical_Shift\",0.0],"
+        "PARAMETER[\"Direction\",1.0],"
+        "UNIT[\"Meter\",1.0]]";
+    auto dbContext = DatabaseContext::create();
+    auto obj = WKTParser().attachDatabaseContext(dbContext).createFromWKT(wkt);
+    auto crs = nn_dynamic_pointer_cast<ProjectedCRS>(obj);
+    ASSERT_TRUE(crs != nullptr);
+    EXPECT_EQ(crs->coordinateSystem()->axisList().size(), 3U);
+
+    EXPECT_EQ(
+        crs->exportToWKT(
+            WKTFormatter::create(WKTFormatter::Convention::WKT1_ESRI, dbContext)
+                .get()),
         wkt);
 }