Fix support of WKT1_GDAL with netCDF rotated pole formulation (#2185)

Contributes to fixing issue raised in
https://lists.osgeo.org/pipermail/gdal-dev/2020-April/052003.html

3 files changed, 83 insertions, 4 deletions

diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp
index 4ae5f8a6..ec515cd7 100644
--- a/src/iso19111/coordinateoperation.cpp
+++ b/src/iso19111/coordinateoperation.cpp
@@ -6380,10 +6380,9 @@ void Conversion::_exportToPROJString(
 
         auto derivedGeographicCRS =
             dynamic_cast<const crs::DerivedGeographicCRS *>(horiz);
-        if (derivedGeographicCRS) {
-            auto baseGeodCRS = derivedGeographicCRS->baseCRS();
+        if (!formatter->getCRSExport() && derivedGeographicCRS) {
             formatter->setOmitProjLongLatIfPossible(true);
-            baseGeodCRS->_exportToPROJString(formatter);
+            derivedGeographicCRS->addAngularUnitConvertAndAxisSwap(formatter);
             formatter->setOmitProjLongLatIfPossible(false);
         }
     }
diff --git a/src/iso19111/io.cpp b/src/iso19111/io.cpp
index 260cf1c7..86f21db7 100644
--- a/src/iso19111/io.cpp
+++ b/src/iso19111/io.cpp
@@ -4514,7 +4514,28 @@ CRSPtr WKTParser::Private::buildCRS(const WKTNodeNNPtr &node) {
     if (ci_equal(name, WKTConstants::PROJCS) ||
         ci_equal(name, WKTConstants::PROJCRS) ||
         ci_equal(name, WKTConstants::PROJECTEDCRS)) {
-        return util::nn_static_pointer_cast<CRS>(buildProjectedCRS(node));
+        auto projCRS =
+            util::nn_static_pointer_cast<CRS>(buildProjectedCRS(node));
+        auto projString = projCRS->getExtensionProj4();
+        if (starts_with(projString, "+proj=ob_tran +o_proj=longlat") ||
+            starts_with(projString, "+proj=ob_tran +o_proj=lonlat") ||
+            starts_with(projString, "+proj=ob_tran +o_proj=latlong") ||
+            starts_with(projString, "+proj=ob_tran +o_proj=latlon")) {
+            // Those are not a projected CRS, but a DerivedGeographic one...
+            if (projString.find(" +type=crs") == std::string::npos) {
+                projString += " +type=crs";
+            }
+            try {
+                auto projObj =
+                    PROJStringParser().createFromPROJString(projString);
+                auto crs = nn_dynamic_pointer_cast<CRS>(projObj);
+                if (crs) {
+                    return crs;
+                }
+            } catch (const io::ParsingException &) {
+            }
+        }
+        return projCRS.as_nullable();
     }
 
     if (ci_equal(name, WKTConstants::VERT_CS) ||
diff --git a/test/unit/test_io.cpp b/test/unit/test_io.cpp
index 1c69bb05..2a262111 100644
--- a/test/unit/test_io.cpp
+++ b/test/unit/test_io.cpp
@@ -3625,6 +3625,65 @@ TEST(wkt_parse, DerivedGeodeticCRS) {
 
 // ---------------------------------------------------------------------------
 
+TEST(wkt_parse, DerivedGeographicCRS_GDAL_PROJ4_EXSTENSION_hack) {
+    auto wkt =
+        "PROJCS[\"unnamed\","
+        "   GEOGCS[\"unknown\","
+        "       DATUM[\"unnamed\","
+        "           SPHEROID[\"Spheroid\",6367470,594.313048347956]],"
+        "       PRIMEM[\"Greenwich\",0],"
+        "       UNIT[\"degree\",0.0174532925199433,"
+        "           AUTHORITY[\"EPSG\",\"9122\"]]],"
+        "   PROJECTION[\"Rotated_pole\"],"
+        "       UNIT[\"metre\",1,AUTHORITY[\"EPSG\",\"9001\"]],"
+        "       AXIS[\"Easting\",EAST],"
+        "       AXIS[\"Northing\",NORTH],"
+        "       EXTENSION[\"PROJ4\",\"+proj=ob_tran +o_proj=longlat +lon_0=18 "
+        "+o_lon_p=0 +o_lat_p=39.25 +a=6367470 +b=6367470 "
+        "+to_meter=0.0174532925199 +wktext\"]]";
+
+    auto obj = WKTParser().createFromWKT(wkt);
+    auto crs = nn_dynamic_pointer_cast<DerivedGeographicCRS>(obj);
+    ASSERT_TRUE(crs != nullptr);
+
+    auto obj2 = PROJStringParser().createFromPROJString(
+        "+proj=ob_tran +o_proj=longlat +lon_0=18 "
+        "+o_lon_p=0 +o_lat_p=39.25 +a=6367470 +b=6367470 "
+        "+to_meter=0.0174532925199 +wktext +type=crs");
+    auto crs2 = nn_dynamic_pointer_cast<DerivedGeographicCRS>(obj2);
+    ASSERT_TRUE(crs2 != nullptr);
+
+    EXPECT_TRUE(
+        crs->isEquivalentTo(crs2.get(), IComparable::Criterion::EQUIVALENT));
+
+    {
+        auto op = CoordinateOperationFactory::create()->createOperation(
+            crs->baseCRS(), NN_NO_CHECK(crs));
+        ASSERT_TRUE(op != nullptr);
+        EXPECT_EQ(op->exportToPROJString(PROJStringFormatter::create().get()),
+                  "+proj=pipeline "
+                  "+step +proj=unitconvert +xy_in=deg +xy_out=rad "
+                  "+step +proj=ob_tran +o_proj=longlat +lon_0=18 +o_lon_p=0 "
+                  "+o_lat_p=39.25 +R=6367470 "
+                  "+step +proj=unitconvert +xy_in=rad +xy_out=deg");
+    }
+
+    {
+        auto op = CoordinateOperationFactory::create()->createOperation(
+            NN_NO_CHECK(crs), crs->baseCRS());
+        ASSERT_TRUE(op != nullptr);
+        EXPECT_EQ(
+            op->exportToPROJString(PROJStringFormatter::create().get()),
+            "+proj=pipeline "
+            "+step +proj=unitconvert +xy_in=deg +xy_out=rad "
+            "+step +inv +proj=ob_tran +o_proj=longlat +lon_0=18 +o_lon_p=0 "
+            "+o_lat_p=39.25 +R=6367470 "
+            "+step +proj=unitconvert +xy_in=rad +xy_out=deg");
+    }
+}
+
+// ---------------------------------------------------------------------------
+
 TEST(wkt_parse, DerivedProjectedCRS) {
     auto wkt =
         "DERIVEDPROJCRS[\"derived projectedCRS\",\n"