OGRSpatialReference::importFromUSGS(): use plain WGS84 datum with WGS…

…84 ellipsoid, and identify EPSG code
rouault · Jun 17, 2024 · 94dc28c · 94dc28c
1 parent 3b62bf2
commit 94dc28c
Show file tree

Hide file tree

Showing 3 changed files with 41 additions and 17 deletions.
diff --git a/autotest/gdrivers/fast.py b/autotest/gdrivers/fast.py
@@ -31,7 +31,7 @@
 
 import gdaltest
 import pytest
-from osgeo import gdal
+from osgeo import gdal, osr
 
 ###############################################################################
 # Verify we have the driver.
@@ -210,20 +210,9 @@ def test_fast_7():
     gt = (676565.09, 5, 0, 5348341.5, 0, -5)
 
     # Expected definition of the projection
-    proj = """PROJCS["UTM Zone 32, Northern Hemisphere",
-    GEOGCS["Unknown datum based upon the WGS 84 ellipsoid",
-        DATUM["Not specified (based on WGS 84 spheroid)",
-            SPHEROID["WGS 84",6378137,298.257223563,
-                AUTHORITY["EPSG","7030"]]],
-        PRIMEM["Greenwich",0],
-        UNIT["degree",0.0174532925199433]],
-    PROJECTION["Transverse_Mercator"],
-    PARAMETER["latitude_of_origin",0],
-    PARAMETER["central_meridian",9],
-    PARAMETER["scale_factor",0.9996],
-    PARAMETER["false_easting",500000],
-    PARAMETER["false_northing",0],
-    UNIT["Meter",1]]"""
+    srs = osr.SpatialReference()
+    srs.ImportFromEPSG(32632)
+    proj = srs.ExportToWkt()
 
     return tst.testOpen(check_gt=gt, check_prj=proj)
 

diff --git a/autotest/osr/osr_usgs.py b/autotest/osr/osr_usgs.py
@@ -112,3 +112,21 @@ def test_osr_usgs_2():
         and gdal.PackedDMSToDec(params[4]) == pytest.approx(-117.4745429, abs=0.0000005)
         and gdal.PackedDMSToDec(params[5]) == pytest.approx(33.76446203, abs=0.0000005)
     ), "Can not import Lambert Conformal Conic projection."
+
+
+###############################################################################
+# Test the osr.SpatialReference.ImportFromUSGS() function with WGS 84
+#
+
+
+def test_osr_usgs_wgs84():
+
+    srs = osr.SpatialReference()
+    srs.ImportFromEPSG(32631)
+    (proj_code, zone, params, datum_code) = srs.ExportToUSGS()
+
+    srs2 = osr.SpatialReference()
+    srs2.ImportFromUSGS(proj_code, zone, params, datum_code)
+
+    assert srs2.IsSame(srs)
+    assert srs2.GetAuthorityCode(None) == "32631"
diff --git a/ogr/ogr_srs_usgs.cpp b/ogr/ogr_srs_usgs.cpp
@@ -734,8 +734,14 @@ OGRErr OGRSpatialReference::importFromUSGS(long iProjSys, long iZone,
         }
         else if (iDatum < NUMBER_OF_ELLIPSOIDS && aoEllips[iDatum])
         {
-            if (OSRGetEllipsoidInfo(aoEllips[iDatum], &pszName, &dfSemiMajor,
-                                    &dfInvFlattening) == OGRERR_NONE)
+            if (aoEllips[iDatum] == 7030)  // WGS 84 ellipsoid
+            {
+                // Assume a WGS 84 datum
+                SetWellKnownGeogCS("WGS84");
+            }
+            else if (OSRGetEllipsoidInfo(aoEllips[iDatum], &pszName,
+                                         &dfSemiMajor,
+                                         &dfInvFlattening) == OGRERR_NONE)
             {
                 SetGeogCS(
                     CPLString().Printf(
@@ -773,6 +779,17 @@ OGRErr OGRSpatialReference::importFromUSGS(long iProjSys, long iZone,
     if (IsLocal() || IsProjected())
         SetLinearUnits(SRS_UL_METER, 1.0);
 
+    if (iDatum < NUMBER_OF_ELLIPSOIDS && aoEllips[iDatum] == 7030)
+    {
+        if (AutoIdentifyEPSG() == OGRERR_NONE)
+        {
+            const char *pszAuthName = GetAuthorityName(nullptr);
+            const char *pszAuthCode = GetAuthorityCode(nullptr);
+            if (pszAuthName && pszAuthCode && EQUAL(pszAuthName, "EPSG"))
+                CPL_IGNORE_RET_VAL(importFromEPSG(atoi(pszAuthCode)));
+        }
+    }
+
     return OGRERR_NONE;
 }