From 55e1428ea624bbf63db856ec261dd5926dc7d076 Mon Sep 17 00:00:00 2001
From: Even Rouault <even.rouault@spatialys.com>
Date: Thu, 9 Jan 2025 19:33:26 +0100
Subject: [PATCH] GDALGCPsToGeoTransform(): return FALSE when invalid
 geotransform is generated

Fixes case of https://github.com/OSGeo/gdal/issues/11618
---
 autotest/gcore/gcps2geotransform.py | 31 +++++++++++++++++++++++++++++
 gcore/gdal_misc.cpp                 |  7 +++++++
 2 files changed, 38 insertions(+)

diff --git a/autotest/gcore/gcps2geotransform.py b/autotest/gcore/gcps2geotransform.py
index db803924d632..70c8971e6fd7 100755
--- a/autotest/gcore/gcps2geotransform.py
+++ b/autotest/gcore/gcps2geotransform.py
@@ -195,3 +195,34 @@ def test_gcps2gt_8():
         1.596921026218e-05,
     )
     gdaltest.check_geotransform(gt, gt_expected, 0.00001)
+
+
+###############################################################################
+# Test case of https://github.com/OSGeo/gdal/issues/11618
+
+
+def test_gcps2gt_broken_hour_glass():
+
+    gt = gdal.GCPsToGeoTransform(
+        _list2gcps(
+            [
+                (0, 0, 0, 0),
+                (0, 10, 0, 10),
+                (10, 0, 10, 10),
+                (10, 10, 10, 0),
+            ]
+        )
+    )
+    assert gt is None
+
+    gt = gdal.GCPsToGeoTransform(
+        _list2gcps(
+            [
+                (0, 0, 0, 0),
+                (0, 10, 10, 0),
+                (10, 0, 0, 10),
+                (10, 10, 10, 10),
+            ]
+        )
+    )
+    assert gt is None
diff --git a/gcore/gdal_misc.cpp b/gcore/gdal_misc.cpp
index c86e42a7eaa8..bcd9a3930837 100644
--- a/gcore/gdal_misc.cpp
+++ b/gcore/gdal_misc.cpp
@@ -3265,6 +3265,13 @@ int CPL_STDCALL GDALGCPsToGeoTransform(int nGCPCount, const GDAL_GCP *pasGCPs,
     GDALComposeGeoTransforms(pl_normalize, gt_normalized, gt1p2);
     GDALComposeGeoTransforms(gt1p2, inv_geo_normalize, padfGeoTransform);
 
+    // "Hour-glass" like shape of GCPs. Cf https://github.com/OSGeo/gdal/issues/11618
+    if (std::abs(padfGeoTransform[1]) <= 1e-15 ||
+        std::abs(padfGeoTransform[5]) <= 1e-15)
+    {
+        return FALSE;
+    }
+
     /* -------------------------------------------------------------------- */
     /*      Now check if any of the input points fit this poorly.           */
     /* -------------------------------------------------------------------- */