Orthographic projection: do not add f=0 to PROJ string if the ellipsoid is a sphere (fixes GDAL PDS4 tests)

2 files changed, 37 insertions, 1 deletions

diff --git a/src/iso19111/coordinateoperation.cpp b/src/iso19111/coordinateoperation.cpp
index 184ff435..3577b6ff 100644
--- a/src/iso19111/coordinateoperation.cpp
+++ b/src/iso19111/coordinateoperation.cpp
@@ -6332,6 +6332,27 @@ void Conversion::_exportToPROJString(
                     }
                     bAxisSpecFound = true;
                 }
+
+                // No need to add explicit f=0 if the ellipsoid is a sphere
+                if (strcmp(mapping->proj_name_aux, "f=0") == 0) {
+                    crs::CRS *horiz = l_sourceCRS.get();
+                    const auto compound =
+                        dynamic_cast<const crs::CompoundCRS *>(horiz);
+                    if (compound) {
+                        const auto &components =
+                            compound->componentReferenceSystems();
+                        if (!components.empty()) {
+                            horiz = components.front().get();
+                        }
+                    }
+
+                    auto geogCRS =
+                        dynamic_cast<const crs::GeographicCRS *>(horiz);
+                    if (geogCRS && geogCRS->ellipsoid()->isSphere()) {
+                        addAux = false;
+                    }
+                }
+
                 if (addAux) {
                     auto kv = split(mapping->proj_name_aux, '=');
                     if (kv.size() == 2) {
diff --git a/test/unit/test_io.cpp b/test/unit/test_io.cpp
index 47e574c5..1be0f2aa 100644
--- a/test/unit/test_io.cpp
+++ b/test/unit/test_io.cpp
@@ -9272,7 +9272,7 @@ TEST(io, projparse_ortho_ellipsoidal) {
 
 // ---------------------------------------------------------------------------
 
-TEST(io, projparse_ortho_spherical) {
+TEST(io, projparse_ortho_spherical_on_ellipsoid) {
     std::string input("+proj=ortho +f=0 +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 "
                       "+ellps=WGS84 +units=m +no_defs +type=crs");
     auto obj = PROJStringParser().createFromPROJString(input);
@@ -9289,6 +9289,21 @@ TEST(io, projparse_ortho_spherical) {
 
 // ---------------------------------------------------------------------------
 
+TEST(io, projparse_ortho_spherical_on_sphere) {
+    std::string input("+proj=ortho +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 "
+                      "+R=6378137 +units=m +no_defs +type=crs");
+    auto obj = PROJStringParser().createFromPROJString(input);
+    auto crs = nn_dynamic_pointer_cast<ProjectedCRS>(obj);
+    ASSERT_TRUE(crs != nullptr);
+    EXPECT_EQ(
+        crs->exportToPROJString(
+            PROJStringFormatter::create(PROJStringFormatter::Convention::PROJ_4)
+                .get()),
+        input);
+}
+
+// ---------------------------------------------------------------------------
+
 TEST(io, projparse_axisswap_unitconvert_longlat_proj) {
     std::string input =
         "+type=crs +proj=pipeline +step +proj=axisswap +order=2,1 +step "