Repository: shapely/shapely (language: python)

Issue:
from_wkt and from_wkb crash on MULTISURFACE input
GEOS 3.13 introduced basic support for CURVE/SURFACE geometries, mainly to be able to parse them from WKB and WKT. However, `shapely.from_wkt` and `shapely.from_wkb` seem to crash on this type of input.

Noticed in a crashing test in pyogrio: https://github.com/geopandas/pyogrio/issues/478

Testing GEOS directly using `geosop` doesn't seem to trigger the problem (see issue above)?

I've had a quick look at the shapely code and don't immediately see where it can go wrong specifically on the curve types... but I'm not very familiar with the code, so probably missed something.

I also wonder what is going to happen if other shapely functions are going to be called on these geometries... but that's to be seen after this is solved :-)

Script to reproduce:

``` python
import shapely

wkt = "MULTISURFACE (CURVEPOLYGON ((-95.29698763468994 63.36578165250125, -94.5991423640996 53.06992070409323, -73.7073995758013 60.149472648492385, -95.29698763468994 63.36578165250125)))"
print(shapely.from_wkt(wkt))
print("here")
```


Discussion / hints:
Also fails for some of the WKBs in the GEOS test suite for nonlinear WKBs added in [geos #1106](https://github.com/libgeos/geos/pull/1106).

It looks like the actual read operation using GEOS is working properly and we're failing downstream of that.  We have a variety of checks against the type ID of the geometry throughout the codebase, and since these type IDs are effectively new to Shapely (GEOS <3.13 they were rejected by GEOS on read from WKT / WKB) it is possible that we are falling through a conditional statement somewhere and trying to do an unsupported operation on a nonlinear type.  Investigating now.

A short term fix may be to detect nonlinear types immediately after read from WKT / WKB and immediately convert them to their linear counterparts until better support for nonlinear types is added throughout Shapely.
The issue is that after creating the GEOS geometry object, we attempt to lookup the type in the shapely type registry, which is currently limited to non-linear types, and so we get an out-of-bounds indexing operation on that lookup.

The fix for that specific issue is to reject any geometry type IDs that are out of bounds.  This means that nonlinear types would still be unsupported in Shapely until someone has a chance to add comprehensive support, which is similar behavior to before but without the crash.  PR coming up.

Fix patch:
diff --git a/CHANGES.txt b/CHANGES.txt
index f6218524f..13553a011 100644
--- a/CHANGES.txt
+++ b/CHANGES.txt
@@ -9,6 +9,10 @@ Bug fixes:
 - Prevent crash when serializing a number > 1e100 to WKT with GEOS < 3.13. (#1907)
 - Fixes GeoJSON serialization of empty points (#2118)
 - Fixes `__geo_interface__` handling of empty points (#2120)
+- Fixes crash when reading nonlinear geometry types (CircularString,
+  CompoundCurve, MultiCurve, CurvePolygon, MultiSurface) from WKB/WKT with
+  GEOS >= 3.13; these types are not yet supported in Shapely and now raise a
+  ``NotImplementedError`` (#2160)
 
 Improvements:
 
diff --git a/src/pygeom.c b/src/pygeom.c
index f72e97da2..34f2eac13 100644
--- a/src/pygeom.c
+++ b/src/pygeom.c
@@ -23,6 +23,16 @@ PyObject* GeometryObject_FromGEOS(GEOSGeometry* ptr, GEOSContextHandle_t ctx) {
   if (type_id == -1) {
     return NULL;
   }
+
+  // Nonlinear types (CircularString, CompoundCurve, MultiCurve, CurvePolygon,
+  // MultiSurface are not currently supported
+  // TODO: this can be removed once these types are added to the type registry
+  if (type_id >= 8) {
+    PyErr_Format(PyExc_NotImplementedError,
+                 "Nonlinear geometry types are not currently supported");
+    return NULL;
+  }
+
   PyObject* type_obj = PyList_GET_ITEM(geom_registry[0], type_id);
   if (type_obj == NULL) {
     return NULL;


Test patch:
diff --git a/shapely/tests/test_io.py b/shapely/tests/test_io.py
index 2113d6ef4..f88e07b6e 100644
--- a/shapely/tests/test_io.py
+++ b/shapely/tests/test_io.py
@@ -239,6 +239,31 @@ def test_from_wkt_empty(wkt):
     assert shapely.to_wkt(geom) == wkt
 
 
+# WKT from https://github.com/libgeos/geos/blob/main/tests/unit/io/WKBReaderTest.cpp
+@pytest.mark.parametrize(
+    "wkt",
+    (
+        "CIRCULARSTRING(1 3,2 4,3 1)",
+        "COMPOUNDCURVE(CIRCULARSTRING(1 3,2 4,3 1),(3 1,0 0))",
+        "CURVEPOLYGON(COMPOUNDCURVE(CIRCULARSTRING(0 0,2 0,2 1,2 3,4 3),(4 3,4 5,1 4,0 0)),CIRCULARSTRING(1.7 1,1.4 0.4,1.6 0.4,1.6 0.5,1.7 1))",  # noqa: E501
+        "MULTICURVE((0 0,5 5),COMPOUNDCURVE((-1 -1,0 0),CIRCULARSTRING(0 0,1 1,2 0)),CIRCULARSTRING(4 0,4 4,8 4))",  # noqa: E501
+        "MULTISURFACE(CURVEPOLYGON(CIRCULARSTRING(0 0,4 0,4 4,0 4,0 0),(1 1,3 3,3 1,1 1)),((10 10,14 12,11 10,10 10),(11 11,11.5 11,11 11.5,11 11)))",  # noqa: E501
+    ),
+)
+def test_from_wkt_nonlinear_unsupported(wkt):
+    if shapely.geos_version >= (3, 13, 0):
+        with pytest.raises(
+            NotImplementedError,
+            match="Nonlinear geometry types are not currently supported",
+        ):
+            shapely.from_wkt(wkt)
+
+    else:
+        # prior to GEOS 3.13 nonlinear types were rejected by GEOS on read from WKT
+        with pytest.raises(shapely.errors.GEOSException, match="Unknown type"):
+            shapely.from_wkt(wkt)
+
+
 def test_from_wkb():
     expected = shapely.points(1, 1)
     actual = shapely.from_wkb(POINT11_WKB)
@@ -400,6 +425,36 @@ def test_from_wkb_empty(geom):
     assert shapely.to_wkb(geom) == wkb
 
 
+# WKB from https://github.com/libgeos/geos/blob/main/tests/unit/io/WKBReaderTest.cpp
+@pytest.mark.parametrize(
+    "wkb",
+    (
+        # "CIRCULARSTRING(1 3,2 4,3 1)",
+        "010800000003000000000000000000F03F0000000000000840000000000000004000000000000010400000000000000840000000000000F03F",
+        # "COMPOUNDCURVE(CIRCULARSTRING(1 3,2 4,3 1),(3 1,0 0))",
+        "01090000200E16000002000000010800000003000000000000000000F03F0000000000000840000000000000004000000000000010400000000000000840000000000000F03F0102000000020000000000000000000840000000000000F03F00000000000000000000000000000000",
+        # "CURVEPOLYGON(COMPOUNDCURVE(CIRCULARSTRING(0 0,2 0,2 1,2 3,4 3),(4 3,4 5,1 4,0 0)),CIRCULARSTRING(1.7 1,1.4 0.4,1.6 0.4,1.6 0.5,1.7 1))",  # noqa: E501
+        "010A0000200E1600000200000001090000000200000001080000000500000000000000000000000000000000000000000000000000004000000000000000000000000000000040000000000000F03F00000000000000400000000000000840000000000000104000000000000008400102000000040000000000000000001040000000000000084000000000000010400000000000001440000000000000F03F000000000000104000000000000000000000000000000000010800000005000000333333333333FB3F000000000000F03F666666666666F63F9A9999999999D93F9A9999999999F93F9A9999999999D93F9A9999999999F93F000000000000E03F333333333333FB3F000000000000F03F",
+        # "MULTICURVE((0 0,5 5),COMPOUNDCURVE((-1 -1,0 0),CIRCULARSTRING(0 0,1 1,2 0)),CIRCULARSTRING(4 0,4 4,8 4))",  # noqa: E501
+        "010B000000030000000102000000020000000000000000000000000000000000000000000000000014400000000000001440010900000002000000010200000002000000000000000000F0BF000000000000F0BF0000000000000000000000000000000001080000000300000000000000000000000000000000000000000000000000F03F000000000000F03F00000000000000400000000000000000010800000003000000000000000000104000000000000000000000000000001040000000000000104000000000000020400000000000001040",
+        # "MULTISURFACE(CURVEPOLYGON(CIRCULARSTRING(0 0,4 0,4 4,0 4,0 0),(1 1,3 3,3 1,1 1)),((10 10,14 12,11 10,10 10),(11 11,11.5 11,11 11.5,11 11)))",  # noqa: E501
+        "010C00000002000000010A000000020000000108000000050000000000000000000000000000000000000000000000000010400000000000000000000000000000104000000000000010400000000000000000000000000000104000000000000000000000000000000000010200000004000000000000000000F03F000000000000F03F000000000000084000000000000008400000000000000840000000000000F03F000000000000F03F000000000000F03F01030000000200000004000000000000000000244000000000000024400000000000002C40000000000000284000000000000026400000000000002440000000000000244000000000000024400400000000000000000026400000000000002640000000000000274000000000000026400000000000002640000000000000274000000000000026400000000000002640",
+    ),
+)
+def test_from_wkb_nonlinear_unsupported(wkb):
+    if shapely.geos_version >= (3, 13, 0):
+        with pytest.raises(
+            NotImplementedError,
+            match="Nonlinear geometry types are not currently supported",
+        ):
+            shapely.from_wkb(wkb)
+
+    else:
+        # prior to GEOS 3.13 nonlinear types were rejected by GEOS on read from WKB
+        with pytest.raises(shapely.errors.GEOSException, match="Unknown WKB type"):
+            shapely.from_wkb(wkb)
+
+
 def test_to_wkt():
     point = shapely.points(1, 1)
     actual = shapely.to_wkt(point)