Add vectorized spherely.polygons() constructor

src/creation.cpp

@@ -13,6 +13,7 @@
 #include <s2geography/geography.h>
 
 #include <memory>
+#include <optional>
 #include <sstream>
 #include <stdexcept>
 #include <utility>
@@ -64,26 +65,23 @@ std::vector<S2Point> make_s2points(const std::vector<V>& points) {
     return std::move(s2points);
 }
 
-// create a S2Loop from coordinates or Point objects.
+void drop_closed_ring_duplicate(std::vector<S2Point>& s2points) {
+    if (!s2points.empty() && s2points.front() == s2points.back()) {
+        s2points.pop_back();
+    }
+}
+
+// create a S2Loop from a pre-built vector of S2Point. Callers are responsible
+// for dropping any closed-ring duplicate vertex before calling.
 //
 // Normalization (to CCW order for identifying the loop interior) and validation
 // are both enabled by default.
 //
-// Additional normalization is made here:
-// - if the input loop is already closed, remove one of the end nodes
-//
 // TODO: add option to skip normalization.
 //
-template <class V>
-std::unique_ptr<S2Loop> make_s2loop(const std::vector<V>& vertices,
-                                    bool check = true,
-                                    bool oriented = false) {
-    auto s2points = make_s2points(vertices);
-
-    if (s2points.front() == s2points.back()) {
-        s2points.pop_back();
-    }
-
+std::unique_ptr<S2Loop> make_s2loop_from_points(const std::vector<S2Point>& s2points,
+                                                bool check = true,
+                                                bool oriented = false) {
     auto loop_ptr = std::make_unique<S2Loop>();
 
     loop_ptr->set_s2debug_override(S2Debug::DISABLE);
@@ -101,13 +99,23 @@ std::unique_ptr<S2Loop> make_s2loop(const std::vector<V>& vertices,
         loop_ptr->Normalize();
     }
 
-    return std::move(loop_ptr);
+    return loop_ptr;
+}
+
+template <class V>
+std::unique_ptr<S2Loop> make_s2loop(const std::vector<V>& vertices,
+                                    bool check = true,
+                                    bool oriented = false) {
+    auto s2points = make_s2points(vertices);
+    drop_closed_ring_duplicate(s2points);
+    return make_s2loop_from_points(s2points, check, oriented);
 }
 
 // create a S2Polygon.
 //
 std::unique_ptr<S2Polygon> make_s2polygon(std::vector<std::unique_ptr<S2Loop>> loops,
-                                          bool oriented = false) {
+                                          bool oriented = false,
+                                          bool check = true) {
     auto polygon_ptr = std::make_unique<S2Polygon>();
     polygon_ptr->set_s2debug_override(S2Debug::DISABLE);
 
@@ -118,7 +126,7 @@ std::unique_ptr<S2Polygon> make_s2polygon(std::vector<std::unique_ptr<S2Loop>> l
     }
 
     // Note: this also checks each loop of the polygon
-    if (!polygon_ptr->IsValid()) {
+    if (check && !polygon_ptr->IsValid()) {
         std::stringstream err;
         S2Error s2err;
         err << "polygon is not valid: ";
@@ -159,6 +167,44 @@ py::array_t<PyObjectGeography> points(const py::array_t<double>& coords) {
     return points;
 }
 
+py::array_t<PyObjectGeography> polygons(const py::array_t<double>& shells,
+                                        bool oriented,
+                                        bool validate) {
+    auto shells_data = shells.unchecked<3>();
+    if (shells_data.shape(2) != 2) {
+        throw std::runtime_error("shells array must be of shape (N, K, 2)");
+    }
+
+    auto npolys = shells_data.shape(0);
+    auto nverts = shells_data.shape(1);
+
+    auto result = py::array_t<PyObjectGeography>(npolys);
+    py::buffer_info buf = result.request();
+    py::object* data = static_cast<py::object*>(buf.ptr);
+
+    // Scratch buffer reused across polygons. Per-ring validity is always
+    // skipped; bad rings still surface via the polygon-level IsValid() in
+    // make_s2polygon when validate=true.
+    std::vector<S2Point> pts;
+    pts.reserve(static_cast<size_t>(nverts));
+
+    for (py::ssize_t i = 0; i < npolys; i++) {
+        pts.clear();
+        for (py::ssize_t j = 0; j < nverts; j++) {
+            pts.push_back(make_s2point(shells_data(i, j, 0), shells_data(i, j, 1)));
+        }
+        drop_closed_ring_duplicate(pts);
+
+        std::vector<std::unique_ptr<S2Loop>> loops;
+        loops.push_back(make_s2loop_from_points(pts, /*check=*/false, oriented));
+
+        data[i] = make_py_geography<s2geog::PolygonGeography>(
+            make_s2polygon(std::move(loops), oriented, validate));
+    }
+
+    return result;
+}
+
 template <class V>
 std::unique_ptr<Geography> create_multipoint(const std::vector<V>& pts) {
     try {
@@ -513,4 +559,51 @@ void init_creation(py::module& m) {
             A array of longitude, latitude coordinate tuples (i.e., with shape (N, 2)).
 
     )pbdoc");
+
+    m.def("polygons",
+          &polygons,
+          py::arg("shells"),
+          py::arg("oriented") = false,
+          py::arg("validate") = true,
+          R"pbdoc(polygons(shells, oriented=False, validate=True)
+
+        Create an array of polygons from a numpy array of ring coordinates.
+
+        This vectorized constructor is functionally equivalent to calling
+        :py:func:`create_polygon` for each shell but avoids the per-polygon
+        Python parsing overhead, making it much faster when building many
+        polygons of uniform shape (e.g. grid cells).
+
+        Limitations compared to :py:func:`create_polygon`: all polygons
+        must have the same number of shell vertices (the numpy ndarray
+        input enforces this uniformity), and holes are not supported —
+        use :py:func:`create_polygon` in a Python loop if either varies
+        across the batch.
+
+        Parameters
+        ----------
+        shells : array_like
+            Array of shape ``(N, K, 2)`` giving ``N`` shell rings, each with
+            ``K`` vertices expressed as ``(longitude, latitude)`` in degrees.
+            Rings may be open (first vertex not repeated as last) or closed;
+            in the latter case the duplicate closing vertex is dropped.
+        oriented : bool, default False
+            Set to True if polygon ring directions are known to be correct
+            (i.e., shell ring vertices are defined counter clockwise).
+            By default (False), each ring is normalized so that the polygon
+            corresponds to the smaller area on the sphere.
+        validate : bool, default True
+            Validate each resulting polygon via ``S2Polygon::IsValid``. Set
+            to False only if the input polygons are known to be valid (e.g.
+            rectilinear grid cells); skipping the check is a meaningful
+            speedup on large batches but an invalid polygon will then be
+            constructed silently.
+
+        Returns
+        -------
+        polygons : ndarray
+            A 1-d object array of shape ``(N,)`` containing POLYGON
+            :class:`~spherely.Geography` objects.
+
+    )pbdoc");
 }

src/spherely.pyi

@@ -212,6 +212,11 @@ def create_collection(geographies: Iterable[Geography]) -> GeometryCollection: .
 def points(
     longitude: npt.ArrayLike, latitude: npt.ArrayLike
 ) -> PointGeography | T_NDArray_Geography: ...
+def polygons(
+    shells: npt.NDArray[np.float64],
+    oriented: bool = False,
+    validate: bool = True,
+) -> T_NDArray_Geography: ...
 
 # Geography utils
 

tests/test_creation.py

@@ -1,3 +1,4 @@
+import numpy as np
 import pytest
 
 import spherely
@@ -387,6 +388,97 @@ def test_multipolygon_invalid_geography() -> None:
         spherely.create_multipolygon([poly, line])
 
 
+def test_polygons_vectorized_basic() -> None:
+    shells = np.array(
+        [
+            [(0, 0), (2, 0), (2, 2), (0, 2)],
+            [(4, 0), (6, 0), (6, 2), (4, 2)],
+            [(10, -5), (12, -5), (12, -3), (10, -3)],
+            [(-30, 40), (-28, 40), (-28, 42), (-30, 42)],
+        ],
+        dtype=np.float64,
+    )
+    result = spherely.polygons(shells)
+
+    assert isinstance(result, np.ndarray)
+    assert result.shape == (4,)
+    assert result.dtype == object
+    assert repr(result[0]).startswith("POLYGON ((0 0")
+    assert repr(result[1]).startswith("POLYGON ((4 0")
+
+    # each vectorized polygon matches its scalar create_polygon equivalent
+    scalar = [spherely.create_polygon(ring.tolist()) for ring in shells]
+    assert [spherely.to_wkt(p) for p in result] == [spherely.to_wkt(p) for p in scalar]
+
+
+def test_polygons_vectorized_closed_ring() -> None:
+    # open and closed input rings should produce the same polygon
+    open_ring = np.array(
+        [[(0, 0), (2, 0), (2, 2), (0, 2)]],
+        dtype=np.float64,
+    )
+    closed_ring = np.array(
+        [[(0, 0), (2, 0), (2, 2), (0, 2), (0, 0)]],
+        dtype=np.float64,
+    )
+    a = spherely.polygons(open_ring)[0]
+    b = spherely.polygons(closed_ring)[0]
+    assert spherely.to_wkt(a) == spherely.to_wkt(b)
+
+
+def test_polygons_vectorized_oriented() -> None:
+    # same pattern as test_polygon_oriented but using the vectorized path
+    shells = np.array(
+        [[(0, 0), (0, 2), (2, 2), (2, 0)]],
+        dtype=np.float64,
+    )
+    (poly_cw,) = spherely.polygons(shells, oriented=True)
+
+    point = spherely.create_point(1, 1)
+    # CW polygon + oriented=True => point at (1, 1) is NOT in the interior
+    assert not spherely.contains(poly_cw, point)
+    assert repr(poly_cw) == "POLYGON ((0 0, 0 2, 2 2, 2 0, 0 0))"
+
+
+def test_polygons_vectorized_shape_errors() -> None:
+    # wrong trailing dimension
+    with pytest.raises(RuntimeError, match="shape"):
+        spherely.polygons(np.zeros((2, 4, 3), dtype=np.float64))
+
+    # wrong number of dimensions (surfaces from pybind11's unchecked view)
+    with pytest.raises(ValueError, match="number of dimensions"):
+        spherely.polygons(np.zeros((4, 2), dtype=np.float64))
+
+
+def test_polygons_vectorized_invalid_ring() -> None:
+    # same validation behavior as scalar create_polygon: bad rings surface as
+    # "polygon is not valid" from the final polygon-level check.
+    shells = np.array(
+        [
+            [(0, 0), (2, 0), (2, 2), (0, 2)],
+            [(0, 0), (0, 2), (0, 2), (2, 0)],  # invalid: duplicate + self-cross
+        ],
+        dtype=np.float64,
+    )
+    with pytest.raises(ValueError, match="polygon is not valid"):
+        spherely.polygons(shells)
+
+
+def test_polygons_vectorized_validate_false_skips_validation() -> None:
+    # same bad shell as above, but validate=False must not raise.
+    shells = np.array(
+        [
+            [(0, 0), (2, 0), (2, 2), (0, 2)],
+            [(0, 0), (0, 2), (0, 2), (2, 0)],
+        ],
+        dtype=np.float64,
+    )
+    result = spherely.polygons(shells, validate=False)
+    # we still get a Geography back for each input; no guarantee it's valid.
+    assert result.shape == (2,)
+    assert all(isinstance(g, spherely.Geography) for g in result)
+
+
 def test_collection() -> None:
     objs = [
         spherely.create_point(0, 0),