structure saas with tools

.venv/lib/python3.10/site-packages/shapely/__init__.py

@@ -0,0 +1,36 @@
+"""Manipulation and analysis of geometric objects in the Cartesian plane."""
+
+from shapely.lib import GEOSException
+from shapely.lib import Geometry
+from shapely.lib import geos_version, geos_version_string
+from shapely.lib import geos_capi_version, geos_capi_version_string
+from shapely.errors import setup_signal_checks
+from shapely._geometry import *
+from shapely.creation import *
+from shapely.constructive import *
+from shapely.predicates import *
+from shapely.measurement import *
+from shapely.set_operations import *
+from shapely.linear import *
+from shapely.coordinates import *
+from shapely.strtree import *
+from shapely.io import *
+from shapely._coverage import *
+
+# Submodule always needs to be imported to ensure Geometry subclasses are registered
+from shapely.geometry import (
+    Point,
+    LineString,
+    Polygon,
+    MultiPoint,
+    MultiLineString,
+    MultiPolygon,
+    GeometryCollection,
+    LinearRing,
+)
+
+from shapely import _version
+
+__version__ = _version.get_versions()["version"]
+
+setup_signal_checks()

.venv/lib/python3.10/site-packages/shapely/_coverage.py

@@ -0,0 +1,168 @@
+import numpy as np
+
+from shapely import Geometry, GeometryType, lib
+from shapely._geometry import get_parts
+from shapely.decorators import multithreading_enabled, requires_geos
+
+__all__ = ["coverage_invalid_edges", "coverage_is_valid", "coverage_simplify"]
+
+
+@requires_geos("3.12.0")
+@multithreading_enabled
+def coverage_is_valid(geometry, gap_width=0.0, **kwargs):
+    """Verify if a coverage is valid.
+
+    The coverage is represented by an array of polygonal geometries with
+    exactly matching edges and no overlap.
+
+    A valid coverage may contain holes (regions of no coverage). However,
+    sometimes it might be desirable to detect narrow gaps as invalidities in
+    the coverage. The `gap_width` parameter allows to specify the maximum
+    width of gaps to detect. When gaps are detected, this function will
+    return False and the `coverage_invalid_edges` function can be used to
+    find the edges of those gaps.
+
+    Geometries that are not Polygon or MultiPolygon are ignored.
+
+    .. versionadded:: 2.1.0
+
+    Parameters
+    ----------
+    geometry : array_like
+        Array of geometries to verify.
+    gap_width : float, default 0.0
+        The maximum width of gaps to detect.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Returns
+    -------
+    bool
+
+    See Also
+    --------
+    coverage_invalid_edges, coverage_simplify
+
+    """
+    geometries = np.asarray(geometry)
+    # we always consider the full array as a single coverage -> ravel the input
+    # to pass a 1D array
+    return lib.coverage_is_valid(geometries.ravel(order="K"), gap_width, **kwargs)
+
+
+@requires_geos("3.12.0")
+@multithreading_enabled
+def coverage_invalid_edges(geometry, gap_width=0.0, **kwargs):
+    """Verify if a coverage is valid and return invalid edges.
+
+    This functions returns linear indicators showing the location of invalid
+    edges (if any) in each polygon in the input array.
+
+    The coverage is represented by an array of polygonal geometries with
+    exactly matching edges and no overlap.
+
+    A valid coverage may contain holes (regions of no coverage). However,
+    sometimes it might be desirable to detect narrow gaps as invalidities in
+    the coverage. The `gap_width` parameter allows to specify the maximum
+    width of gaps to detect. When gaps are detected, the `coverage_is_valid`
+    function will return False and this function can be used to find the
+    edges of those gaps.
+
+    Geometries that are not Polygon or MultiPolygon are ignored.
+
+    .. versionadded:: 2.1.0
+
+    Parameters
+    ----------
+    geometry : array_like
+        Array of geometries to verify.
+    gap_width : float, default 0.0
+        The maximum width of gaps to detect.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Returns
+    -------
+    numpy.ndarray | shapely.Geometry
+
+    See Also
+    --------
+    coverage_is_valid, coverage_simplify
+
+    """
+    geometries = np.asarray(geometry)
+    # we always consider the full array as a single coverage -> ravel the input
+    # to pass a 1D array
+    return lib.coverage_invalid_edges(geometries.ravel(order="K"), gap_width, **kwargs)
+
+
+@requires_geos("3.12.0")
+@multithreading_enabled
+def coverage_simplify(geometry, tolerance, *, simplify_boundary=True):
+    """Return a simplified version of an input geometry using coverage simplification.
+
+    Assumes that the geometry forms a polygonal coverage. Under this assumption, the
+    function simplifies the edges using the Visvalingam-Whyatt algorithm, while
+    preserving a valid coverage. In the most simplified case, polygons are reduced to
+    triangles.
+
+    A collection of valid polygons is considered a coverage if the polygons are:
+
+    * **Non-overlapping** - polygons do not overlap (their interiors do not intersect)
+    * **Edge-Matched** - vertices along shared edges are identical
+
+    The function allows simplification of all edges including the outer boundaries of
+    the coverage or simplification of only the inner (shared) edges.
+
+    If there are other geometry types than Polygons or MultiPolygons present,
+    the function will raise an error.
+
+    If the geometry is polygonal but does not form a valid coverage due to overlaps,
+    it will be simplified but it may result in invalid topology.
+
+    .. versionadded:: 2.1.0
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+    tolerance : float or array_like
+        The degree of simplification roughly equal to the square root of the area
+        of triangles that will be removed.
+    simplify_boundary : bool, optional
+        By default (True), simplifies both internal edges of the coverage as well
+        as its boundary. If set to False, only simplifies internal edges.
+
+    Returns
+    -------
+    numpy.ndarray | shapely.Geometry
+
+    See Also
+    --------
+    coverage_is_valid, coverage_invalid_edges
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Polygon
+    >>> poly = Polygon([(0, 0), (20, 0), (20, 10), (10, 5), (0, 10), (0, 0)])
+    >>> shapely.coverage_simplify(poly, tolerance=2)
+    <POLYGON ((0 0, 20 0, 20 10, 10 5, 0 10, 0 0))>
+    """
+    scalar = False
+    if isinstance(geometry, Geometry):
+        scalar = True
+
+    geometries = np.asarray(geometry)
+    shape = geometries.shape
+    geometries = geometries.ravel()
+
+    # create_collection acts on the inner axis
+    collections = lib.create_collection(
+        geometries, np.intc(GeometryType.GEOMETRYCOLLECTION)
+    )
+
+    simplified = lib.coverage_simplify(collections, tolerance, simplify_boundary)
+    parts = get_parts(simplified).reshape(shape)
+    if scalar:
+        return parts.item()
+    return parts

.venv/lib/python3.10/site-packages/shapely/_enum.py

@@ -0,0 +1,23 @@
+from enum import IntEnum
+
+
+class ParamEnum(IntEnum):
+    """Wraps IntEnum to provide validation of a requested item.
+
+    Intended for enums used for function parameters.
+
+    Use enum.get_value(item) for this behavior instead of builtin enum[item].
+    """
+
+    @classmethod
+    def get_value(cls, item):
+        """Validate item and raise a ValueError with valid options if not present."""
+        try:
+            return cls[item].value
+        except KeyError:
+            valid_options = {e.name for e in cls}
+            raise ValueError(
+                "'{}' is not a valid option, must be one of '{}'".format(
+                    item, "', '".join(valid_options)
+                )
+            )

.venv/lib/python3.10/site-packages/shapely/_geos.pxd

@@ -0,0 +1,57 @@
+"""Provides a wrapper for GEOS types and functions.
+
+Note: GEOS functions in Cython must be called using the get_geos_handle context
+manager.
+
+Examples
+--------
+    with get_geos_handle() as geos_handle:
+        SomeGEOSFunc(geos_handle, ...<other params>)
+
+"""
+
+cdef extern from "geos_c.h":
+    # Types
+    ctypedef void *GEOSContextHandle_t
+    ctypedef struct GEOSGeometry
+    ctypedef struct GEOSCoordSequence
+    ctypedef void (*GEOSMessageHandler_r)(const char *message, void *userdata)
+
+    # GEOS Context & Messaging
+    GEOSContextHandle_t GEOS_init_r() nogil
+    void GEOS_finish_r(GEOSContextHandle_t handle) nogil
+    void GEOSContext_setErrorMessageHandler_r(GEOSContextHandle_t handle, GEOSMessageHandler_r ef, void* userData) nogil
+    void GEOSContext_setNoticeMessageHandler_r(GEOSContextHandle_t handle, GEOSMessageHandler_r nf, void* userData) nogil
+
+    # Geometry functions
+    const GEOSGeometry* GEOSGetGeometryN_r(GEOSContextHandle_t handle, const GEOSGeometry* g, int n) nogil
+    const GEOSGeometry* GEOSGetExteriorRing_r(GEOSContextHandle_t handle, const GEOSGeometry* g) nogil
+    const GEOSGeometry* GEOSGetInteriorRingN_r(GEOSContextHandle_t handle, const GEOSGeometry* g, int n) nogil
+    int GEOSGeomTypeId_r(GEOSContextHandle_t handle, GEOSGeometry* g) nogil
+
+    # Geometry creation / destruction
+    GEOSGeometry* GEOSGeom_clone_r(GEOSContextHandle_t handle, const GEOSGeometry* g) nogil
+    GEOSGeometry* GEOSGeom_createPoint_r(GEOSContextHandle_t handle, GEOSCoordSequence* s) nogil
+    GEOSGeometry* GEOSGeom_createLineString_r(GEOSContextHandle_t handle, GEOSCoordSequence* s) nogil
+    GEOSGeometry* GEOSGeom_createLinearRing_r(GEOSContextHandle_t handle, GEOSCoordSequence* s) nogil
+    GEOSGeometry* GEOSGeom_createEmptyPolygon_r(GEOSContextHandle_t handle) nogil
+    GEOSGeometry* GEOSGeom_createPolygon_r(GEOSContextHandle_t handle, GEOSGeometry* shell, GEOSGeometry** holes, unsigned int nholes) nogil
+    GEOSGeometry* GEOSGeom_createCollection_r(GEOSContextHandle_t handle, int type, GEOSGeometry** geoms, unsigned int ngeoms) nogil
+    void GEOSGeom_destroy_r(GEOSContextHandle_t handle, GEOSGeometry* g) nogil
+
+    # Coordinate sequences
+    const GEOSCoordSequence* GEOSGeom_getCoordSeq_r(GEOSContextHandle_t handle, const GEOSGeometry* g)
+    GEOSCoordSequence* GEOSCoordSeq_clone_r(GEOSContextHandle_t handle, const GEOSCoordSequence* s)
+    GEOSCoordSequence* GEOSCoordSeq_create_r(GEOSContextHandle_t handle, unsigned int size, unsigned int dims) nogil
+    void GEOSCoordSeq_destroy_r(GEOSContextHandle_t handle, GEOSCoordSequence* s) nogil
+    int GEOSCoordSeq_setX_r(GEOSContextHandle_t handle, GEOSCoordSequence* s, unsigned int idx, double val) nogil
+    int GEOSCoordSeq_setY_r(GEOSContextHandle_t handle, GEOSCoordSequence* s, unsigned int idx, double val) nogil
+    int GEOSCoordSeq_setZ_r(GEOSContextHandle_t handle, GEOSCoordSequence* s, unsigned int idx, double val) nogil
+    int GEOSCoordSeq_getSize_r(GEOSContextHandle_t handle, GEOSCoordSequence* s, unsigned int* size) nogil
+
+
+cdef class get_geos_handle:
+    cdef GEOSContextHandle_t handle
+    cdef char* last_error
+    cdef char* last_warning
+    cdef GEOSContextHandle_t __enter__(self)

.venv/lib/python3.10/site-packages/shapely/_pygeos_api.pxd

@@ -0,0 +1,58 @@
+"""
+Provides a wrapper for the shapely.lib C API for use in Cython.
+Internally, the shapely C extension uses a PyCapsule to provide run-time access
+to function pointers within the C API.
+
+To use these functions, you must first call the following function in each Cython module:
+`import_shapely_c_api()`
+
+This uses a macro to dynamically load the functions from pointers in the PyCapsule.
+Each C function in shapely.lib exposed in the C API must be specially-wrapped to enable
+this capability.
+
+Segfaults will occur if the C API is not imported properly.
+"""
+
+cimport numpy as np
+from cpython.ref cimport PyObject
+
+from shapely._geos cimport GEOSContextHandle_t, GEOSCoordSequence, GEOSGeometry
+
+
+cdef extern from "c_api.h":
+    cdef enum ShapelyErrorCode:
+        PGERR_SUCCESS,
+        PGERR_NOT_A_GEOMETRY,
+        PGERR_GEOS_EXCEPTION,
+        PGERR_NO_MALLOC,
+        PGERR_GEOMETRY_TYPE,
+        PGERR_MULTIPOINT_WITH_POINT_EMPTY,
+        PGERR_COORD_OUT_OF_BOUNDS,
+        PGERR_EMPTY_GEOMETRY,
+        PGERR_GEOJSON_EMPTY_POINT,
+        PGERR_LINEARRING_NCOORDS,
+        PGERR_NAN_COORD,
+        PGWARN_INVALID_WKB,
+        PGWARN_INVALID_WKT,
+        PGWARN_INVALID_GEOJSON,
+        PGERR_PYSIGNAL
+
+    cpdef enum ShapelyHandleNan:
+        SHAPELY_HANDLE_NAN_ALLOW,
+        SHAPELY_HANDLE_NAN_SKIP,
+        SHAPELY_HANDLE_NANS_ERROR
+
+    # shapely.lib C API loader; returns -1 on error
+    # MUST be called before calling other C API functions
+    int import_shapely_c_api() except -1
+
+    # C functions provided by the shapely.lib C API
+    # Note: GeometryObjects are always managed as Python objects
+    # in Cython to avoid memory leaks, not PyObject* (even though
+    # they are declared that way in the header file).
+    object PyGEOS_CreateGeometry(GEOSGeometry *ptr, GEOSContextHandle_t ctx)
+    char PyGEOS_GetGEOSGeometry(PyObject *obj, GEOSGeometry **out) nogil
+    int PyGEOS_CoordSeq_FromBuffer(
+        GEOSContextHandle_t ctx, const double* buf, unsigned int size,
+        unsigned int dims, char is_ring, int handle_nan,
+        GEOSCoordSequence** coord_seq) nogil

.venv/lib/python3.10/site-packages/shapely/_ragged_array.py

@@ -0,0 +1,475 @@
+"""Provides a conversion to / from a ragged array representation of geometries.
+
+A ragged (or "jagged") array is an irregular array of arrays of which each
+element can have a different length. As a result, such an array cannot be
+represented as a standard, rectangular nD array.
+The coordinates of geometries can be represented as arrays of arrays of
+coordinate pairs (possibly multiple levels of nesting, depending on the
+geometry type).
+
+Geometries, as a ragged array of coordinates, can be efficiently represented
+as contiguous arrays of coordinates provided that there is another data
+structure that keeps track of which range of coordinate values corresponds
+to a given geometry. This can be done using offsets, counts, or indices.
+
+This module currently implements offsets into the coordinates array. This
+is the ragged array representation defined by the the Apache Arrow project
+as "variable size list array" (https://arrow.apache.org/docs/format/Columnar.html#variable-size-list-layout).
+See for example https://cfconventions.org/Data/cf-conventions/cf-conventions-1.9/cf-conventions.html#representations-features
+for different options.
+
+The exact usage of the Arrow list array with varying degrees of nesting for the
+different geometry types is defined by the GeoArrow project:
+https://github.com/geoarrow/geoarrow
+
+"""
+
+import numpy as np
+
+from shapely import creation, geos_version
+from shapely._geometry import (
+    GeometryType,
+    get_parts,
+    get_rings,
+    get_type_id,
+)
+from shapely._geometry_helpers import (
+    _from_ragged_array_multi_linear,
+    _from_ragged_array_multipolygon,
+)
+from shapely.coordinates import get_coordinates
+from shapely.predicates import is_empty, is_missing
+
+__all__ = ["from_ragged_array", "to_ragged_array"]
+
+_geos_ge_312 = geos_version >= (3, 12, 0)
+
+
+# # GEOS -> coords/offset arrays (to_ragged_array)
+
+
+def _get_arrays_point(arr, include_z, include_m):
+    # only one array of coordinates
+    coords = get_coordinates(arr, include_z=include_z, include_m=include_m)
+
+    # empty points are represented by NaNs
+    # + missing geometries should also be present with some value
+    empties = is_empty(arr) | is_missing(arr)
+    if empties.any():
+        indices = np.nonzero(empties)[0]
+        indices = indices - np.arange(len(indices))
+        coords = np.insert(coords, indices, np.nan, axis=0)
+
+    return coords, ()
+
+
+def _indices_to_offsets(indices, n):
+    # default to int32 offsets if possible (to prefer the non-large arrow list variants)
+    # n_coords is the length of the array the indices are poin
+    if len(indices) > 2147483647:
+        dtype = np.int64
+    else:
+        dtype = np.int32
+
+    offsets = np.insert(np.bincount(indices).cumsum(dtype=dtype), 0, 0)
+
+    if len(offsets) != n + 1:
+        # last geometries might be empty or missing
+        offsets = np.pad(
+            offsets,
+            (0, n + 1 - len(offsets)),
+            "constant",
+            constant_values=offsets[-1],
+        )
+    return offsets
+
+
+def _get_arrays_multipoint(arr, include_z, include_m):
+    # explode/flatten the MultiPoints
+    _, part_indices = get_parts(arr, return_index=True)
+    # the offsets into the multipoint parts
+    offsets = _indices_to_offsets(part_indices, len(arr))
+
+    # only one array of coordinates
+    coords = get_coordinates(arr, include_z=include_z, include_m=include_m)
+
+    return coords, (offsets,)
+
+
+def _get_arrays_linestring(arr, include_z, include_m):
+    # the coords and offsets into the coordinates of the linestrings
+    coords, indices = get_coordinates(
+        arr, return_index=True, include_z=include_z, include_m=include_m
+    )
+    offsets = _indices_to_offsets(indices, len(arr))
+
+    return coords, (offsets,)
+
+
+def _get_arrays_multilinestring(arr, include_z, include_m):
+    # explode/flatten the MultiLineStrings
+    arr_flat, part_indices = get_parts(arr, return_index=True)
+    # the offsets into the multilinestring parts
+    offsets2 = _indices_to_offsets(part_indices, len(arr))
+
+    # the coords and offsets into the coordinates of the linestrings
+    coords, indices = get_coordinates(
+        arr_flat, return_index=True, include_z=include_z, include_m=include_m
+    )
+    offsets1 = _indices_to_offsets(indices, len(arr_flat))
+
+    return coords, (offsets1, offsets2)
+
+
+def _get_arrays_polygon(arr, include_z, include_m):
+    # explode/flatten the Polygons into Rings
+    arr_flat, ring_indices = get_rings(arr, return_index=True)
+    # the offsets into the exterior/interior rings of the multipolygon parts
+    offsets2 = _indices_to_offsets(ring_indices, len(arr))
+
+    # the coords and offsets into the coordinates of the rings
+    coords, indices = get_coordinates(
+        arr_flat, return_index=True, include_z=include_z, include_m=include_m
+    )
+    offsets1 = _indices_to_offsets(indices, len(arr_flat))
+
+    return coords, (offsets1, offsets2)
+
+
+def _get_arrays_multipolygon(arr, include_z, include_m):
+    # explode/flatten the MultiPolygons
+    arr_flat, part_indices = get_parts(arr, return_index=True)
+    # the offsets into the multipolygon parts
+    offsets3 = _indices_to_offsets(part_indices, len(arr))
+
+    # explode/flatten the Polygons into Rings
+    arr_flat2, ring_indices = get_rings(arr_flat, return_index=True)
+    # the offsets into the exterior/interior rings of the multipolygon parts
+    offsets2 = _indices_to_offsets(ring_indices, len(arr_flat))
+
+    # the coords and offsets into the coordinates of the rings
+    coords, indices = get_coordinates(
+        arr_flat2, return_index=True, include_z=include_z, include_m=include_m
+    )
+    offsets1 = _indices_to_offsets(indices, len(arr_flat2))
+
+    return coords, (offsets1, offsets2, offsets3)
+
+
+def to_ragged_array(geometries, include_z=None, include_m=None):
+    """Convert geometries to a ragged array representation.
+
+    This function converts an array of geometries to a ragged array
+    (i.e. irregular array of arrays) of coordinates, represented in memory
+    using a single contiguous array of the coordinates, and
+    up to 3 offset arrays that keep track where each sub-array
+    starts and ends.
+
+    This follows the in-memory layout of the variable size list arrays defined
+    by Apache Arrow, as specified for geometries by the GeoArrow project:
+    https://github.com/geoarrow/geoarrow.
+
+    Parameters
+    ----------
+    geometries : array_like
+        Array of geometries (1-dimensional).
+    include_z, include_m : bool, default None
+        If both are False, return XY (2D) geometries.
+        If both are True, return XYZM (4D) geometries.
+        If either is True, return either XYZ or XYM (3D) geometries.
+        If a geometry has no Z or M dimension, extra coordinate data will be NaN.
+        By default, will infer the dimensionality from the
+        input geometries. Note that this inference can be unreliable with
+        empty geometries (for a guaranteed result, it is recommended to
+        specify the keyword).
+
+        .. versionadded:: 2.1.0
+            The ``include_m`` parameter was added to support XYM (3D) and
+            XYZM (4D) geometries available with GEOS 3.12.0 or later.
+            With older GEOS versions, M dimension coordinates will be NaN.
+
+    Returns
+    -------
+    tuple of (geometry_type, coords, offsets)
+        geometry_type : GeometryType
+            The type of the input geometries (required information for
+            roundtrip).
+        coords : np.ndarray
+            Contiguous array of shape (n, 2), (n, 3), or (n, 4) of all
+            coordinates of all input geometries.
+        offsets: tuple of np.ndarray
+            Offset arrays that make it possible to reconstruct the
+            geometries from the flat coordinates array. The number of
+            offset arrays depends on the geometry type. See
+            https://github.com/geoarrow/geoarrow/blob/main/format.md
+            for details.
+            Uses int32 dtype offsets if possible, otherwise int64 for
+            large inputs (coordinates > 32GB).
+
+    Notes
+    -----
+    Mixed singular and multi geometry types of the same basic type are
+    allowed (e.g., Point and MultiPoint) and all singular types will be
+    treated as multi types.
+    GeometryCollections and other mixed geometry types are not supported.
+
+    See Also
+    --------
+    from_ragged_array
+
+    Examples
+    --------
+    Consider a Polygon with one hole (interior ring):
+
+    >>> import shapely
+    >>> from shapely import Polygon
+    >>> polygon = Polygon(
+    ...     [(0, 0), (10, 0), (10, 10), (0, 10)],
+    ...     holes=[[(2, 2), (3, 2), (2, 3)]]
+    ... )
+    >>> polygon
+    <POLYGON ((0 0, 10 0, 10 10, 0 10, 0 0), (2 2, 3 2, 2 3, 2 2))>
+
+    This polygon can be thought of as a list of rings (first ring is the
+    exterior ring, subsequent rings are the interior rings), and each ring
+    as a list of coordinate pairs. This is very similar to how GeoJSON
+    represents the coordinates:
+
+    >>> import json
+    >>> json.loads(shapely.to_geojson(polygon))["coordinates"]
+    [[[0.0, 0.0], [10.0, 0.0], [10.0, 10.0], [0.0, 10.0], [0.0, 0.0]],
+     [[2.0, 2.0], [3.0, 2.0], [2.0, 3.0], [2.0, 2.0]]]
+
+    This function will return a similar list of lists of lists, but
+    using a single contiguous array of coordinates, and multiple arrays of
+    offsets:
+
+    >>> geometry_type, coords, offsets = shapely.to_ragged_array([polygon])
+    >>> geometry_type
+    <GeometryType.POLYGON: 3>
+    >>> coords
+    array([[ 0.,  0.],
+           [10.,  0.],
+           [10., 10.],
+           [ 0., 10.],
+           [ 0.,  0.],
+           [ 2.,  2.],
+           [ 3.,  2.],
+           [ 2.,  3.],
+           [ 2.,  2.]])
+
+    >>> offsets
+    (array([0, 5, 9], dtype=int32), array([0, 2], dtype=int32))
+
+    As an example how to interpret the offsets: the i-th ring in the
+    coordinates is represented by ``offsets[0][i]`` to ``offsets[0][i+1]``:
+
+    >>> exterior_ring_start, exterior_ring_end = offsets[0][0], offsets[0][1]
+    >>> coords[exterior_ring_start:exterior_ring_end]
+    array([[ 0.,  0.],
+           [10.,  0.],
+           [10., 10.],
+           [ 0., 10.],
+           [ 0.,  0.]])
+
+    """
+    from shapely import has_m, has_z  # avoid circular import
+
+    geometries = np.asarray(geometries)
+    if include_z is None:
+        include_z = np.any(has_z(geometries[~is_empty(geometries)]))
+    if include_m is None:
+        if _geos_ge_312:
+            include_m = np.any(has_m(geometries[~is_empty(geometries)]))
+        else:
+            include_m = False
+
+    geom_types = np.unique(get_type_id(geometries))
+    # ignore missing values (type of -1)
+    geom_types = geom_types[geom_types >= 0]
+
+    get_arrays_args = geometries, include_z, include_m
+    if len(geom_types) == 1:
+        typ = GeometryType(geom_types[0])
+        if typ == GeometryType.POINT:
+            coords, offsets = _get_arrays_point(*get_arrays_args)
+        elif typ == GeometryType.LINESTRING:
+            coords, offsets = _get_arrays_linestring(*get_arrays_args)
+        elif typ == GeometryType.POLYGON:
+            coords, offsets = _get_arrays_polygon(*get_arrays_args)
+        elif typ == GeometryType.MULTIPOINT:
+            coords, offsets = _get_arrays_multipoint(*get_arrays_args)
+        elif typ == GeometryType.MULTILINESTRING:
+            coords, offsets = _get_arrays_multilinestring(*get_arrays_args)
+        elif typ == GeometryType.MULTIPOLYGON:
+            coords, offsets = _get_arrays_multipolygon(*get_arrays_args)
+        else:
+            raise ValueError(f"Geometry type {typ.name} is not supported")
+
+    elif len(geom_types) == 2:
+        if set(geom_types) == {GeometryType.POINT, GeometryType.MULTIPOINT}:
+            typ = GeometryType.MULTIPOINT
+            coords, offsets = _get_arrays_multipoint(*get_arrays_args)
+        elif set(geom_types) == {GeometryType.LINESTRING, GeometryType.MULTILINESTRING}:
+            typ = GeometryType.MULTILINESTRING
+            coords, offsets = _get_arrays_multilinestring(*get_arrays_args)
+        elif set(geom_types) == {GeometryType.POLYGON, GeometryType.MULTIPOLYGON}:
+            typ = GeometryType.MULTIPOLYGON
+            coords, offsets = _get_arrays_multipolygon(*get_arrays_args)
+        else:
+            raise ValueError(
+                "Geometry type combination is not supported "
+                f"({[GeometryType(t).name for t in geom_types]})"
+            )
+    else:
+        raise ValueError(
+            "Geometry type combination is not supported "
+            f"({[GeometryType(t).name for t in geom_types]})"
+        )
+
+    return typ, coords, offsets
+
+
+# # coords/offset arrays -> GEOS (from_ragged_array)
+
+
+def _point_from_flatcoords(coords):
+    result = creation.points(coords)
+
+    # Older versions of GEOS (<= 3.9) don't automatically convert NaNs
+    # to empty points -> do manually
+    empties = np.isnan(coords).all(axis=1)
+    if empties.any():
+        result[empties] = creation.empty(1, geom_type=GeometryType.POINT).item()
+
+    return result
+
+
+def _multipoint_from_flatcoords(coords, offsets):
+    # recreate points
+    if len(offsets):
+        coords = coords[offsets[0] :]
+    points = creation.points(coords)
+
+    # recreate multipoints
+    multipoint_parts = np.diff(offsets)
+    multipoint_indices = np.repeat(np.arange(len(multipoint_parts)), multipoint_parts)
+
+    result = np.empty(len(offsets) - 1, dtype=object)
+    result = creation.multipoints(points, indices=multipoint_indices, out=result)
+    result[multipoint_parts == 0] = creation.empty(
+        1, geom_type=GeometryType.MULTIPOINT
+    ).item()
+
+    return result
+
+
+def _linestring_from_flatcoords(coords, offsets):
+    # recreate linestrings
+    if len(offsets):
+        coords = coords[offsets[0] :]
+    linestring_n = np.diff(offsets)
+    linestring_indices = np.repeat(np.arange(len(linestring_n)), linestring_n)
+
+    result = np.empty(len(offsets) - 1, dtype=object)
+    result = creation.linestrings(coords, indices=linestring_indices, out=result)
+    result[linestring_n == 0] = creation.empty(
+        1, geom_type=GeometryType.LINESTRING
+    ).item()
+    return result
+
+
+def _multilinestrings_from_flatcoords(coords, offsets1, offsets2):
+    # ensure correct dtypes
+    offsets1 = np.asarray(offsets1, dtype="int64")
+    offsets2 = np.asarray(offsets2, dtype="int64")
+
+    # recreate multilinestrings
+    result = _from_ragged_array_multi_linear(
+        coords, offsets1, offsets2, geometry_type=GeometryType.MULTILINESTRING
+    )
+    return result
+
+
+def _polygon_from_flatcoords(coords, offsets1, offsets2):
+    # ensure correct dtypes
+    offsets1 = np.asarray(offsets1, dtype="int64")
+    offsets2 = np.asarray(offsets2, dtype="int64")
+
+    # recreate polygons
+    result = _from_ragged_array_multi_linear(
+        coords, offsets1, offsets2, geometry_type=GeometryType.POLYGON
+    )
+    return result
+
+
+def _multipolygons_from_flatcoords(coords, offsets1, offsets2, offsets3):
+    # ensure correct dtypes
+    offsets1 = np.asarray(offsets1, dtype="int64")
+    offsets2 = np.asarray(offsets2, dtype="int64")
+    offsets3 = np.asarray(offsets3, dtype="int64")
+
+    # recreate multipolygons
+    result = _from_ragged_array_multipolygon(coords, offsets1, offsets2, offsets3)
+    return result
+
+
+def from_ragged_array(geometry_type, coords, offsets=None):
+    """Create geometries from a contiguous array of coordinates and offset arrays.
+
+    This function creates geometries from the ragged array representation
+    as returned by ``to_ragged_array``.
+
+    This follows the in-memory layout of the variable size list arrays defined
+    by Apache Arrow, as specified for geometries by the GeoArrow project:
+    https://github.com/geoarrow/geoarrow.
+
+    See :func:`to_ragged_array` for more details.
+
+    Parameters
+    ----------
+    geometry_type : GeometryType
+        The type of geometry to create.
+    coords : np.ndarray
+        Contiguous array of shape (n, 2) or (n, 3) of all coordinates
+        for the geometries.
+    offsets: tuple of np.ndarray
+        Offset arrays that allow to reconstruct the geometries based on the
+        flat coordinates array. The number of offset arrays depends on the
+        geometry type. See
+        https://github.com/geoarrow/geoarrow/blob/main/format.md for details.
+
+    Returns
+    -------
+    np.ndarray
+        Array of geometries (1-dimensional).
+
+    See Also
+    --------
+    to_ragged_array
+
+    """
+    coords = np.asarray(coords, dtype="float64")
+
+    if geometry_type == GeometryType.POINT:
+        if not (offsets is None or len(offsets) == 0):
+            raise ValueError("'offsets' should not be provided for geometry type Point")
+        return _point_from_flatcoords(coords)
+
+    if offsets is None:
+        raise ValueError(
+            "'offsets' must be provided for any geometry type except for Point"
+        )
+
+    if geometry_type == GeometryType.LINESTRING:
+        return _linestring_from_flatcoords(coords, *offsets)
+    elif geometry_type == GeometryType.POLYGON:
+        return _polygon_from_flatcoords(coords, *offsets)
+    elif geometry_type == GeometryType.MULTIPOINT:
+        return _multipoint_from_flatcoords(coords, *offsets)
+    elif geometry_type == GeometryType.MULTILINESTRING:
+        return _multilinestrings_from_flatcoords(coords, *offsets)
+    elif geometry_type == GeometryType.MULTIPOLYGON:
+        return _multipolygons_from_flatcoords(coords, *offsets)
+    else:
+        raise ValueError(f"Geometry type {geometry_type.name} is not supported")

.venv/lib/python3.10/site-packages/shapely/_version.py

@@ -0,0 +1,21 @@
+
+# This file was generated by 'versioneer.py' (0.28) from
+# revision-control system data, or from the parent directory name of an
+# unpacked source archive. Distribution tarballs contain a pre-generated copy
+# of this file.
+
+import json
+
+version_json = '''
+{
+ "date": "2025-04-03T10:55:05+0200",
+ "dirty": false,
+ "error": null,
+ "full-revisionid": "4940c6405ac9ef2d77c9e9990954b68294c3c399",
+ "version": "2.1.0"
+}
+'''  # END VERSION_JSON
+
+
+def get_versions():
+    return json.loads(version_json)

.venv/lib/python3.10/site-packages/shapely/affinity.py

@@ -0,0 +1,266 @@
+"""Affine transforms, both in general and specific, named transforms."""
+
+from math import cos, pi, sin, tan
+
+import numpy as np
+
+import shapely
+
+__all__ = ["affine_transform", "rotate", "scale", "skew", "translate"]
+
+
+def affine_transform(geom, matrix):
+    r"""Return a transformed geometry using an affine transformation matrix.
+
+    The coefficient matrix is provided as a list or tuple with 6 or 12 items
+    for 2D or 3D transformations, respectively.
+
+    For 2D affine transformations, the 6 parameter matrix is::
+
+        [a, b, d, e, xoff, yoff]
+
+    which represents the augmented matrix::
+
+        [x']   / a  b xoff \ [x]
+        [y'] = | d  e yoff | [y]
+        [1 ]   \ 0  0   1  / [1]
+
+    or the equations for the transformed coordinates::
+
+        x' = a * x + b * y + xoff
+        y' = d * x + e * y + yoff
+
+    For 3D affine transformations, the 12 parameter matrix is::
+
+        [a, b, c, d, e, f, g, h, i, xoff, yoff, zoff]
+
+    which represents the augmented matrix::
+
+        [x']   / a  b  c xoff \ [x]
+        [y'] = | d  e  f yoff | [y]
+        [z']   | g  h  i zoff | [z]
+        [1 ]   \ 0  0  0   1  / [1]
+
+    or the equations for the transformed coordinates::
+
+        x' = a * x + b * y + c * z + xoff
+        y' = d * x + e * y + f * z + yoff
+        z' = g * x + h * y + i * z + zoff
+    """
+    if len(matrix) == 6:
+        ndim = 2
+        a, b, d, e, xoff, yoff = matrix
+        if geom.has_z:
+            ndim = 3
+            i = 1.0
+            c = f = g = h = zoff = 0.0
+    elif len(matrix) == 12:
+        ndim = 3
+        a, b, c, d, e, f, g, h, i, xoff, yoff, zoff = matrix
+        if not geom.has_z:
+            ndim = 2
+    else:
+        raise ValueError("'matrix' expects either 6 or 12 coefficients")
+
+    # if ndim == 2:
+    #     A = np.array([[a, b], [d, e]], dtype=float)
+    #     off = np.array([xoff, yoff], dtype=float)
+    # else:
+    #     A = np.array([[a, b, c], [d, e, f], [g, h, i]], dtype=float)
+    #     off = np.array([xoff, yoff, zoff], dtype=float)
+
+    def _affine_coords(coords):
+        # These are equivalent, but unfortunately not robust
+        #   result = np.matmul(coords, A.T) + off
+        #   result = np.matmul(A, coords.T).T + off
+        # Therefore, manual matrix multiplication is needed
+        if ndim == 2:
+            x, y = coords.T
+            xp = a * x + b * y + xoff
+            yp = d * x + e * y + yoff
+            result = np.stack([xp, yp]).T
+        elif ndim == 3:
+            x, y, z = coords.T
+            xp = a * x + b * y + c * z + xoff
+            yp = d * x + e * y + f * z + yoff
+            zp = g * x + h * y + i * z + zoff
+            result = np.stack([xp, yp, zp]).T
+        return result
+
+    return shapely.transform(geom, _affine_coords, include_z=ndim == 3)
+
+
+def interpret_origin(geom, origin, ndim):
+    """Return interpreted coordinate tuple for origin parameter.
+
+    This is a helper function for other transform functions.
+
+    The point of origin can be a keyword 'center' for the 2D bounding box
+    center, 'centroid' for the geometry's 2D centroid, a Point object or a
+    coordinate tuple (x0, y0, z0).
+    """
+    # get coordinate tuple from 'origin' from keyword or Point type
+    if origin == "center":
+        # bounding box center
+        minx, miny, maxx, maxy = geom.bounds
+        origin = ((maxx + minx) / 2.0, (maxy + miny) / 2.0)
+    elif origin == "centroid":
+        origin = geom.centroid.coords[0]
+    elif isinstance(origin, str):
+        raise ValueError(f"'origin' keyword {origin!r} is not recognized")
+    elif getattr(origin, "geom_type", None) == "Point":
+        origin = origin.coords[0]
+
+    # origin should now be tuple-like
+    if len(origin) not in (2, 3):
+        raise ValueError("Expected number of items in 'origin' to be either 2 or 3")
+    if ndim == 2:
+        return origin[0:2]
+    else:  # 3D coordinate
+        if len(origin) == 2:
+            return origin + (0.0,)
+        else:
+            return origin
+
+
+def rotate(geom, angle, origin="center", use_radians=False):
+    r"""Return a rotated geometry on a 2D plane.
+
+    The angle of rotation can be specified in either degrees (default) or
+    radians by setting ``use_radians=True``. Positive angles are
+    counter-clockwise and negative are clockwise rotations.
+
+    The point of origin can be a keyword 'center' for the bounding box
+    center (default), 'centroid' for the geometry's centroid, a Point object
+    or a coordinate tuple (x0, y0).
+
+    The affine transformation matrix for 2D rotation is:
+
+      / cos(r) -sin(r) xoff \
+      | sin(r)  cos(r) yoff |
+      \   0       0      1  /
+
+    where the offsets are calculated from the origin Point(x0, y0):
+
+        xoff = x0 - x0 * cos(r) + y0 * sin(r)
+        yoff = y0 - x0 * sin(r) - y0 * cos(r)
+    """
+    if geom.is_empty:
+        return geom
+    if not use_radians:  # convert from degrees
+        angle = angle * pi / 180.0
+    cosp = cos(angle)
+    sinp = sin(angle)
+    if abs(cosp) < 2.5e-16:
+        cosp = 0.0
+    if abs(sinp) < 2.5e-16:
+        sinp = 0.0
+    x0, y0 = interpret_origin(geom, origin, 2)
+
+    # fmt: off
+    matrix = (cosp, -sinp, 0.0,
+              sinp, cosp, 0.0,
+              0.0, 0.0, 1.0,
+              x0 - x0 * cosp + y0 * sinp, y0 - x0 * sinp - y0 * cosp, 0.0)
+    # fmt: on
+    return affine_transform(geom, matrix)
+
+
+def scale(geom, xfact=1.0, yfact=1.0, zfact=1.0, origin="center"):
+    r"""Return a scaled geometry, scaled by factors along each dimension.
+
+    The point of origin can be a keyword 'center' for the 2D bounding box
+    center (default), 'centroid' for the geometry's 2D centroid, a Point
+    object or a coordinate tuple (x0, y0, z0).
+
+    Negative scale factors will mirror or reflect coordinates.
+
+    The general 3D affine transformation matrix for scaling is:
+
+        / xfact  0    0   xoff \
+        |   0  yfact  0   yoff |
+        |   0    0  zfact zoff |
+        \   0    0    0     1  /
+
+    where the offsets are calculated from the origin Point(x0, y0, z0):
+
+        xoff = x0 - x0 * xfact
+        yoff = y0 - y0 * yfact
+        zoff = z0 - z0 * zfact
+    """
+    if geom.is_empty:
+        return geom
+    x0, y0, z0 = interpret_origin(geom, origin, 3)
+
+    # fmt: off
+    matrix = (xfact, 0.0, 0.0,
+              0.0, yfact, 0.0,
+              0.0, 0.0, zfact,
+              x0 - x0 * xfact, y0 - y0 * yfact, z0 - z0 * zfact)
+    # fmt: on
+    return affine_transform(geom, matrix)
+
+
+def skew(geom, xs=0.0, ys=0.0, origin="center", use_radians=False):
+    r"""Return a skewed geometry, sheared by angles along x and y dimensions.
+
+    The shear angle can be specified in either degrees (default) or radians
+    by setting ``use_radians=True``.
+
+    The point of origin can be a keyword 'center' for the bounding box
+    center (default), 'centroid' for the geometry's centroid, a Point object
+    or a coordinate tuple (x0, y0).
+
+    The general 2D affine transformation matrix for skewing is:
+
+        /   1    tan(xs) xoff \
+        | tan(ys)  1     yoff |
+        \   0      0       1  /
+
+    where the offsets are calculated from the origin Point(x0, y0):
+
+        xoff = -y0 * tan(xs)
+        yoff = -x0 * tan(ys)
+    """
+    if geom.is_empty:
+        return geom
+    if not use_radians:  # convert from degrees
+        xs = xs * pi / 180.0
+        ys = ys * pi / 180.0
+    tanx = tan(xs)
+    tany = tan(ys)
+    if abs(tanx) < 2.5e-16:
+        tanx = 0.0
+    if abs(tany) < 2.5e-16:
+        tany = 0.0
+    x0, y0 = interpret_origin(geom, origin, 2)
+
+    # fmt: off
+    matrix = (1.0, tanx, 0.0,
+              tany, 1.0, 0.0,
+              0.0, 0.0, 1.0,
+              -y0 * tanx, -x0 * tany, 0.0)
+    # fmt: on
+    return affine_transform(geom, matrix)
+
+
+def translate(geom, xoff=0.0, yoff=0.0, zoff=0.0):
+    r"""Return a translated geometry shifted by offsets along each dimension.
+
+    The general 3D affine transformation matrix for translation is:
+
+        / 1  0  0 xoff \
+        | 0  1  0 yoff |
+        | 0  0  1 zoff |
+        \ 0  0  0   1  /
+    """
+    if geom.is_empty:
+        return geom
+
+    # fmt: off
+    matrix = (1.0, 0.0, 0.0,
+              0.0, 1.0, 0.0,
+              0.0, 0.0, 1.0,
+              xoff, yoff, zoff)
+    # fmt: on
+    return affine_transform(geom, matrix)

.venv/lib/python3.10/site-packages/shapely/algorithms/__init__.py

@@ -0,0 +1 @@
+"""Algorithms implemented in Shapely."""

.venv/lib/python3.10/site-packages/shapely/algorithms/_oriented_envelope.py

@@ -0,0 +1,56 @@
+import math
+from itertools import islice
+
+import numpy as np
+
+import shapely
+from shapely.affinity import affine_transform
+
+
+def _oriented_envelope_min_area(geometry, **kwargs):
+    """Compute the oriented envelope (minimum rotated rectangle).
+
+    This is a fallback implementation for GEOS < 3.12 to have the correct
+    minimum area behaviour.
+    """
+    if geometry is None:
+        return None
+    if geometry.is_empty:
+        return shapely.from_wkt("POLYGON EMPTY")
+
+    # first compute the convex hull
+    hull = geometry.convex_hull
+    try:
+        coords = hull.exterior.coords
+    except AttributeError:  # may be a Point or a LineString
+        return hull
+    # generate the edge vectors between the convex hull's coords
+    edges = (
+        (pt2[0] - pt1[0], pt2[1] - pt1[1])
+        for pt1, pt2 in zip(coords, islice(coords, 1, None))
+    )
+
+    def _transformed_rects():
+        for dx, dy in edges:
+            # compute the normalized direction vector of the edge
+            # vector.
+            length = math.sqrt(dx**2 + dy**2)
+            ux, uy = dx / length, dy / length
+            # compute the normalized perpendicular vector
+            vx, vy = -uy, ux
+            # transform hull from the original coordinate system to
+            # the coordinate system defined by the edge and compute
+            # the axes-parallel bounding rectangle.
+            transf_rect = affine_transform(hull, (ux, uy, vx, vy, 0, 0)).envelope
+            # yield the transformed rectangle and a matrix to
+            # transform it back to the original coordinate system.
+            yield (transf_rect, (ux, vx, uy, vy, 0, 0))
+
+    # check for the minimum area rectangle and return it
+    transf_rect, inv_matrix = min(_transformed_rects(), key=lambda r: r[0].area)
+    return affine_transform(transf_rect, inv_matrix)
+
+
+_oriented_envelope_min_area_vectorized = np.frompyfunc(
+    _oriented_envelope_min_area, 1, 1
+)

.venv/lib/python3.10/site-packages/shapely/algorithms/cga.py

@@ -0,0 +1,50 @@
+"""Shapely CGA algorithms."""
+
+import numpy as np
+
+import shapely
+
+
+def signed_area(ring):
+    """Return the signed area enclosed by a ring in linear time.
+
+    Algorithm used: https://web.archive.org/web/20080209143651/http://cgafaq.info:80/wiki/Polygon_Area
+    """
+    coords = np.array(ring.coords)[:, :2]
+    xs, ys = np.vstack([coords, coords[1]]).T
+    return np.sum(xs[1:-1] * (ys[2:] - ys[:-2])) / 2.0
+
+
+def _reverse_conditioned(rings, condition):
+    """Return a copy of the rings potentially reversed depending on `condition`."""
+    condition = np.asarray(condition)
+    if np.all(condition):
+        rings = shapely.reverse(rings)
+    elif np.any(condition):
+        rings = np.array(rings)
+        rings[condition] = shapely.reverse(rings[condition])
+    return rings
+
+
+def _orient_polygon(geometry, exterior_cw=False):
+    if geometry is None:
+        return None
+    if geometry.geom_type in ["MultiPolygon", "GeometryCollection"]:
+        return geometry.__class__(
+            [_orient_polygon(geom, exterior_cw) for geom in geometry.geoms]
+        )
+    # elif geometry.geom_type in ["LinearRing"]:
+    #     return reverse_conditioned(geometry, is_ccw(geometry) != ccw)
+    elif geometry.geom_type == "Polygon":
+        rings = np.array([geometry.exterior, *geometry.interiors])
+        reverse_condition = shapely.is_ccw(rings)
+        reverse_condition[0] = not reverse_condition[0]
+        if exterior_cw:
+            reverse_condition = np.logical_not(reverse_condition)
+        if np.any(reverse_condition):
+            rings = _reverse_conditioned(rings, reverse_condition)
+            return geometry.__class__(rings[0], rings[1:])
+    return geometry
+
+
+_orient_polygons_vectorized = np.frompyfunc(_orient_polygon, nin=2, nout=1)

.venv/lib/python3.10/site-packages/shapely/algorithms/polylabel.py

@@ -0,0 +1,43 @@
+"""Provides functions for finding the pole of inaccessibility for a given polygon."""
+
+from shapely._geometry import get_point
+from shapely.constructive import maximum_inscribed_circle
+
+
+def polylabel(polygon, tolerance=1.0):
+    """Find pole of inaccessibility for a given polygon.
+
+    Based on Vladimir Agafonkin's https://github.com/mapbox/polylabel
+
+    Parameters
+    ----------
+    polygon : shapely.geometry.Polygon
+        Polygon for which to find the pole of inaccessibility.
+    tolerance : int or float, optional
+        `tolerance` represents the highest resolution in units of the
+        input geometry that will be considered for a solution. (default
+        value is 1.0).
+
+    Returns
+    -------
+    shapely.geometry.Point
+        A point representing the pole of inaccessibility for the given input
+        polygon.
+
+    Raises
+    ------
+    shapely.errors.TopologicalError
+        If the input polygon is not a valid geometry.
+
+    Examples
+    --------
+    >>> from shapely.ops import polylabel
+    >>> from shapely import LineString
+    >>> polygon = LineString([(0, 0), (50, 200), (100, 100), (20, 50),
+    ... (-100, -20), (-150, -200)]).buffer(100)
+    >>> polylabel(polygon, tolerance=0.001)
+    <POINT (59.733 111.33)>
+
+    """
+    line = maximum_inscribed_circle(polygon, tolerance)
+    return get_point(line, 0)

.venv/lib/python3.10/site-packages/shapely/conftest.py

@@ -0,0 +1,70 @@
+"""Pytest and scipy-doctest configuration for Shapely."""
+
+import numpy
+import pytest
+
+from shapely import geos_version_string
+
+try:
+    from scipy_doctest.conftest import dt_config
+
+    HAVE_SCPDT = True
+except ModuleNotFoundError:
+    HAVE_SCPDT = False
+
+shapely20_todo = pytest.mark.xfail(
+    strict=True, reason="Not yet implemented for Shapely 2.0"
+)
+shapely20_wontfix = pytest.mark.xfail(strict=True, reason="Will fail for Shapely 2.0")
+
+
+def pytest_report_header(config):
+    """Header for pytest."""
+    vers = [
+        f"GEOS version: {geos_version_string}",
+        f"NumPy version: {numpy.__version__}",
+    ]
+    return "\n".join(vers)
+
+
+if HAVE_SCPDT:
+    import doctest
+    import warnings
+    from contextlib import contextmanager
+
+    @contextmanager
+    def warnings_errors_and_rng(test=None):
+        """Filter out some warnings."""
+        depr_msgs = "|".join(
+            [
+                # https://github.com/pyproj4/pyproj/issues/1468
+                "Conversion of an array with ndim",
+            ]
+        )
+
+        runtime_msgs = "|".join(
+            [
+                # https://github.com/libgeos/geos/pull/1226
+                "invalid value encountered in coverage_union",
+            ]
+        )
+
+        with warnings.catch_warnings():
+            if depr_msgs:
+                warnings.filterwarnings("ignore", depr_msgs, DeprecationWarning)
+            if runtime_msgs:
+                warnings.filterwarnings("ignore", runtime_msgs, RuntimeWarning)
+            yield
+
+    # find and check doctests under this context manager
+    dt_config.user_context_mgr = warnings_errors_and_rng
+
+    # relax all NumPy scalar type repr, e.g. `np.int32(0)` matches `0`
+    dt_config.strict_check = False
+
+    dt_config.optionflags = doctest.NORMALIZE_WHITESPACE | doctest.ELLIPSIS
+
+    # ignores are for things fail doctest collection (optionals etc)
+    dt_config.pytest_extra_ignore = [
+        "shapely/geos.py",
+    ]

.venv/lib/python3.10/site-packages/shapely/coordinates.py

@@ -0,0 +1,335 @@
+"""Methods that operate on the coordinates of geometries."""
+
+import numpy as np
+
+import shapely
+from shapely import lib
+from shapely.decorators import deprecate_positional
+
+__all__ = ["count_coordinates", "get_coordinates", "set_coordinates", "transform"]
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0: only supported XY and XYZ geometries
+#   transform(geometry, transformation, include_z=False)
+# shapely 2.1: shows deprecation warning about positional 'include_z' arg
+#   transform(geometry, transformation, include_z=False, *, interleaved=True)
+# shapely 2.2(?): enforce keyword-only arguments after 'transformation'
+#   transform(geometry, transformation, *, include_z=False, interleaved=True)
+
+
+@deprecate_positional(["include_z"], category=DeprecationWarning)
+def transform(
+    geometry,
+    transformation,
+    include_z: bool | None = False,
+    *,
+    interleaved: bool = True,
+):
+    """Apply a function to the coordinates of a geometry.
+
+    With the default of ``include_z=False``, all returned geometries will be
+    two-dimensional; the third dimension will be discarded, if present.
+    When specifying ``include_z=True``, the returned geometries preserve
+    the dimensionality of the respective input geometries.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries to transform.
+    transformation : function
+        A function that transforms a (N, 2) or (N, 3) ndarray of float64 to
+        another (N, 2) or (N, 3) ndarray of float64.
+        The function may not change N.
+    include_z : bool, optional, default False
+        If False, always return 2D geometries.
+        If True, the data being passed to the
+        transformation function will include the third dimension
+        (if a geometry has no third dimension, the z-coordinates
+        will be NaN). If None, will infer the dimensionality per
+        input geometry using ``has_z``, which may result in 2 calls to
+        the transformation function. Note that this inference
+        can be unreliable with empty geometries or NaN coordinates: for a
+        guaranteed result, it is recommended to specify ``include_z`` explicitly.
+    interleaved : bool, default True
+        If set to False, the transformation function should accept 2 or 3 separate
+        one-dimensional arrays (x, y and optional z) instead of a single
+        two-dimensional array.
+
+        .. versionadded:: 2.1.0
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``include_z`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    has_z
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point
+    >>> shapely.transform(Point(0, 0), lambda x: x + 1)
+    <POINT (1 1)>
+    >>> shapely.transform(LineString([(2, 2), (4, 4)]), lambda x: x * [2, 3])
+    <LINESTRING (4 6, 8 12)>
+    >>> shapely.transform(None, lambda x: x) is None
+    True
+    >>> shapely.transform([Point(0, 0), None], lambda x: x).tolist()
+    [<POINT (0 0)>, None]
+
+    The presence of a third dimension can be automatically detected, or
+    controlled explicitly:
+
+    >>> shapely.transform(Point(0, 0, 0), lambda x: x + 1)
+    <POINT (1 1)>
+    >>> shapely.transform(Point(0, 0, 0), lambda x: x + 1, include_z=True)
+    <POINT Z (1 1 1)>
+    >>> shapely.transform(Point(0, 0, 0), lambda x: x + 1, include_z=None)
+    <POINT Z (1 1 1)>
+
+    With interleaved=False, the call signature of the transformation is different:
+
+    >>> shapely.transform(LineString([(1, 2), (3, 4)]), lambda x, y: (x + 1, y), \
+interleaved=False)
+    <LINESTRING (2 2, 4 4)>
+
+    Or with a z coordinate:
+
+    >>> shapely.transform(Point(0, 0, 0), lambda x, y, z: (x + 1, y, z + 2), \
+interleaved=False, include_z=True)
+    <POINT Z (1 0 2)>
+
+    Using pyproj >= 2.1, the following example will reproject Shapely geometries
+    from EPSG 4326 to EPSG 32618:
+
+    >>> from pyproj import Transformer
+    >>> transformer = Transformer.from_crs(4326, 32618, always_xy=True)
+    >>> shapely.transform(Point(-75, 50), transformer.transform, interleaved=False)
+    <POINT (500000 5538630.703)>
+
+    """
+    geometry_arr = np.array(geometry, dtype=np.object_)  # makes a copy
+    if include_z is None:
+        has_z = shapely.has_z(geometry_arr)
+        result = np.empty_like(geometry_arr)
+        result[has_z] = transform(
+            geometry_arr[has_z], transformation, include_z=True, interleaved=interleaved
+        )
+        result[~has_z] = transform(
+            geometry_arr[~has_z],
+            transformation,
+            include_z=False,
+            interleaved=interleaved,
+        )
+    else:
+        # TODO: expose include_m
+        include_m = False
+        coordinates = lib.get_coordinates(geometry_arr, include_z, include_m, False)
+        if interleaved:
+            new_coordinates = transformation(coordinates)
+        else:
+            new_coordinates = np.asarray(
+                transformation(*coordinates.T), dtype=np.float64
+            ).T
+        # check the array to yield understandable error messages
+        if not isinstance(new_coordinates, np.ndarray) or new_coordinates.ndim != 2:
+            raise ValueError(
+                "The provided transformation did not return a two-dimensional numpy "
+                "array"
+            )
+        if new_coordinates.dtype != np.float64:
+            raise ValueError(
+                "The provided transformation returned an array with an unexpected "
+                f"dtype ({new_coordinates.dtype})"
+            )
+        if new_coordinates.shape != coordinates.shape:
+            # if the shape is too small we will get a segfault
+            raise ValueError(
+                "The provided transformation returned an array with an unexpected "
+                f"shape ({new_coordinates.shape})"
+            )
+        result = lib.set_coordinates(geometry_arr, new_coordinates)
+    if result.ndim == 0 and not isinstance(geometry, np.ndarray):
+        return result.item()
+    return result
+
+
+def count_coordinates(geometry):
+    """Count the number of coordinate pairs in a geometry array.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries to count the coordinates of.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point
+    >>> shapely.count_coordinates(Point(0, 0))
+    1
+    >>> shapely.count_coordinates(LineString([(2, 2), (4, 2)]))
+    2
+    >>> shapely.count_coordinates(None)
+    0
+    >>> shapely.count_coordinates([Point(0, 0), None])
+    1
+
+    """
+    return lib.count_coordinates(np.asarray(geometry, dtype=np.object_))
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0: only supported XY and XYZ geometries
+#   get_coordinates(geometry, include_z=False, return_index=False)
+# shapely 2.1: shows deprecation warning about positional 'include_z' and 'return_index'
+#   get_coordinates(geometry, include_z=False, return_index=False, *, include_m=False)
+# shapely 2.2(?): enforce keyword-only arguments after 'geometry'
+#   get_coordinates(geometry, *, include_z=False, include_m=False, return_index=False)
+
+
+@deprecate_positional(["include_z", "return_index"], category=DeprecationWarning)
+def get_coordinates(geometry, include_z=False, return_index=False, *, include_m=False):
+    """Get coordinates from a geometry array as an array of floats.
+
+    The shape of the returned array is (N, 2), with N being the number of
+    coordinate pairs. The shape of the data may also be (N, 3) or (N, 4),
+    depending on ``include_z`` and ``include_m`` options.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries to get the coordinates of.
+    include_z, include_m : bool, default False
+        If both are False, return XY (2D) geometries.
+        If both are True, return XYZM (4D) geometries.
+        If either are True, return XYZ or XYM (3D) geometries.
+        If a geometry has no Z or M dimension, extra coordinate data will be NaN.
+
+        .. versionadded:: 2.1.0
+            The ``include_m`` parameter was added to support XYM (3D) and
+            XYZM (4D) geometries available with GEOS 3.12.0 or later.
+            With older GEOS versions, M dimension coordinates will be NaN.
+
+    return_index : bool, default False
+        If True, also return the index of each returned geometry as a separate
+        ndarray of integers. For multidimensional arrays, this indexes into the
+        flattened array (in C contiguous order).
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``include_z`` or ``return_index`` are
+        specified as positional arguments. In a future release, these will
+        need to be specified as keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point
+    >>> shapely.get_coordinates(Point(1, 2)).tolist()
+    [[1.0, 2.0]]
+    >>> shapely.get_coordinates(LineString([(2, 2), (4, 4)])).tolist()
+    [[2.0, 2.0], [4.0, 4.0]]
+    >>> shapely.get_coordinates(None)
+    array([], shape=(0, 2), dtype=float64)
+
+    By default the third dimension is ignored:
+
+    >>> shapely.get_coordinates(Point(1, 2, 3)).tolist()
+    [[1.0, 2.0]]
+    >>> shapely.get_coordinates(Point(1, 2, 3), include_z=True).tolist()
+    [[1.0, 2.0, 3.0]]
+
+    If geometries don't have Z or M dimension, these values will be NaN:
+
+    >>> pt = Point(1, 2)
+    >>> shapely.get_coordinates(pt, include_z=True).tolist()
+    [[1.0, 2.0, nan]]
+    >>> shapely.get_coordinates(pt, include_z=True, include_m=True).tolist()
+    [[1.0, 2.0, nan, nan]]
+
+    When ``return_index=True``, indexes are returned also:
+
+    >>> geometries = [LineString([(2, 2), (4, 4)]), Point(0, 0)]
+    >>> coordinates, index = shapely.get_coordinates(geometries, return_index=True)
+    >>> coordinates.tolist(), index.tolist()
+    ([[2.0, 2.0], [4.0, 4.0], [0.0, 0.0]], [0, 0, 1])
+
+    """
+    return lib.get_coordinates(
+        np.asarray(geometry, dtype=np.object_), include_z, include_m, return_index
+    )
+
+
+def set_coordinates(geometry, coordinates):
+    """Adapts the coordinates of a geometry array in-place.
+
+    If the coordinates array has shape (N, 2), all returned geometries
+    will be two-dimensional, and the third dimension will be discarded,
+    if present. If the coordinates array has shape (N, 3), the returned
+    geometries preserve the dimensionality of the input geometries.
+
+    .. warning::
+
+        The geometry array is modified in-place! If you do not want to
+        modify the original array, you can do
+        ``set_coordinates(arr.copy(), newcoords)``.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries to set the coordinates of.
+    coordinates: array_like
+        An array of coordinates to set.
+
+    See Also
+    --------
+    transform : Returns a copy of a geometry array with a function applied to its
+        coordinates.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point
+    >>> shapely.set_coordinates(Point(0, 0), [[1, 1]])
+    <POINT (1 1)>
+    >>> shapely.set_coordinates(
+    ...     [Point(0, 0), LineString([(0, 0), (0, 0)])],
+    ...     [[1, 2], [3, 4], [5, 6]]
+    ... ).tolist()
+    [<POINT (1 2)>, <LINESTRING (3 4, 5 6)>]
+    >>> shapely.set_coordinates([None, Point(0, 0)], [[1, 2]]).tolist()
+    [None, <POINT (1 2)>]
+
+    Third dimension of input geometry is discarded if coordinates array does
+    not include one:
+
+    >>> shapely.set_coordinates(Point(0, 0, 0), [[1, 1]])
+    <POINT (1 1)>
+    >>> shapely.set_coordinates(Point(0, 0, 0), [[1, 1, 1]])
+    <POINT Z (1 1 1)>
+
+    """
+    geometry_arr = np.asarray(geometry, dtype=np.object_)
+    coordinates = np.atleast_2d(np.asarray(coordinates)).astype(np.float64)
+    if coordinates.ndim != 2:
+        raise ValueError(
+            f"The coordinate array should have dimension of 2 (has {coordinates.ndim})"
+        )
+    n_coords = lib.count_coordinates(geometry_arr)
+    if (coordinates.shape[0] != n_coords) or (coordinates.shape[1] not in {2, 3}):
+        raise ValueError(
+            f"The coordinate array has an invalid shape {coordinates.shape}"
+        )
+    lib.set_coordinates(geometry_arr, coordinates)
+    if geometry_arr.ndim == 0 and not isinstance(geometry, np.ndarray):
+        return geometry_arr.item()
+    return geometry_arr

.venv/lib/python3.10/site-packages/shapely/coords.py

@@ -0,0 +1,119 @@
+"""Coordinate sequence utilities."""
+
+from array import array
+
+
+class CoordinateSequence:
+    """Access to coordinate tuples from the parent geometry's coordinate sequence.
+
+    Examples
+    --------
+    >>> from shapely.wkt import loads
+    >>> g = loads('POINT (0.0 0.0)')
+    >>> list(g.coords)
+    [(0.0, 0.0)]
+    >>> g = loads('POINT M (1 2 4)')
+    >>> g.coords[:]
+    [(1.0, 2.0, 4.0)]
+
+    """
+
+    def __init__(self, coords):
+        """Initialize the CoordinateSequence.
+
+        Parameters
+        ----------
+        coords : array
+            The coordinate array.
+
+        """
+        self._coords = coords
+
+    def __len__(self):
+        """Return the length of the CoordinateSequence.
+
+        Returns
+        -------
+        int
+            The length of the CoordinateSequence.
+
+        """
+        return self._coords.shape[0]
+
+    def __iter__(self):
+        """Iterate over the CoordinateSequence."""
+        for i in range(self.__len__()):
+            yield tuple(self._coords[i].tolist())
+
+    def __getitem__(self, key):
+        """Get the item at the specified index or slice.
+
+        Parameters
+        ----------
+        key : int or slice
+            The index or slice.
+
+        Returns
+        -------
+        tuple or list
+            The item at the specified index or slice.
+
+        """
+        m = self.__len__()
+        if isinstance(key, int):
+            if key + m < 0 or key >= m:
+                raise IndexError("index out of range")
+            if key < 0:
+                i = m + key
+            else:
+                i = key
+            return tuple(self._coords[i].tolist())
+        elif isinstance(key, slice):
+            res = []
+            start, stop, stride = key.indices(m)
+            for i in range(start, stop, stride):
+                res.append(tuple(self._coords[i].tolist()))
+            return res
+        else:
+            raise TypeError("key must be an index or slice")
+
+    def __array__(self, dtype=None, copy=None):
+        """Return a copy of the coordinate array.
+
+        Parameters
+        ----------
+        dtype : data-type, optional
+            The desired data-type for the array.
+        copy : bool, optional
+            If None (default) or True, a copy of the array is always returned.
+            If False, a ValueError is raised as this is not supported.
+
+        Returns
+        -------
+        array
+            The coordinate array.
+
+        Raises
+        ------
+        ValueError
+            If `copy=False` is specified.
+
+        """
+        if copy is False:
+            raise ValueError("`copy=False` isn't supported. A copy is always created.")
+        elif copy is True:
+            return self._coords.copy()
+        else:
+            return self._coords
+
+    @property
+    def xy(self):
+        """X and Y arrays."""
+        m = self.__len__()
+        x = array("d")
+        y = array("d")
+        for i in range(m):
+            xy = self._coords[i].tolist()
+            x.append(xy[0])
+            y.append(xy[1])
+        return x, y

.venv/lib/python3.10/site-packages/shapely/creation.py

@@ -0,0 +1,831 @@
+"""Methods to create geometries."""
+
+import numpy as np
+
+from shapely import Geometry, GeometryType, lib
+from shapely._enum import ParamEnum
+from shapely._geometry_helpers import collections_1d, simple_geometries_1d
+from shapely.decorators import deprecate_positional, multithreading_enabled
+from shapely.io import from_wkt
+
+__all__ = [
+    "box",
+    "destroy_prepared",
+    "empty",
+    "geometrycollections",
+    "linearrings",
+    "linestrings",
+    "multilinestrings",
+    "multipoints",
+    "multipolygons",
+    "points",
+    "polygons",
+    "prepare",
+]
+
+
+class HandleNaN(ParamEnum):
+    allow = 0
+    skip = 1
+    error = 2
+
+
+def _xyz_to_coords(x, y, z):
+    if y is None:
+        return x
+    if z is None:
+        coords = np.broadcast_arrays(x, y)
+    else:
+        coords = np.broadcast_arrays(x, y, z)
+    return np.stack(coords, axis=-1)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   points(coords, y=None, z=None, indices=None, out=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'indices' arg
+#   points(coords, y=None, z=None, indices=None, *, handle_nan=HandleNaN.allow, out=None, **kwargs)  # noqa: E501
+# shapely 2.2(?): enforce keyword-only arguments after 'z'
+#   points(coords, y=None, z=None, *, indices=None, handle_nan=HandleNaN.allow, out=None, **kwargs)  # noqa: E501
+
+
+@deprecate_positional(["indices"], category=DeprecationWarning)
+@multithreading_enabled
+def points(
+    coords,
+    y=None,
+    z=None,
+    indices=None,
+    *,
+    handle_nan=HandleNaN.allow,
+    out=None,
+    **kwargs,
+):
+    """Create an array of points.
+
+    Parameters
+    ----------
+    coords : array_like
+        An array of coordinate tuples (2- or 3-dimensional) or, if ``y`` is
+        provided, an array of x coordinates.
+    y : array_like, optional
+        An array of y coordinates.
+    z : array_like, optional
+        An array of z coordinates.
+    indices : array_like, optional
+        Indices into the target array where input coordinates belong. If
+        provided, the coords should be 2D with shape (N, 2) or (N, 3) and
+        indices should be an array of shape (N,) with integers in increasing
+        order. Missing indices result in a ValueError unless ``out`` is
+        provided, in which case the original value in ``out`` is kept.
+    handle_nan : shapely.HandleNaN or {'allow', 'skip', 'error'}, default 'allow'
+        Specifies what to do when a NaN or Inf is encountered in the coordinates:
+
+        - 'allow': the geometries are created with NaN or Inf coordinates.
+          Note that this can result in unexpected behaviour in subsequent
+          operations, and generally it is discouraged to have non-finite
+          coordinate values. One can use this option if you know all
+          coordinates are finite and want to avoid the overhead of checking
+          for this.
+        - 'skip': if any of x, y or z values are NaN or Inf, an empty point
+          will be created.
+        - 'error': if any NaN or Inf is detected in the coordinates, a ValueError
+          is raised. This option ensures that the created geometries have all
+          finite coordinate values.
+
+        .. versionadded:: 2.1.0
+    out : ndarray, optional
+        An array (with dtype object) to output the geometries into.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+        Ignored if ``indices`` is provided.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> shapely.points([[0, 1], [4, 5]]).tolist()
+    [<POINT (0 1)>, <POINT (4 5)>]
+    >>> shapely.points([0, 1, 2])
+    <POINT Z (0 1 2)>
+
+    Notes
+    -----
+    - GEOS 3.10, 3.11 and 3.12 automatically converts POINT (nan nan) to POINT EMPTY.
+    - GEOS 3.10 and 3.11 will transform a 3D point to 2D if its Z coordinate is NaN.
+    - Usage of the ``y`` and ``z`` arguments will prevents lazy evaluation in
+      ``dask``. Instead provide the coordinates as an array with shape
+      ``(..., 2)`` or ``(..., 3)`` using only the ``coords`` argument.
+
+    """
+    coords = _xyz_to_coords(coords, y, z)
+    if isinstance(handle_nan, str):
+        handle_nan = HandleNaN.get_value(handle_nan)
+    if indices is None:
+        return lib.points(coords, np.intc(handle_nan), out=out, **kwargs)
+    else:
+        return simple_geometries_1d(
+            coords, indices, GeometryType.POINT, handle_nan=handle_nan, out=out
+        )
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   linestrings(coords, y=None, z=None, indices=None, out=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'indices' arg
+#   linestrings(coords, y=None, z=None, indices=None, *, handle_nan=HandleNaN.allow, out=None, **kwargs)  # noqa: E501
+# shapely 2.2(?): enforce keyword-only arguments after 'z'
+#   linestrings(coords, y=None, z=None, *, indices=None, handle_nan=HandleNaN.allow, out=None, **kwargs)  # noqa: E501
+
+
+@deprecate_positional(["indices"], category=DeprecationWarning)
+@multithreading_enabled
+def linestrings(
+    coords,
+    y=None,
+    z=None,
+    indices=None,
+    *,
+    handle_nan=HandleNaN.allow,
+    out=None,
+    **kwargs,
+):
+    """Create an array of linestrings.
+
+    This function will raise an exception if a linestring contains less than
+    two points.
+
+    Parameters
+    ----------
+    coords : array_like
+        An array of lists of coordinate tuples (2- or 3-dimensional) or, if ``y``
+        is provided, an array of lists of x coordinates.
+    y : array_like, optional
+        An array of y coordinates.
+    z : array_like, optional
+        An array of z coordinates.
+    indices : array_like, optional
+        Indices into the target array where input coordinates belong. If
+        provided, the coords should be 2D with shape (N, 2) or (N, 3) and
+        indices should be an array of shape (N,) with integers in increasing
+        order. Missing indices result in a ValueError unless ``out`` is
+        provided, in which case the original value in ``out`` is kept.
+    handle_nan : shapely.HandleNaN or {'allow', 'skip', 'error'}, default 'allow'
+        Specifies what to do when a NaN or Inf is encountered in the coordinates:
+
+        - 'allow': the geometries are created with NaN or Inf coordinates.
+          Note that this can result in unexpected behaviour in subsequent
+          operations, and generally it is discouraged to have non-finite
+          coordinate values. One can use this option if you know all
+          coordinates are finite and want to avoid the overhead of checking
+          for this.
+        - 'skip': the coordinate pairs where any of x, y or z values are
+          NaN or Inf are ignored. If this results in ignoring all coordinates
+          for one geometry, an empty geometry is created.
+        - 'error': if any NaN or Inf is detected in the coordinates, a ValueError
+          is raised. This option ensures that the created geometries have all
+          finite coordinate values.
+
+        .. versionadded:: 2.1.0
+
+    out : ndarray, optional
+        An array (with dtype object) to output the geometries into.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+        Ignored if ``indices`` is provided.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> shapely.linestrings([[[0, 1], [4, 5]], [[2, 3], [5, 6]]]).tolist()
+    [<LINESTRING (0 1, 4 5)>, <LINESTRING (2 3, 5 6)>]
+    >>> shapely.linestrings(
+    ...     [[0, 1], [4, 5], [2, 3], [5, 6], [7, 8]],
+    ...     indices=[0, 0, 1, 1, 1]
+    ... ).tolist()
+    [<LINESTRING (0 1, 4 5)>, <LINESTRING (2 3, 5 6, 7 8)>]
+
+    Notes
+    -----
+    - Usage of the ``y`` and ``z`` arguments will prevents lazy evaluation in
+      ``dask``. Instead provide the coordinates as a ``(..., 2)`` or
+      ``(..., 3)`` array using only ``coords``.
+
+    """
+    coords = _xyz_to_coords(coords, y, z)
+    if isinstance(handle_nan, str):
+        handle_nan = HandleNaN.get_value(handle_nan)
+    if indices is None:
+        return lib.linestrings(coords, np.intc(handle_nan), out=out, **kwargs)
+    else:
+        return simple_geometries_1d(
+            coords, indices, GeometryType.LINESTRING, handle_nan=handle_nan, out=out
+        )
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   linearrings(coords, y=None, z=None, indices=None, out=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'indices' arg
+#   linearrings(coords, y=None, z=None, indices=None, *, handle_nan=HandleNaN.allow, out=None, **kwargs)  # noqa: E501
+# shapely 2.2(?): enforce keyword-only arguments after 'z'
+#   linearrings(coords, y=None, z=None, *, indices=None, handle_nan=HandleNaN.allow, out=None, **kwargs)  # noqa: E501
+
+
+@deprecate_positional(["indices"], category=DeprecationWarning)
+@multithreading_enabled
+def linearrings(
+    coords,
+    y=None,
+    z=None,
+    indices=None,
+    *,
+    handle_nan=HandleNaN.allow,
+    out=None,
+    **kwargs,
+):
+    """Create an array of linearrings.
+
+    If the provided coords do not constitute a closed linestring, or if there
+    are only 3 provided coords, the first
+    coordinate is duplicated at the end to close the ring. This function will
+    raise an exception if a linearring contains less than three points or if
+    the terminal coordinates contain NaN (not-a-number).
+
+    Parameters
+    ----------
+    coords : array_like
+        An array of lists of coordinate tuples (2- or 3-dimensional) or, if ``y``
+        is provided, an array of lists of x coordinates
+    y : array_like, optional
+        An array of y coordinates.
+    z : array_like, optional
+        An array of z coordinates.
+    indices : array_like, optional
+        Indices into the target array where input coordinates belong. If
+        provided, the coords should be 2D with shape (N, 2) or (N, 3) and
+        indices should be an array of shape (N,) with integers in increasing
+        order. Missing indices result in a ValueError unless ``out`` is
+        provided, in which case the original value in ``out`` is kept.
+    handle_nan : shapely.HandleNaN or {'allow', 'skip', 'error'}, default 'allow'
+        Specifies what to do when a NaN or Inf is encountered in the coordinates:
+
+        - 'allow': the geometries are created with NaN or Inf coordinates.
+          Note that this can result in unexpected behaviour in subsequent
+          operations, and generally it is discouraged to have non-finite
+          coordinate values. One can use this option if you know all
+          coordinates are finite and want to avoid the overhead of checking
+          for this.
+        - 'skip': the coordinate pairs where any of x, y or z values are
+          NaN or Inf are ignored. If this results in ignoring all coordinates
+          for one geometry, an empty geometry is created.
+        - 'error': if any NaN or Inf is detected in the coordinates, a ValueError
+          is raised. This option ensures that the created geometries have all
+          finite coordinate values.
+
+        .. versionadded:: 2.1.0
+
+    out : ndarray, optional
+        An array (with dtype object) to output the geometries into.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+        Ignored if ``indices`` is provided.
+
+    See Also
+    --------
+    linestrings
+
+    Examples
+    --------
+    >>> import shapely
+    >>> shapely.linearrings([[0, 0], [0, 1], [1, 1], [0, 0]])
+    <LINEARRING (0 0, 0 1, 1 1, 0 0)>
+    >>> shapely.linearrings([[0, 0], [0, 1], [1, 1]])
+    <LINEARRING (0 0, 0 1, 1 1, 0 0)>
+
+    Notes
+    -----
+    - Usage of the ``y`` and ``z`` arguments will prevents lazy evaluation in
+      ``dask``. Instead provide the coordinates as a ``(..., 2)`` or
+      ``(..., 3)`` array using only ``coords``.
+
+    """
+    coords = _xyz_to_coords(coords, y, z)
+    if isinstance(handle_nan, str):
+        handle_nan = HandleNaN.get_value(handle_nan)
+    if indices is None:
+        return lib.linearrings(coords, np.intc(handle_nan), out=out, **kwargs)
+    else:
+        return simple_geometries_1d(
+            coords, indices, GeometryType.LINEARRING, handle_nan=handle_nan, out=out
+        )
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   polygons(geometries, holes=None, indices=None, out=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'indices' arg
+#   polygons(geometries, holes=None, indices=None, *, out=None, **kwargs)
+# shapely 2.2(?): enforce keyword-only arguments after 'holes'
+#   polygons(geometries, holes=None, *, indices=None, out=None, **kwargs)
+
+
+@deprecate_positional(["indices"], category=DeprecationWarning)
+@multithreading_enabled
+def polygons(geometries, holes=None, indices=None, *, out=None, **kwargs):
+    """Create an array of polygons.
+
+    Parameters
+    ----------
+    geometries : array_like
+        An array of linearrings or coordinates (see linearrings).
+        Unless ``indices`` are given (see description below), this
+        include the outer shells only. The ``holes`` argument should be used
+        to create polygons with holes.
+    holes : array_like, optional
+        An array of lists of linearrings that constitute holes for each shell.
+        Not to be used in combination with ``indices``.
+    indices : array_like, optional
+        Indices into the target array where input geometries belong. If
+        provided, the holes are expected to be present inside ``geometries``;
+        the first geometry for each index is the outer shell
+        and all subsequent geometries in that index are the holes.
+        Both geometries and indices should be 1D and have matching sizes.
+        Indices should be in increasing order. Missing indices result in a
+        ValueError unless ``out`` is  provided, in which case the original value
+        in ``out`` is kept.
+    out : ndarray, optional
+        An array (with dtype object) to output the geometries into.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+        Ignored if ``indices`` is provided.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``indices`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    Examples
+    --------
+    >>> import shapely
+
+    Polygons are constructed from rings:
+
+    >>> ring_1 = shapely.linearrings([[0, 0], [0, 10], [10, 10], [10, 0]])
+    >>> ring_2 = shapely.linearrings([[2, 6], [2, 7], [3, 7], [3, 6]])
+    >>> shapely.polygons([ring_1, ring_2])[0]
+    <POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0))>
+    >>> shapely.polygons([ring_1, ring_2])[1]
+    <POLYGON ((2 6, 2 7, 3 7, 3 6, 2 6))>
+
+    Or from coordinates directly:
+
+    >>> shapely.polygons([[0, 0], [0, 10], [10, 10], [10, 0]])
+    <POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0))>
+
+    Adding holes can be done using the ``holes`` keyword argument:
+
+    >>> shapely.polygons(ring_1, holes=[ring_2])
+    <POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0), (2 6, 2 7, 3 7, 3 6, 2 6))>
+
+    Or using the ``indices`` argument:
+
+    >>> shapely.polygons([ring_1, ring_2], indices=[0, 1])[0]
+    <POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0))>
+    >>> shapely.polygons([ring_1, ring_2], indices=[0, 1])[1]
+    <POLYGON ((2 6, 2 7, 3 7, 3 6, 2 6))>
+    >>> shapely.polygons([ring_1, ring_2], indices=[0, 0])[0]
+    <POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0), (2 6, 2 7, 3 7, 3 6, 2 6))>
+
+    Missing input values (``None``) are skipped and may result in an
+    empty polygon:
+
+    >>> shapely.polygons(None)
+    <POLYGON EMPTY>
+    >>> shapely.polygons(ring_1, holes=[None])
+    <POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0))>
+    >>> shapely.polygons([ring_1, None], indices=[0, 0])[0]
+    <POLYGON ((0 0, 0 10, 10 10, 10 0, 0 0))>
+
+    """
+    geometries = np.asarray(geometries)
+    if not isinstance(geometries, Geometry) and np.issubdtype(
+        geometries.dtype, np.number
+    ):
+        geometries = linearrings(geometries)
+
+    if indices is not None:
+        if holes is not None:
+            raise TypeError("Cannot specify separate holes array when using indices.")
+        return collections_1d(geometries, indices, GeometryType.POLYGON, out=out)
+
+    if holes is None:
+        # no holes provided: initialize an empty holes array matching shells
+        shape = geometries.shape + (0,) if isinstance(geometries, np.ndarray) else (0,)
+        holes = np.empty(shape, dtype=object)
+    else:
+        holes = np.asarray(holes)
+        # convert holes coordinates into linearrings
+        if np.issubdtype(holes.dtype, np.number):
+            holes = linearrings(holes)
+
+    return lib.polygons(geometries, holes, out=out, **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   box(xmin, ymin, xmax, ymax, ccw=True, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'ccw' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'ymax'
+#   box(xmin, ymin, xmax, ymax, *, ccw=True, **kwargs)
+
+
+@deprecate_positional(["ccw"], category=DeprecationWarning)
+@multithreading_enabled
+def box(xmin, ymin, xmax, ymax, ccw=True, **kwargs):
+    """Create box polygons.
+
+    Parameters
+    ----------
+    xmin : float or array_like
+        Float or array of minimum x coordinates.
+    ymin : float or array_like
+        Float or array of minimum y coordinates.
+    xmax : float or array_like
+        Float or array of maximum x coordinates.
+    ymax : float or array_like
+        Float or array of maximum y coordinates.
+    ccw : bool, default True
+        If True, box will be created in counterclockwise direction starting
+        from bottom right coordinate (xmax, ymin).
+        If False, box will be created in clockwise direction starting from
+        bottom left coordinate (xmin, ymin).
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``ccw`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> shapely.box(0, 0, 1, 1)
+    <POLYGON ((1 0, 1 1, 0 1, 0 0, 1 0))>
+    >>> shapely.box(0, 0, 1, 1, ccw=False)
+    <POLYGON ((0 0, 0 1, 1 1, 1 0, 0 0))>
+
+    """
+    return lib.box(xmin, ymin, xmax, ymax, ccw, **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   multipoints(geometries, indices=None, out=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'indices' arg
+#   multipoints(geometries, indices=None, *, out=None, **kwargs)
+# shapely 2.2(?): enforce keyword-only arguments after 'indices'
+#   multipoints(geometries, *, indices=None, out=None, **kwargs)
+
+
+@deprecate_positional(["indices"], category=DeprecationWarning)
+@multithreading_enabled
+def multipoints(geometries, indices=None, *, out=None, **kwargs):
+    """Create multipoints from arrays of points.
+
+    Parameters
+    ----------
+    geometries : array_like
+        An array of points or coordinates (see points).
+    indices : array_like, optional
+        Indices into the target array where input geometries belong. If
+        provided, both geometries and indices should be 1D and have matching
+        sizes. Indices should be in increasing order. Missing indices result
+        in a ValueError unless ``out`` is  provided, in which case the original
+        value in ``out`` is kept.
+    out : ndarray, optional
+        An array (with dtype object) to output the geometries into.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+        Ignored if ``indices`` is provided.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``indices`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    Examples
+    --------
+    >>> import shapely
+
+    Multipoints are constructed from points:
+
+    >>> point_1 = shapely.points([1, 1])
+    >>> point_2 = shapely.points([2, 2])
+    >>> shapely.multipoints([point_1, point_2])
+    <MULTIPOINT ((1 1), (2 2))>
+    >>> shapely.multipoints([[point_1, point_2], [point_2, None]]).tolist()
+    [<MULTIPOINT ((1 1), (2 2))>, <MULTIPOINT ((2 2))>]
+
+    Or from coordinates directly:
+
+    >>> shapely.multipoints([[0, 0], [2, 2], [3, 3]])
+    <MULTIPOINT ((0 0), (2 2), (3 3))>
+
+    Multiple multipoints of different sizes can be constructed efficiently using the
+    ``indices`` keyword argument:
+
+    >>> shapely.multipoints([point_1, point_2, point_2], indices=[0, 0, 1]).tolist()
+    [<MULTIPOINT ((1 1), (2 2))>, <MULTIPOINT ((2 2))>]
+
+    Missing input values (``None``) are skipped and may result in an
+    empty multipoint:
+
+    >>> shapely.multipoints([None])
+    <MULTIPOINT EMPTY>
+    >>> shapely.multipoints([point_1, None], indices=[0, 0]).tolist()
+    [<MULTIPOINT ((1 1))>]
+    >>> shapely.multipoints([point_1, None], indices=[0, 1]).tolist()
+    [<MULTIPOINT ((1 1))>, <MULTIPOINT EMPTY>]
+
+    """
+    typ = GeometryType.MULTIPOINT
+    geometries = np.asarray(geometries)
+    if not isinstance(geometries, Geometry) and np.issubdtype(
+        geometries.dtype, np.number
+    ):
+        geometries = points(geometries)
+    if indices is None:
+        return lib.create_collection(geometries, np.intc(typ), out=out, **kwargs)
+    else:
+        return collections_1d(geometries, indices, typ, out=out)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   multilinestrings(geometries, indices=None, out=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'indices' arg
+#   multilinestrings(geometries, indices=None, *, out=None, **kwargs)
+# shapely 2.2(?): enforce keyword-only arguments after 'indices'
+#   multilinestrings(geometries, *, indices=None, out=None, **kwargs)
+
+
+@deprecate_positional(["indices"], category=DeprecationWarning)
+@multithreading_enabled
+def multilinestrings(geometries, indices=None, *, out=None, **kwargs):
+    """Create multilinestrings from arrays of linestrings.
+
+    Parameters
+    ----------
+    geometries : array_like
+        An array of linestrings or coordinates (see linestrings).
+    indices : array_like, optional
+        Indices into the target array where input geometries belong. If
+        provided, both geometries and indices should be 1D and have matching
+        sizes. Indices should be in increasing order. Missing indices result
+        in a ValueError unless ``out`` is  provided, in which case the original
+        value in ``out`` is kept.
+    out : ndarray, optional
+        An array (with dtype object) to output the geometries into.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+        Ignored if ``indices`` is provided.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``indices`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    multipoints
+
+    """
+    typ = GeometryType.MULTILINESTRING
+    geometries = np.asarray(geometries)
+    if not isinstance(geometries, Geometry) and np.issubdtype(
+        geometries.dtype, np.number
+    ):
+        geometries = linestrings(geometries)
+
+    if indices is None:
+        return lib.create_collection(geometries, np.intc(typ), out=out, **kwargs)
+    else:
+        return collections_1d(geometries, indices, typ, out=out)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   multipolygons(geometries, indices=None, out=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'indices' arg
+#   multipolygons(geometries, indices=None, *, out=None, **kwargs)
+# shapely 2.2(?): enforce keyword-only arguments after 'indices'
+#   multipolygons(geometries, *, indices=None, out=None, **kwargs)
+
+
+@deprecate_positional(["indices"], category=DeprecationWarning)
+@multithreading_enabled
+def multipolygons(geometries, indices=None, *, out=None, **kwargs):
+    """Create multipolygons from arrays of polygons.
+
+    Parameters
+    ----------
+    geometries : array_like
+        An array of polygons or coordinates (see polygons).
+    indices : array_like, optional
+        Indices into the target array where input geometries belong. If
+        provided, both geometries and indices should be 1D and have matching
+        sizes. Indices should be in increasing order. Missing indices result
+        in a ValueError unless ``out`` is  provided, in which case the original
+        value in ``out`` is kept.
+    out : ndarray, optional
+        An array (with dtype object) to output the geometries into.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+        Ignored if ``indices`` is provided.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``indices`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    multipoints
+
+    """
+    typ = GeometryType.MULTIPOLYGON
+    geometries = np.asarray(geometries)
+    if not isinstance(geometries, Geometry) and np.issubdtype(
+        geometries.dtype, np.number
+    ):
+        geometries = polygons(geometries)
+    if indices is None:
+        return lib.create_collection(geometries, np.intc(typ), out=out, **kwargs)
+    else:
+        return collections_1d(geometries, indices, typ, out=out)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   geometrycollections(geometries, indices=None, out=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'indices' arg
+#   geometrycollections(geometries, indices=None, *, out=None, **kwargs)
+# shapely 2.2(?): enforce keyword-only arguments after 'indices'
+#   geometrycollections(geometries, *, indices=None, out=None, **kwargs)
+
+
+@deprecate_positional(["indices"], category=DeprecationWarning)
+@multithreading_enabled
+def geometrycollections(geometries, indices=None, out=None, **kwargs):
+    """Create geometrycollections from arrays of geometries.
+
+    Parameters
+    ----------
+    geometries : array_like
+        An array of geometries.
+    indices : array_like, optional
+        Indices into the target array where input geometries belong. If
+        provided, both geometries and indices should be 1D and have matching
+        sizes. Indices should be in increasing order. Missing indices result
+        in a ValueError unless ``out`` is  provided, in which case the original
+        value in ``out`` is kept.
+    out : ndarray, optional
+        An array (with dtype object) to output the geometries into.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+        Ignored if ``indices`` is provided.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``indices`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    multipoints
+
+    """
+    typ = GeometryType.GEOMETRYCOLLECTION
+    if indices is None:
+        return lib.create_collection(geometries, np.intc(typ), out=out, **kwargs)
+    else:
+        return collections_1d(geometries, indices, typ, out=out)
+
+
+def prepare(geometry, **kwargs):
+    """Prepare a geometry, improving performance of other operations.
+
+    A prepared geometry is a normal geometry with added information such as an
+    index on the line segments. This improves the performance of the following
+    operations: contains, contains_properly, covered_by, covers, crosses,
+    disjoint, intersects, overlaps, touches, and within.
+
+    Note that if a prepared geometry is modified, the newly created Geometry
+    object is not prepared. In that case, ``prepare`` should be called again.
+
+    This function does not recompute previously prepared geometries;
+    it is efficient to call this function on an array that partially contains
+    prepared geometries.
+
+    This function does not return any values; geometries are modified in place.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometries are changed in place
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    See Also
+    --------
+    is_prepared : Identify whether a geometry is prepared already.
+    destroy_prepared : Destroy the prepared part of a geometry.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Point
+    >>> poly = shapely.buffer(Point(1.0, 1.0), 1)
+    >>> shapely.prepare(poly)
+    >>> shapely.contains_properly(poly, [Point(0.0, 0.0), Point(0.5, 0.5)]).tolist()
+    [False, True]
+
+    """
+    lib.prepare(geometry, **kwargs)
+
+
+def destroy_prepared(geometry, **kwargs):
+    """Destroy the prepared part of a geometry, freeing up memory.
+
+    Note that the prepared geometry will always be cleaned up if the geometry itself
+    is dereferenced. This function needs only be called in very specific circumstances,
+    such as freeing up memory without losing the geometries, or benchmarking.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometries are changed in-place
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    See Also
+    --------
+    prepare
+
+    """
+    lib.destroy_prepared(geometry, **kwargs)
+
+
+def empty(shape, geom_type=None, order="C"):
+    """Create a geometry array prefilled with None or with empty geometries.
+
+    Parameters
+    ----------
+    shape : int or tuple of int
+        Shape of the empty array, e.g., ``(2, 3)`` or ``2``.
+    geom_type : shapely.GeometryType, optional
+        The desired geometry type in case the array should be prefilled
+        with empty geometries. Default ``None``.
+    order : {'C', 'F'}, optional, default: 'C'
+        Whether to store multi-dimensional data in row-major
+        (C-style) or column-major (Fortran-style) order in
+        memory.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> shapely.empty((2, 3)).tolist()
+    [[None, None, None], [None, None, None]]
+    >>> shapely.empty(2, geom_type=shapely.GeometryType.POINT).tolist()
+    [<POINT EMPTY>, <POINT EMPTY>]
+
+    """
+    if geom_type is None:
+        return np.empty(shape, dtype=object, order=order)
+
+    geom_type = GeometryType(geom_type)  # cast int to GeometryType
+    if geom_type is GeometryType.MISSING:
+        return np.empty(shape, dtype=object, order=order)
+
+    fill_value = from_wkt(geom_type.name + " EMPTY")
+    return np.full(shape, fill_value, dtype=object, order=order)

.venv/lib/python3.10/site-packages/shapely/decorators.py

@@ -0,0 +1,133 @@
+"""Decorators for Shapely functions."""
+
+import inspect
+import os
+import warnings
+from functools import wraps
+
+import numpy as np
+
+from shapely import lib
+from shapely.errors import UnsupportedGEOSVersionError
+
+
+class requires_geos:
+    """Decorator to require a minimum GEOS version."""
+
+    def __init__(self, version):
+        """Create a decorator that requires a minimum GEOS version."""
+        if version.count(".") != 2:
+            raise ValueError("Version must be <major>.<minor>.<patch> format")
+        self.version = tuple(int(x) for x in version.split("."))
+
+    def __call__(self, func):
+        """Return the wrapped function."""
+        is_compatible = lib.geos_version >= self.version
+        is_doc_build = os.environ.get("SPHINX_DOC_BUILD") == "1"  # set in docs/conf.py
+        if is_compatible and not is_doc_build:
+            return func  # return directly, do not change the docstring
+
+        msg = "'{}' requires at least GEOS {}.{}.{}.".format(
+            func.__name__, *self.version
+        )
+        if is_compatible:
+
+            @wraps(func)
+            def wrapped(*args, **kwargs):
+                return func(*args, **kwargs)
+
+        else:
+
+            @wraps(func)
+            def wrapped(*args, **kwargs):
+                raise UnsupportedGEOSVersionError(msg)
+
+        doc = wrapped.__doc__
+        if doc:
+            # Insert the message at the first double newline
+            position = doc.find("\n\n") + 2
+            # Figure out the indentation level
+            indent = 0
+            while True:
+                if doc[position + indent] == " ":
+                    indent += 1
+                else:
+                    break
+            wrapped.__doc__ = doc.replace(
+                "\n\n", "\n\n{}.. note:: {}\n\n".format(" " * indent, msg), 1
+            )
+
+        return wrapped
+
+
+def multithreading_enabled(func):
+    """Enable multithreading.
+
+    To do this, the writable flags of object type ndarrays are set to False.
+
+    NB: multithreading also requires the GIL to be released, which is done in
+    the C extension (ufuncs.c).
+    """
+
+    @wraps(func)
+    def wrapped(*args, **kwargs):
+        array_args = [
+            arg for arg in args if isinstance(arg, np.ndarray) and arg.dtype == object
+        ] + [
+            arg
+            for name, arg in kwargs.items()
+            if name not in {"where", "out"}
+            and isinstance(arg, np.ndarray)
+            and arg.dtype == object
+        ]
+        old_flags = [arr.flags.writeable for arr in array_args]
+        try:
+            for arr in array_args:
+                arr.flags.writeable = False
+            return func(*args, **kwargs)
+        finally:
+            for arr, old_flag in zip(array_args, old_flags):
+                arr.flags.writeable = old_flag
+
+    return wrapped
+
+
+def deprecate_positional(should_be_kwargs, category=DeprecationWarning):
+    """Show warning if positional arguments are used that should be keyword."""
+
+    def decorator(func):
+        @wraps(func)
+        def wrapper(*args, **kwargs):
+            # call the function first, to make sure the signature matches
+            ret_value = func(*args, **kwargs)
+
+            # check signature to see which positional args were used
+            sig = inspect.signature(func)
+            args_bind = sig.bind_partial(*args)
+            warn_args = [
+                f"`{arg}`"
+                for arg in args_bind.arguments.keys()
+                if arg in should_be_kwargs
+            ]
+            if warn_args:
+                if len(warn_args) == 1:
+                    plr = ""
+                    isare = "is"
+                    args = warn_args[0]
+                else:
+                    plr = "s"
+                    isare = "are"
+                    if len(warn_args) < 3:
+                        args = " and ".join(warn_args)
+                    else:
+                        args = ", ".join(warn_args[:-1]) + ", and " + warn_args[-1]
+                msg = (
+                    f"positional argument{plr} {args} for `{func.__name__}` "
+                    f"{isare} deprecated.  Please use keyword argument{plr} instead."
+                )
+                warnings.warn(msg, category=category, stacklevel=2)
+            return ret_value
+
+        return wrapper
+
+    return decorator

.venv/lib/python3.10/site-packages/shapely/errors.py

@@ -0,0 +1,79 @@
+"""Shapely errors."""
+
+import threading
+
+from shapely.lib import GEOSException, ShapelyError, _setup_signal_checks  # noqa: F401
+
+
+def setup_signal_checks(interval=10000):
+    """Enable Python signal checks in the ufunc inner loops.
+
+    Doing so allows termination (using CTRL+C) of operations on large arrays of
+    vectors.
+
+    Parameters
+    ----------
+    interval : int, default 10000
+        Check for interrupts every x iterations. The higher the number, the
+        slower shapely will respond to a signal. However, at low values there
+        will be a negative effect on performance. The default of 10000 does not
+        have any measureable effects on performance.
+
+    Notes
+    -----
+    For more information on signals consult the Python docs:
+
+    https://docs.python.org/3/library/signal.html
+
+    """
+    if interval <= 0:
+        raise ValueError("Signal checks interval must be greater than zero.")
+
+    _setup_signal_checks(interval, threading.main_thread().ident)
+
+
+class UnsupportedGEOSVersionError(ShapelyError):
+    """Raised when the GEOS library version does not support a certain operation."""
+
+
+class DimensionError(ShapelyError):
+    """An error in the number of coordinate dimensions."""
+
+
+class TopologicalError(ShapelyError):
+    """A geometry is invalid or topologically incorrect."""
+
+
+class ShapelyDeprecationWarning(FutureWarning):
+    """Warning for features that will be removed or changed in a future release."""
+
+
+class EmptyPartError(ShapelyError):
+    """An error signifying an empty part was encountered when creating a multi-part."""
+
+
+class GeometryTypeError(ShapelyError):
+    """An error raised when the geometry has an unrecognized or inappropriate type."""
+
+
+def __getattr__(name):
+    import warnings
+
+    # Alias Shapely 1.8 error classes to ShapelyError with deprecation warning
+    if name in [
+        "ReadingError",
+        "WKBReadingError",
+        "WKTReadingError",
+        "PredicateError",
+        "InvalidGeometryError",
+    ]:
+        warnings.warn(
+            f"{name} is deprecated and will be removed in a future version. "
+            "Use ShapelyError instead (functions previously raising {name} "
+            "will now raise a ShapelyError instead).",
+            FutureWarning,
+            stacklevel=2,
+        )
+        return ShapelyError
+
+    raise AttributeError(f"module 'shapely.errors' has no attribute '{name}'")

.venv/lib/python3.10/site-packages/shapely/geometry/__init__.py

@@ -0,0 +1,27 @@
+"""Geometry classes and factories."""
+
+from shapely.geometry.base import CAP_STYLE, JOIN_STYLE
+from shapely.geometry.collection import GeometryCollection
+from shapely.geometry.geo import box, mapping, shape
+from shapely.geometry.linestring import LineString
+from shapely.geometry.multilinestring import MultiLineString
+from shapely.geometry.multipoint import MultiPoint
+from shapely.geometry.multipolygon import MultiPolygon
+from shapely.geometry.point import Point
+from shapely.geometry.polygon import LinearRing, Polygon
+
+__all__ = [
+    "CAP_STYLE",
+    "JOIN_STYLE",
+    "GeometryCollection",
+    "LineString",
+    "LinearRing",
+    "MultiLineString",
+    "MultiPoint",
+    "MultiPolygon",
+    "Point",
+    "Polygon",
+    "box",
+    "mapping",
+    "shape",
+]

.venv/lib/python3.10/site-packages/shapely/geometry/collection.py

@@ -0,0 +1,58 @@
+"""Multi-part collections of geometries."""
+
+import shapely
+from shapely.geometry.base import BaseGeometry, BaseMultipartGeometry
+
+
+class GeometryCollection(BaseMultipartGeometry):
+    """Collection of one or more geometries that can be of different types.
+
+    Parameters
+    ----------
+    geoms : list
+        A list of shapely geometry instances, which may be of varying geometry
+        types.
+
+    Attributes
+    ----------
+    geoms : sequence
+        A sequence of Shapely geometry instances
+
+    Examples
+    --------
+    Create a GeometryCollection with a Point and a LineString
+
+    >>> from shapely import GeometryCollection, LineString, Point
+    >>> p = Point(51, -1)
+    >>> l = LineString([(52, -1), (49, 2)])
+    >>> gc = GeometryCollection([p, l])
+
+    """
+
+    __slots__ = []
+
+    def __new__(self, geoms=None):
+        """Create a new GeometryCollection."""
+        if isinstance(geoms, BaseGeometry):
+            # TODO(shapely-2.0) do we actually want to split Multi-part geometries?
+            # this is needed for the split() tests
+            if hasattr(geoms, "geoms"):
+                geoms = geoms.geoms
+            else:
+                geoms = [geoms]
+        elif geoms is None or len(geoms) == 0:
+            # TODO better empty constructor
+            return shapely.from_wkt("GEOMETRYCOLLECTION EMPTY")
+
+        return shapely.geometrycollections(geoms)
+
+    @property
+    def __geo_interface__(self):
+        """Return a GeoJSON-like mapping of the geometry collection."""
+        geometries = []
+        for geom in self.geoms:
+            geometries.append(geom.__geo_interface__)
+        return dict(type="GeometryCollection", geometries=geometries)
+
+
+shapely.lib.registry[7] = GeometryCollection

.venv/lib/python3.10/site-packages/shapely/geometry/geo.py

@@ -0,0 +1,143 @@
+"""Geometry factories based on the geo interface."""
+
+import numpy as np
+
+from shapely.errors import GeometryTypeError
+from shapely.geometry.collection import GeometryCollection
+from shapely.geometry.linestring import LineString
+from shapely.geometry.multilinestring import MultiLineString
+from shapely.geometry.multipoint import MultiPoint
+from shapely.geometry.multipolygon import MultiPolygon
+from shapely.geometry.point import Point
+from shapely.geometry.polygon import LinearRing, Polygon
+
+
+def _is_coordinates_empty(coordinates):
+    """Identify if coordinates or subset of coordinates are empty."""
+    if coordinates is None:
+        return True
+
+    if isinstance(coordinates, (list, tuple, np.ndarray)):
+        if len(coordinates) == 0:
+            return True
+        return all(map(_is_coordinates_empty, coordinates))
+    else:
+        return False
+
+
+def _empty_shape_for_no_coordinates(geom_type):
+    """Return empty counterpart for geom_type."""
+    if geom_type == "point":
+        return Point()
+    elif geom_type == "multipoint":
+        return MultiPoint()
+    elif geom_type == "linestring":
+        return LineString()
+    elif geom_type == "multilinestring":
+        return MultiLineString()
+    elif geom_type == "polygon":
+        return Polygon()
+    elif geom_type == "multipolygon":
+        return MultiPolygon()
+    else:
+        raise GeometryTypeError(f"Unknown geometry type: {geom_type!r}")
+
+
+def box(minx, miny, maxx, maxy, ccw=True):
+    """Return a rectangular polygon with configurable normal vector."""
+    coords = [(maxx, miny), (maxx, maxy), (minx, maxy), (minx, miny)]
+    if not ccw:
+        coords = coords[::-1]
+    return Polygon(coords)
+
+
+def shape(context):
+    """Return a new, independent geometry with coordinates copied from the context.
+
+    Changes to the original context will not be reflected in the geometry
+    object.
+
+    Parameters
+    ----------
+    context :
+        a GeoJSON-like dict, which provides a "type" member describing the type
+        of the geometry and "coordinates" member providing a list of coordinates,
+        or an object which implements __geo_interface__.
+
+    Returns
+    -------
+    Geometry object
+
+    Examples
+    --------
+    Create a Point from GeoJSON, and then create a copy using __geo_interface__.
+
+    >>> from shapely.geometry import shape
+    >>> context = {'type': 'Point', 'coordinates': [0, 1]}
+    >>> geom = shape(context)
+    >>> geom.geom_type == 'Point'
+    True
+    >>> geom.wkt
+    'POINT (0 1)'
+    >>> geom2 = shape(geom)
+    >>> geom == geom2
+    True
+
+    """
+    if hasattr(context, "__geo_interface__"):
+        ob = context.__geo_interface__
+    else:
+        ob = context
+    geom_type = ob.get("type").lower()
+
+    if geom_type == "feature":
+        # GeoJSON features must have a 'geometry' field.
+        ob = ob["geometry"]
+        geom_type = ob.get("type").lower()
+
+    if "coordinates" in ob and _is_coordinates_empty(ob["coordinates"]):
+        return _empty_shape_for_no_coordinates(geom_type)
+    elif geom_type == "point":
+        return Point(ob["coordinates"])
+    elif geom_type == "linestring":
+        return LineString(ob["coordinates"])
+    elif geom_type == "linearring":
+        return LinearRing(ob["coordinates"])
+    elif geom_type == "polygon":
+        return Polygon(ob["coordinates"][0], ob["coordinates"][1:])
+    elif geom_type == "multipoint":
+        return MultiPoint(ob["coordinates"])
+    elif geom_type == "multilinestring":
+        return MultiLineString(ob["coordinates"])
+    elif geom_type == "multipolygon":
+        return MultiPolygon([[c[0], c[1:]] for c in ob["coordinates"]])
+    elif geom_type == "geometrycollection":
+        geoms = [shape(g) for g in ob.get("geometries", [])]
+        return GeometryCollection(geoms)
+    else:
+        raise GeometryTypeError(f"Unknown geometry type: {geom_type!r}")
+
+
+def mapping(ob):
+    """Return a GeoJSON-like mapping.
+
+    Input should be a Geometry or an object which implements __geo_interface__.
+
+    Parameters
+    ----------
+    ob : geometry or object
+        An object which implements __geo_interface__.
+
+    Returns
+    -------
+    dict
+
+    Examples
+    --------
+    >>> from shapely.geometry import mapping, Point
+    >>> pt = Point(0, 0)
+    >>> mapping(pt)
+    {'type': 'Point', 'coordinates': (0.0, 0.0)}
+
+    """
+    return ob.__geo_interface__

.venv/lib/python3.10/site-packages/shapely/geometry/linestring.py

@@ -0,0 +1,211 @@
+"""Line strings and related utilities."""
+
+import numpy as np
+
+import shapely
+from shapely.decorators import deprecate_positional
+from shapely.geometry.base import JOIN_STYLE, BaseGeometry
+from shapely.geometry.point import Point
+
+__all__ = ["LineString"]
+
+
+class LineString(BaseGeometry):
+    """A geometry type composed of one or more line segments.
+
+    A LineString is a one-dimensional feature and has a non-zero length but
+    zero area. It may approximate a curve and need not be straight. A LineString may
+    be closed.
+
+    Parameters
+    ----------
+    coordinates : sequence
+        A sequence of (x, y, [,z]) numeric coordinate pairs or triples, or
+        an array-like with shape (N, 2) or (N, 3).
+        Also can be a sequence of Point objects, or combination of both.
+
+    Examples
+    --------
+    Create a LineString with two segments
+
+    >>> from shapely import LineString
+    >>> a = LineString([[0, 0], [1, 0], [1, 1]])
+    >>> a.length
+    2.0
+
+    """
+
+    __slots__ = []
+
+    def __new__(self, coordinates=None):
+        """Create a new LineString geometry."""
+        if coordinates is None:
+            # empty geometry
+            # TODO better constructor
+            return shapely.from_wkt("LINESTRING EMPTY")
+        elif isinstance(coordinates, LineString):
+            if type(coordinates) is LineString:
+                # return original objects since geometries are immutable
+                return coordinates
+            else:
+                # LinearRing
+                # TODO convert LinearRing to LineString more directly
+                coordinates = coordinates.coords
+        else:
+            if hasattr(coordinates, "__array__"):
+                coordinates = np.asarray(coordinates)
+            if isinstance(coordinates, np.ndarray) and np.issubdtype(
+                coordinates.dtype, np.number
+            ):
+                pass
+            else:
+                # check coordinates on points
+                def _coords(o):
+                    if isinstance(o, Point):
+                        return o.coords[0]
+                    else:
+                        return [float(c) for c in o]
+
+                coordinates = [_coords(o) for o in coordinates]
+
+        if len(coordinates) == 0:
+            # empty geometry
+            # TODO better constructor + should shapely.linestrings handle this?
+            return shapely.from_wkt("LINESTRING EMPTY")
+
+        geom = shapely.linestrings(coordinates)
+        if not isinstance(geom, LineString):
+            raise ValueError("Invalid values passed to LineString constructor")
+        return geom
+
+    @property
+    def __geo_interface__(self):
+        """Return a GeoJSON-like mapping of the LineString geometry."""
+        return {"type": "LineString", "coordinates": tuple(self.coords)}
+
+    def svg(self, scale_factor=1.0, stroke_color=None, opacity=None):
+        """Return SVG polyline element for the LineString geometry.
+
+        Parameters
+        ----------
+        scale_factor : float
+            Multiplication factor for the SVG stroke-width.  Default is 1.
+        stroke_color : str, optional
+            Hex string for stroke color. Default is to use "#66cc99" if
+            geometry is valid, and "#ff3333" if invalid.
+        opacity : float
+            Float number between 0 and 1 for color opacity. Default value is 0.8
+
+        """
+        if self.is_empty:
+            return "<g />"
+        if stroke_color is None:
+            stroke_color = "#66cc99" if self.is_valid else "#ff3333"
+        if opacity is None:
+            opacity = 0.8
+        pnt_format = " ".join(["{},{}".format(*c) for c in self.coords])
+        return (
+            f'<polyline fill="none" stroke="{stroke_color}" '
+            f'stroke-width="{2.0 * scale_factor}" '
+            f'points="{pnt_format}" opacity="{opacity}" />'
+        )
+
+    @property
+    def xy(self):
+        """Separate arrays of X and Y coordinate values.
+
+        Examples
+        --------
+        >>> from shapely import LineString
+        >>> x, y = LineString([(0, 0), (1, 1)]).xy
+        >>> list(x)
+        [0.0, 1.0]
+        >>> list(y)
+        [0.0, 1.0]
+
+        """
+        return self.coords.xy
+
+    # Note: future plan is to change this signature over a few releases:
+    # shapely 2.0:
+    #   offset_curve(self, distance, quad_segs=16, ...)
+    # shapely 2.1: shows deprecation warning about positional 'quad_segs', etc.
+    #   same signature as 2.0
+    # shapely 2.2(?): enforce keyword-only arguments after 'distance'
+    #   offset_curve(self, distance, *, quad_segs=16, ...)
+
+    @deprecate_positional(
+        ["quad_segs", "join_style", "mitre_limit"], category=DeprecationWarning
+    )
+    def offset_curve(
+        self,
+        distance,
+        quad_segs=16,
+        join_style=JOIN_STYLE.round,
+        mitre_limit=5.0,
+    ):
+        """Return a (Multi)LineString at a distance from the object.
+
+        The side, left or right, is determined by the sign of the `distance`
+        parameter (negative for right side offset, positive for left side
+        offset). The resolution of the buffer around each vertex of the object
+        increases by increasing the `quad_segs` keyword parameter.
+
+        The join style is for outside corners between line segments. Accepted
+        values are JOIN_STYLE.round (1), JOIN_STYLE.mitre (2), and
+        JOIN_STYLE.bevel (3).
+
+        The mitre ratio limit is used for very sharp corners. It is the ratio
+        of the distance from the corner to the end of the mitred offset corner.
+        When two line segments meet at a sharp angle, a miter join will extend
+        far beyond the original geometry. To prevent unreasonable geometry, the
+        mitre limit allows controlling the maximum length of the join corner.
+        Corners with a ratio which exceed the limit will be beveled.
+
+        Note: the behaviour regarding orientation of the resulting line
+        depends on the GEOS version. With GEOS < 3.11, the line retains the
+        same direction for a left offset (positive distance) or has reverse
+        direction for a right offset (negative distance), and this behaviour
+        was documented as such in previous Shapely versions. Starting with
+        GEOS 3.11, the function tries to preserve the orientation of the
+        original line.
+        """
+        if mitre_limit == 0.0:
+            raise ValueError("Cannot compute offset from zero-length line segment")
+        elif not np.isfinite(distance):
+            raise ValueError("offset_curve distance must be finite")
+        return shapely.offset_curve(
+            self,
+            distance,
+            quad_segs=quad_segs,
+            join_style=join_style,
+            mitre_limit=mitre_limit,
+        )
+
+    def parallel_offset(
+        self,
+        distance,
+        side="right",
+        resolution=16,
+        join_style=JOIN_STYLE.round,
+        mitre_limit=5.0,
+    ):
+        """Alternative method to :meth:`offset_curve` method.
+
+        Older alternative method to the :meth:`offset_curve` method, but uses
+        ``resolution`` instead of ``quad_segs`` and a ``side`` keyword
+        ('left' or 'right') instead of sign of the distance. This method is
+        kept for backwards compatibility for now, but is is recommended to
+        use :meth:`offset_curve` instead.
+        """
+        if side == "right":
+            distance *= -1
+        return self.offset_curve(
+            distance,
+            quad_segs=resolution,
+            join_style=join_style,
+            mitre_limit=mitre_limit,
+        )
+
+
+shapely.lib.registry[1] = LineString

.venv/lib/python3.10/site-packages/shapely/geometry/multilinestring.py

@@ -0,0 +1,95 @@
+"""Collections of linestrings and related utilities."""
+
+import shapely
+from shapely.errors import EmptyPartError
+from shapely.geometry import linestring
+from shapely.geometry.base import BaseMultipartGeometry
+
+__all__ = ["MultiLineString"]
+
+
+class MultiLineString(BaseMultipartGeometry):
+    """A collection of one or more LineStrings.
+
+    A MultiLineString has non-zero length and zero area.
+
+    Parameters
+    ----------
+    lines : sequence
+        A sequence LineStrings, or a sequence of line-like coordinate
+        sequences or array-likes (see accepted input for LineString).
+
+    Attributes
+    ----------
+    geoms : sequence
+        A sequence of LineStrings
+
+    Examples
+    --------
+    Construct a MultiLineString containing two LineStrings.
+
+    >>> from shapely import MultiLineString
+    >>> lines = MultiLineString([[[0, 0], [1, 2]], [[4, 4], [5, 6]]])
+
+    """
+
+    __slots__ = []
+
+    def __new__(self, lines=None):
+        """Create a new MultiLineString geometry."""
+        if not lines:
+            # allow creation of empty multilinestrings, to support unpickling
+            # TODO better empty constructor
+            return shapely.from_wkt("MULTILINESTRING EMPTY")
+        elif isinstance(lines, MultiLineString):
+            return lines
+
+        lines = getattr(lines, "geoms", lines)
+        subs = []
+        for item in lines:
+            line = linestring.LineString(item)
+            if line.is_empty:
+                raise EmptyPartError(
+                    "Can't create MultiLineString with empty component"
+                )
+            subs.append(line)
+
+        if len(lines) == 0:
+            return shapely.from_wkt("MULTILINESTRING EMPTY")
+
+        return shapely.multilinestrings(subs)
+
+    @property
+    def __geo_interface__(self):
+        """Return a GeoJSON-like mapping interface for this MultiLineString."""
+        return {
+            "type": "MultiLineString",
+            "coordinates": tuple(tuple(c for c in g.coords) for g in self.geoms),
+        }
+
+    def svg(self, scale_factor=1.0, stroke_color=None, opacity=None):
+        """Return a group of SVG polyline elements for the LineString geometry.
+
+        Parameters
+        ----------
+        scale_factor : float
+            Multiplication factor for the SVG stroke-width.  Default is 1.
+        stroke_color : str, optional
+            Hex string for stroke color. Default is to use "#66cc99" if
+            geometry is valid, and "#ff3333" if invalid.
+        opacity : float
+            Float number between 0 and 1 for color opacity. Default value is 0.8
+
+        """
+        if self.is_empty:
+            return "<g />"
+        if stroke_color is None:
+            stroke_color = "#66cc99" if self.is_valid else "#ff3333"
+        return (
+            "<g>"
+            + "".join(p.svg(scale_factor, stroke_color, opacity) for p in self.geoms)
+            + "</g>"
+        )
+
+
+shapely.lib.registry[5] = MultiLineString

.venv/lib/python3.10/site-packages/shapely/geometry/multipoint.py

@@ -0,0 +1,104 @@
+"""Collections of points and related utilities."""
+
+import numpy as np
+
+import shapely
+from shapely.errors import EmptyPartError
+from shapely.geometry import point
+from shapely.geometry.base import BaseMultipartGeometry
+
+__all__ = ["MultiPoint"]
+
+
+class MultiPoint(BaseMultipartGeometry):
+    """A collection of one or more Points.
+
+    A MultiPoint has zero area and zero length.
+
+    Parameters
+    ----------
+    points : sequence
+        A sequence of Points, or a sequence of (x, y [,z]) numeric coordinate
+        pairs or triples, or an array-like of shape (N, 2) or (N, 3).
+
+    Attributes
+    ----------
+    geoms : sequence
+        A sequence of Points
+
+    Examples
+    --------
+    Construct a MultiPoint containing two Points
+
+    >>> from shapely import MultiPoint, Point
+    >>> ob = MultiPoint([[0.0, 0.0], [1.0, 2.0]])
+    >>> len(ob.geoms)
+    2
+    >>> type(ob.geoms[0]) == Point
+    True
+
+    """
+
+    __slots__ = []
+
+    def __new__(self, points=None):
+        """Create a new MultiPoint geometry."""
+        if points is None:
+            # allow creation of empty multipoints, to support unpickling
+            # TODO better empty constructor
+            return shapely.from_wkt("MULTIPOINT EMPTY")
+        elif isinstance(points, MultiPoint):
+            return points
+        elif len(points) == 0:
+            return shapely.from_wkt("MULTIPOINT EMPTY")
+
+        if isinstance(points, np.ndarray) and np.issubdtype(points.dtype, np.number):
+            subs = shapely.points(points)
+            if not subs.ndim == 1:
+                raise ValueError("Invalid values passed to MultiPoint constructor")
+            if shapely.is_empty(subs).any():
+                raise EmptyPartError("Can't create MultiPoint with empty component")
+        else:
+            subs = []
+            for item in points:
+                p = point.Point(item)
+                if p.is_empty:
+                    raise EmptyPartError("Can't create MultiPoint with empty component")
+                subs.append(p)
+
+        return shapely.multipoints(subs)
+
+    @property
+    def __geo_interface__(self):
+        """Return a GeoJSON-like mapping interface for this MultiPoint."""
+        return {
+            "type": "MultiPoint",
+            "coordinates": tuple(g.coords[0] for g in self.geoms),
+        }
+
+    def svg(self, scale_factor=1.0, fill_color=None, opacity=None):
+        """Return a group of SVG circle elements for the MultiPoint geometry.
+
+        Parameters
+        ----------
+        scale_factor : float
+            Multiplication factor for the SVG circle diameters.  Default is 1.
+        fill_color : str, optional
+            Hex string for fill color. Default is to use "#66cc99" if
+            geometry is valid, and "#ff3333" if invalid.
+        opacity : float
+            Float number between 0 and 1 for color opacity. Default value is 0.6
+
+        """
+        if self.is_empty:
+            return "<g />"
+        if fill_color is None:
+            fill_color = "#66cc99" if self.is_valid else "#ff3333"
+        return (
+            "<g>"
+            + "".join(p.svg(scale_factor, fill_color, opacity) for p in self.geoms)
+            + "</g>"
+        )
+
+
+shapely.lib.registry[4] = MultiPoint

.venv/lib/python3.10/site-packages/shapely/geometry/multipolygon.py

@@ -0,0 +1,125 @@
+"""Collections of polygons and related utilities."""
+
+import shapely
+from shapely.geometry import polygon
+from shapely.geometry.base import BaseMultipartGeometry
+
+__all__ = ["MultiPolygon"]
+
+
+class MultiPolygon(BaseMultipartGeometry):
+    """A collection of one or more Polygons.
+
+    If component polygons overlap the collection is invalid and some
+    operations on it may fail.
+
+    Parameters
+    ----------
+    polygons : sequence
+        A sequence of Polygons, or a sequence of (shell, holes) tuples
+        where shell is the sequence representation of a linear ring
+        (see LinearRing) and holes is a sequence of such linear rings.
+
+    Attributes
+    ----------
+    geoms : sequence
+        A sequence of `Polygon` instances
+
+    Examples
+    --------
+    Construct a MultiPolygon from a sequence of coordinate tuples
+
+    >>> from shapely import MultiPolygon, Polygon
+    >>> ob = MultiPolygon([
+    ...     (
+    ...     ((0.0, 0.0), (0.0, 1.0), (1.0, 1.0), (1.0, 0.0)),
+    ...     [((0.1,0.1), (0.1,0.2), (0.2,0.2), (0.2,0.1))]
+    ...     )
+    ... ])
+    >>> len(ob.geoms)
+    1
+    >>> type(ob.geoms[0]) == Polygon
+    True
+
+    """
+
+    __slots__ = []
+
+    def __new__(self, polygons=None):
+        """Create a new MultiPolygon geometry."""
+        if polygons is None:
+            # allow creation of empty multipolygons, to support unpickling
+            # TODO better empty constructor
+            return shapely.from_wkt("MULTIPOLYGON EMPTY")
+        elif isinstance(polygons, MultiPolygon):
+            return polygons
+
+        polygons = getattr(polygons, "geoms", polygons)
+        # remove None and empty polygons from list of Polygons
+        polygons = [p for p in polygons if p]
+
+        L = len(polygons)
+
+        # Bail immediately if we have no input points.
+        if L == 0:
+            return shapely.from_wkt("MULTIPOLYGON EMPTY")
+
+        # This function does not accept sequences of MultiPolygons: there is
+        # no implicit flattening.
+        if any(isinstance(p, MultiPolygon) for p in polygons):
+            raise ValueError("Sequences of multi-polygons are not valid arguments")
+
+        subs = []
+        for i in range(L):
+            ob = polygons[i]
+            if not isinstance(ob, polygon.Polygon):
+                shell = ob[0]
+                if len(ob) > 1:
+                    holes = ob[1]
+                else:
+                    holes = None
+                p = polygon.Polygon(shell, holes)
+            else:
+                p = polygon.Polygon(ob)
+            subs.append(p)
+
+        return shapely.multipolygons(subs)
+
+    @property
+    def __geo_interface__(self):
+        """Return a GeoJSON-like mapping of the MultiPolygon geometry."""
+        allcoords = []
+        for geom in self.geoms:
+            coords = []
+            coords.append(tuple(geom.exterior.coords))
+            for hole in geom.interiors:
+                coords.append(tuple(hole.coords))
+            allcoords.append(tuple(coords))
+        return {"type": "MultiPolygon", "coordinates": allcoords}
+
+    def svg(self, scale_factor=1.0, fill_color=None, opacity=None):
+        """Return group of SVG path elements for the MultiPolygon geometry.
+
+        Parameters
+        ----------
+        scale_factor : float
+            Multiplication factor for the SVG stroke-width.  Default is 1.
+        fill_color : str, optional
+            Hex string for fill color. Default is to use "#66cc99" if
+            geometry is valid, and "#ff3333" if invalid.
+        opacity : float
+            Float number between 0 and 1 for color opacity. Default value is 0.6
+
+        """
+        if self.is_empty:
+            return "<g />"
+        if fill_color is None:
+            fill_color = "#66cc99" if self.is_valid else "#ff3333"
+        return (
+            "<g>"
+            + "".join(p.svg(scale_factor, fill_color, opacity) for p in self.geoms)
+            + "</g>"
+        )
+
+
+shapely.lib.registry[6] = MultiPolygon

.venv/lib/python3.10/site-packages/shapely/geometry/point.py

@@ -0,0 +1,166 @@
+"""Points and related utilities."""
+
+import numpy as np
+
+import shapely
+from shapely.errors import DimensionError
+from shapely.geometry.base import BaseGeometry
+
+__all__ = ["Point"]
+
+
+class Point(BaseGeometry):
+    """A geometry type that represents a single coordinate.
+
+    Each coordinate has x, y and possibly z and/or m values.
+
+    A point is a zero-dimensional feature and has zero length and zero area.
+
+    Parameters
+    ----------
+    args : float, or sequence of floats
+        The coordinates can either be passed as a single parameter, or as
+        individual float values using multiple parameters:
+
+        1) 1 parameter: a sequence or array-like of with 2 or 3 values.
+        2) 2 or 3 parameters (float): x, y, and possibly z.
+
+    Attributes
+    ----------
+    x, y, z, m : float
+        Coordinate values
+
+    Examples
+    --------
+    Constructing the Point using separate parameters for x and y:
+
+    >>> from shapely import Point
+    >>> p = Point(1.0, -1.0)
+
+    Constructing the Point using a list of x, y coordinates:
+
+    >>> p = Point([1.0, -1.0])
+    >>> print(p)
+    POINT (1 -1)
+    >>> p.y
+    -1.0
+    >>> p.x
+    1.0
+
+    """
+
+    __slots__ = []
+
+    def __new__(self, *args):
+        """Create a new Point geometry."""
+        if len(args) == 0:
+            # empty geometry
+            # TODO better constructor
+            return shapely.from_wkt("POINT EMPTY")
+        elif len(args) > 3:
+            raise TypeError(f"Point() takes at most 3 arguments ({len(args)} given)")
+        elif len(args) == 1:
+            coords = args[0]
+            if isinstance(coords, Point):
+                return coords
+
+            # Accept either (x, y) or [(x, y)]
+            if not hasattr(coords, "__getitem__"):  # generators
+                coords = list(coords)
+            coords = np.asarray(coords).squeeze()
+        else:
+            # 2 or 3 args
+            coords = np.array(args).squeeze()
+
+        if coords.ndim > 1:
+            raise ValueError(
+                f"Point() takes only scalar or 1-size vector arguments, got {args}"
+            )
+        if not np.issubdtype(coords.dtype, np.number):
+            coords = [float(c) for c in coords]
+        geom = shapely.points(coords)
+        if not isinstance(geom, Point):
+            raise ValueError("Invalid values passed to Point constructor")
+        return geom
+
+    # Coordinate getters and setters
+
+    @property
+    def x(self):
+        """Return x coordinate."""
+        return float(shapely.get_x(self))
+
+    @property
+    def y(self):
+        """Return y coordinate."""
+        return float(shapely.get_y(self))
+
+    @property
+    def z(self):
+        """Return z coordinate."""
+        z = shapely.get_z(self)
+        if np.isnan(z) and not shapely.has_z(self):
+            raise DimensionError("This point has no z coordinate.")
+        return float(z)
+
+    @property
+    def m(self):
+        """Return m coordinate.
+
+        .. versionadded:: 2.1.0
+           Also requires GEOS 3.12.0 or later.
+        """
+        if not shapely.has_m(self):
+            raise DimensionError("This point has no m coordinate.")
+        return float(shapely.get_m(self))
+
+    @property
+    def __geo_interface__(self):
+        """Return a GeoJSON-like mapping of the Point geometry."""
+        coords = self.coords
+        return {"type": "Point", "coordinates": coords[0] if len(coords) > 0 else ()}
+
+    def svg(self, scale_factor=1.0, fill_color=None, opacity=None):
+        """Return SVG circle element for the Point geometry.
+
+        Parameters
+        ----------
+        scale_factor : float
+            Multiplication factor for the SVG circle diameter.  Default is 1.
+        fill_color : str, optional
+            Hex string for fill color. Default is to use "#66cc99" if
+            geometry is valid, and "#ff3333" if invalid.
+        opacity : float
+            Float number between 0 and 1 for color opacity. Default value is 0.6
+
+        """
+        if self.is_empty:
+            return "<g />"
+        if fill_color is None:
+            fill_color = "#66cc99" if self.is_valid else "#ff3333"
+        if opacity is None:
+            opacity = 0.6
+        return (
+            f'<circle cx="{self.x}" cy="{self.y}" r="{3.0 * scale_factor}" '
+            f'stroke="#555555" stroke-width="{1.0 * scale_factor}" fill="{fill_color}" '
+            f'opacity="{opacity}" />'
+        )
+
+    @property
+    def xy(self):
+        """Separate arrays of X and Y coordinate values.
+
+        Examples
+        --------
+        >>> from shapely import Point
+        >>> x, y = Point(0, 0).xy
+        >>> list(x)
+        [0.0]
+        >>> list(y)
+        [0.0]
+
+        """
+        return self.coords.xy
+
+
+shapely.lib.registry[0] = Point

.venv/lib/python3.10/site-packages/shapely/geometry/polygon.py

@@ -0,0 +1,345 @@
+"""Polygons and their linear ring components."""
+
+import numpy as np
+
+import shapely
+from shapely import _geometry_helpers
+from shapely.algorithms.cga import signed_area  # noqa
+from shapely.errors import TopologicalError
+from shapely.geometry.base import BaseGeometry
+from shapely.geometry.linestring import LineString
+from shapely.geometry.point import Point
+
+__all__ = ["LinearRing", "Polygon", "orient"]
+
+
+def _unpickle_linearring(wkb):
+    linestring = shapely.from_wkb(wkb)
+    srid = shapely.get_srid(linestring)
+    linearring = _geometry_helpers.linestring_to_linearring(linestring)
+    if srid:
+        linearring = shapely.set_srid(linearring, srid)
+    return linearring
+
+
+class LinearRing(LineString):
+    """Geometry type composed of one or more line segments that forms a closed loop.
+
+    A LinearRing is a closed, one-dimensional feature.
+    A LinearRing that crosses itself or touches itself at a single point is
+    invalid and operations on it may fail.
+
+    Parameters
+    ----------
+    coordinates : sequence
+        A sequence of (x, y [,z]) numeric coordinate pairs or triples, or
+        an array-like with shape (N, 2) or (N, 3).
+        Also can be a sequence of Point objects.
+
+    Notes
+    -----
+    Rings are automatically closed. There is no need to specify a final
+    coordinate pair identical to the first.
+
+    Examples
+    --------
+    Construct a square ring.
+
+    >>> from shapely import LinearRing
+    >>> ring = LinearRing( ((0, 0), (0, 1), (1 ,1 ), (1 , 0)) )
+    >>> ring.is_closed
+    True
+    >>> list(ring.coords)
+    [(0.0, 0.0), (0.0, 1.0), (1.0, 1.0), (1.0, 0.0), (0.0, 0.0)]
+    >>> ring.length
+    4.0
+
+    """
+
+    __slots__ = []
+
+    def __new__(self, coordinates=None):
+        """Create a new LinearRing geometry."""
+        if coordinates is None:
+            # empty geometry
+            # TODO better way?
+            return shapely.from_wkt("LINEARRING EMPTY")
+        elif isinstance(coordinates, LineString):
+            if type(coordinates) is LinearRing:
+                # return original objects since geometries are immutable
+                return coordinates
+            elif not coordinates.is_valid:
+                raise TopologicalError("An input LineString must be valid.")
+            else:
+                # LineString
+                # TODO convert LineString to LinearRing more directly?
+                coordinates = coordinates.coords
+
+        else:
+            if hasattr(coordinates, "__array__"):
+                coordinates = np.asarray(coordinates)
+            if isinstance(coordinates, np.ndarray) and np.issubdtype(
+                coordinates.dtype, np.number
+            ):
+                pass
+            else:
+                # check coordinates on points
+                def _coords(o):
+                    if isinstance(o, Point):
+                        return o.coords[0]
+                    else:
+                        return [float(c) for c in o]
+
+                coordinates = np.array([_coords(o) for o in coordinates])
+                if not np.issubdtype(coordinates.dtype, np.number):
+                    # conversion of coords to 2D array failed, this might be due
+                    # to inconsistent coordinate dimensionality
+                    raise ValueError("Inconsistent coordinate dimensionality")
+
+        if len(coordinates) == 0:
+            # empty geometry
+            # TODO better constructor + should shapely.linearrings handle this?
+            return shapely.from_wkt("LINEARRING EMPTY")
+
+        geom = shapely.linearrings(coordinates)
+        if not isinstance(geom, LinearRing):
+            raise ValueError("Invalid values passed to LinearRing constructor")
+        return geom
+
+    @property
+    def __geo_interface__(self):
+        """Return a GeoJSON-like mapping of the LinearRing geometry."""
+        return {"type": "LinearRing", "coordinates": tuple(self.coords)}
+
+    def __reduce__(self):
+        """Pickle support.
+
+        WKB doesn't differentiate between LineString and LinearRing so we
+        need to move the coordinate sequence into the correct geometry type
+        """
+        return (_unpickle_linearring, (shapely.to_wkb(self, include_srid=True),))
+
+    @property
+    def is_ccw(self):
+        """True if the ring is oriented counter clock-wise."""
+        return bool(shapely.is_ccw(self))
+
+    @property
+    def is_simple(self):
+        """True if the geometry is simple.
+
+        Simple means that any self-intersections are only at boundary points.
+        """
+        return bool(shapely.is_simple(self))
+
+
+shapely.lib.registry[2] = LinearRing
+
+
+class InteriorRingSequence:
+    _parent = None
+    _ndim = None
+    _index = 0
+    _length = 0
+
+    def __init__(self, parent):
+        self._parent = parent
+        self._ndim = parent._ndim
+
+    def __iter__(self):
+        self._index = 0
+        self._length = self.__len__()
+        return self
+
+    def __next__(self):
+        if self._index < self._length:
+            ring = self._get_ring(self._index)
+            self._index += 1
+            return ring
+        else:
+            raise StopIteration
+
+    def __len__(self):
+        return shapely.get_num_interior_rings(self._parent)
+
+    def __getitem__(self, key):
+        m = self.__len__()
+        if isinstance(key, int):
+            if key + m < 0 or key >= m:
+                raise IndexError("index out of range")
+            if key < 0:
+                i = m + key
+            else:
+                i = key
+            return self._get_ring(i)
+        elif isinstance(key, slice):
+            res = []
+            start, stop, stride = key.indices(m)
+            for i in range(start, stop, stride):
+                res.append(self._get_ring(i))
+            return res
+        else:
+            raise TypeError("key must be an index or slice")
+
+    def _get_ring(self, i):
+        return shapely.get_interior_ring(self._parent, i)
+
+
+class Polygon(BaseGeometry):
+    """A geometry type representing an area that is enclosed by a linear ring.
+
+    A polygon is a two-dimensional feature and has a non-zero area. It may
+    have one or more negative-space "holes" which are also bounded by linear
+    rings. If any rings cross each other, the feature is invalid and
+    operations on it may fail.
+
+    Parameters
+    ----------
+    shell : sequence
+        A sequence of (x, y [,z]) numeric coordinate pairs or triples, or
+        an array-like with shape (N, 2) or (N, 3).
+        Also can be a sequence of Point objects.
+    holes : sequence
+        A sequence of objects which satisfy the same requirements as the
+        shell parameters above
+
+    Attributes
+    ----------
+    exterior : LinearRing
+        The ring which bounds the positive space of the polygon.
+    interiors : sequence
+        A sequence of rings which bound all existing holes.
+
+    Examples
+    --------
+    Create a square polygon with no holes
+
+    >>> from shapely import Polygon
+    >>> coords = ((0., 0.), (0., 1.), (1., 1.), (1., 0.), (0., 0.))
+    >>> polygon = Polygon(coords)
+    >>> polygon.area
+    1.0
+
+    """
+
+    __slots__ = []
+
+    def __new__(self, shell=None, holes=None):
+        """Create a new Polygon geometry."""
+        if shell is None:
+            # empty geometry
+            # TODO better way?
+            return shapely.from_wkt("POLYGON EMPTY")
+        elif isinstance(shell, Polygon):
+            # return original objects since geometries are immutable
+            return shell
+        else:
+            shell = LinearRing(shell)
+
+        if holes is not None:
+            if len(holes) == 0:
+                # shapely constructor cannot handle holes=[]
+                holes = None
+            else:
+                holes = [LinearRing(ring) for ring in holes]
+
+        geom = shapely.polygons(shell, holes=holes)
+        if not isinstance(geom, Polygon):
+            raise ValueError("Invalid values passed to Polygon constructor")
+        return geom
+
+    @property
+    def exterior(self):
+        """Return the exterior ring of the polygon."""
+        return shapely.get_exterior_ring(self)
+
+    @property
+    def interiors(self):
+        """Return the sequence of interior rings of the polygon."""
+        if self.is_empty:
+            return []
+        return InteriorRingSequence(self)
+
+    @property
+    def coords(self):
+        """Not implemented for polygons."""
+        raise NotImplementedError(
+            "Component rings have coordinate sequences, but the polygon does not"
+        )
+
+    @property
+    def __geo_interface__(self):
+        """Return a GeoJSON-like mapping of the Polygon geometry."""
+        if self.exterior == LinearRing():
+            coords = []
+        else:
+            coords = [tuple(self.exterior.coords)]
+            for hole in self.interiors:
+                coords.append(tuple(hole.coords))
+        return {"type": "Polygon", "coordinates": tuple(coords)}
+
+    def svg(self, scale_factor=1.0, fill_color=None, opacity=None):
+        """Return SVG path element for the Polygon geometry.
+
+        Parameters
+        ----------
+        scale_factor : float
+            Multiplication factor for the SVG stroke-width.  Default is 1.
+        fill_color : str, optional
+            Hex string for fill color. Default is to use "#66cc99" if
+            geometry is valid, and "#ff3333" if invalid.
+        opacity : float
+            Float number between 0 and 1 for color opacity. Default value is 0.6
+
+        """
+        if self.is_empty:
+            return "<g />"
+        if fill_color is None:
+            fill_color = "#66cc99" if self.is_valid else "#ff3333"
+        if opacity is None:
+            opacity = 0.6
+        exterior_coords = [["{},{}".format(*c) for c in self.exterior.coords]]
+        interior_coords = [
+            ["{},{}".format(*c) for c in interior.coords] for interior in self.interiors
+        ]
+        path = " ".join(
+            [
+                "M {} L {} z".format(coords[0], " L ".join(coords[1:]))
+                for coords in exterior_coords + interior_coords
+            ]
+        )
+        return (
+            f'<path fill-rule="evenodd" fill="{fill_color}" stroke="#555555" '
+            f'stroke-width="{2.0 * scale_factor}" opacity="{opacity}" d="{path}" />'
+        )
+
+    @classmethod
+    def from_bounds(cls, xmin, ymin, xmax, ymax):
+        """Construct a `Polygon()` from spatial bounds."""
+        return cls([(xmin, ymin), (xmin, ymax), (xmax, ymax), (xmax, ymin)])
+
+
+shapely.lib.registry[3] = Polygon
+
+
+def orient(polygon, sign=1.0):
+    """Return an oriented polygon.
+
+    It is recommended to use :func:`shapely.orient_polygons` instead.
+
+    Parameters
+    ----------
+    polygon : shapely.Polygon
+    sign : float, default 1.
+        The sign of the result's signed area.
+        A non-negative sign means that the coordinates of the geometry's exterior
+        rings will be oriented counter-clockwise.
+
+    Returns
+    -------
+    Geometry or array_like
+
+    Refer to :func:`shapely.orient_polygons` for full documentation.
+
+    """
+    return shapely.orient_polygons(polygon, exterior_cw=sign < 0.0)

.venv/lib/python3.10/site-packages/shapely/geos.py

@@ -0,0 +1,17 @@
+"""Proxies for libgeos, GEOS-specific exceptions, and utilities."""
+
+import warnings
+
+import shapely
+
+warnings.warn(
+    "The 'shapely.geos' module is deprecated, and will be removed in a future version. "
+    "All attributes of 'shapely.geos' are available directly from the top-level "
+    "'shapely' namespace (since shapely 2.0.0).",
+    DeprecationWarning,
+    stacklevel=2,
+)
+
+geos_version_string = shapely.geos_capi_version_string
+geos_version = shapely.geos_version
+geos_capi_version = shapely.geos_capi_version

.venv/lib/python3.10/site-packages/shapely/io.py

@@ -0,0 +1,407 @@
+"""Input/output functions for Shapely geometries."""
+
+import numpy as np
+
+from shapely import geos_version, lib
+from shapely._enum import ParamEnum
+
+# include ragged array functions here for reference documentation purpose
+from shapely._ragged_array import from_ragged_array, to_ragged_array
+from shapely.decorators import requires_geos
+from shapely.errors import UnsupportedGEOSVersionError
+
+__all__ = [
+    "from_geojson",
+    "from_ragged_array",
+    "from_wkb",
+    "from_wkt",
+    "to_geojson",
+    "to_ragged_array",
+    "to_wkb",
+    "to_wkt",
+]
+
+
+# Allowed options for handling WKB/WKT decoding errors
+# Note: cannot use standard constructor since "raise" is a keyword
+DecodingErrorOptions = ParamEnum(
+    "DecodingErrorOptions", {"ignore": 0, "warn": 1, "raise": 2, "fix": 3}
+)
+
+WKBFlavorOptions = ParamEnum("WKBFlavorOptions", {"extended": 1, "iso": 2})
+
+
+def to_wkt(
+    geometry,
+    rounding_precision=6,
+    trim=True,
+    output_dimension=None,
+    old_3d=False,
+    **kwargs,
+):
+    """Convert to the Well-Known Text (WKT) representation of a Geometry.
+
+    The Well-known Text format is defined in the `OGC Simple Features
+    Specification for SQL <https://www.opengeospatial.org/standards/sfs>`__.
+
+    The following limitations apply to WKT serialization:
+
+    - only simple empty geometries can be 3D, empty collections are always 2D
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries to convert to WKT.
+    rounding_precision : int, default 6
+        The rounding precision when writing the WKT string. Set to a value of
+        -1 to indicate the full precision.
+    trim : bool, default True
+        If True, trim unnecessary decimals (trailing zeros).
+    output_dimension : int, default None
+        The output dimension for the WKT string. Supported values are 2, 3 and
+        4 for GEOS 3.12+. Default None will automatically choose 3 or 4,
+        depending on the version of GEOS.
+        Specifying 3 means that up to 3 dimensions will be written but 2D
+        geometries will still be represented as 2D in the WKT string.
+    old_3d : bool, default False
+        Enable old style 3D/4D WKT generation. By default, new style 3D/4D WKT
+        (ie. "POINT Z (10 20 30)") is returned, but with ``old_3d=True``
+        the WKT will be formatted in the style "POINT (10 20 30)".
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Point
+    >>> shapely.to_wkt(Point(0, 0))
+    'POINT (0 0)'
+    >>> shapely.to_wkt(Point(0, 0), rounding_precision=3, trim=False)
+    'POINT (0.000 0.000)'
+    >>> shapely.to_wkt(Point(0, 0), rounding_precision=-1, trim=False)
+    'POINT (0.0000000000000000 0.0000000000000000)'
+    >>> shapely.to_wkt(Point(1, 2, 3), trim=True)
+    'POINT Z (1 2 3)'
+    >>> shapely.to_wkt(Point(1, 2, 3), trim=True, output_dimension=2)
+    'POINT (1 2)'
+    >>> shapely.to_wkt(Point(1, 2, 3), trim=True, old_3d=True)
+    'POINT (1 2 3)'
+
+    Notes
+    -----
+    The defaults differ from the default of some GEOS versions. To mimic this for
+    versions before GEOS 3.12, use::
+
+        shapely.to_wkt(geometry, rounding_precision=-1, trim=False, output_dimension=2)
+
+    """
+    if not np.isscalar(rounding_precision):
+        raise TypeError("rounding_precision only accepts scalar values")
+    if not np.isscalar(trim):
+        raise TypeError("trim only accepts scalar values")
+    if output_dimension is None:
+        output_dimension = 3 if geos_version < (3, 12, 0) else 4
+    elif not np.isscalar(output_dimension):
+        raise TypeError("output_dimension only accepts scalar values")
+    if not np.isscalar(old_3d):
+        raise TypeError("old_3d only accepts scalar values")
+
+    return lib.to_wkt(
+        geometry,
+        np.intc(rounding_precision),
+        np.bool_(trim),
+        np.intc(output_dimension),
+        np.bool_(old_3d),
+        **kwargs,
+    )
+
+
+def to_wkb(
+    geometry,
+    hex=False,
+    output_dimension=None,
+    byte_order=-1,
+    include_srid=False,
+    flavor="extended",
+    **kwargs,
+):
+    r"""Convert to the Well-Known Binary (WKB) representation of a Geometry.
+
+    The Well-Known Binary format is defined in the `OGC Simple Features
+    Specification for SQL <https://www.opengeospatial.org/standards/sfs>`__.
+
+    The following limitations apply to WKB serialization:
+
+    - linearrings will be converted to linestrings
+    - a point with only NaN coordinates is converted to an empty point
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries to convert to WKB.
+    hex : bool, default False
+        If true, export the WKB as a hexadecimal string. The default is to
+        return a binary bytes object.
+    output_dimension : int, default None
+        The output dimension for the WKB. Supported values are 2, 3 and 4 for
+        GEOS 3.12+. Default None will automatically choose 3 or 4, depending on
+        the version of GEOS.
+        Specifying 3 means that up to 3 dimensions will be written but 2D
+        geometries will still be represented as 2D in the WKB representation.
+    byte_order : int, default -1
+        Defaults to native machine byte order (-1). Use 0 to force big endian
+        and 1 for little endian.
+    include_srid : bool, default False
+        If True, the SRID is be included in WKB (this is an extension
+        to the OGC WKB specification). Not allowed when flavor is "iso".
+    flavor : {"iso", "extended"}, default "extended"
+        Which flavor of WKB will be returned. The flavor determines how
+        extra dimensionality is encoded with the type number, and whether
+        SRID can be included in the WKB. ISO flavor is "more standard" for
+        3D output, and does not support SRID embedding.
+        Both flavors are equivalent when ``output_dimension=2`` (or with 2D
+        geometries) and ``include_srid=False``.
+        The `from_wkb` function can read both flavors.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Point
+    >>> point = Point(1, 1)
+    >>> shapely.to_wkb(point, byte_order=1)
+    b'\x01\x01\x00\x00\x00\x00\x00\x00\x00\x00\x00\xf0?\x00\x00\x00\x00\x00\x00\xf0?'
+    >>> shapely.to_wkb(point, hex=True, byte_order=1)
+    '0101000000000000000000F03F000000000000F03F'
+
+    """
+    if not np.isscalar(hex):
+        raise TypeError("hex only accepts scalar values")
+    if output_dimension is None:
+        output_dimension = 3 if geos_version < (3, 12, 0) else 4
+    elif not np.isscalar(output_dimension):
+        raise TypeError("output_dimension only accepts scalar values")
+    if not np.isscalar(byte_order):
+        raise TypeError("byte_order only accepts scalar values")
+    if not np.isscalar(include_srid):
+        raise TypeError("include_srid only accepts scalar values")
+    if not np.isscalar(flavor):
+        raise TypeError("flavor only accepts scalar values")
+    if lib.geos_version < (3, 10, 0) and flavor == "iso":
+        raise UnsupportedGEOSVersionError(
+            'The "iso" option requires at least GEOS 3.10.0'
+        )
+    if flavor == "iso" and include_srid:
+        raise ValueError('flavor="iso" and include_srid=True cannot be used together')
+    flavor = WKBFlavorOptions.get_value(flavor)
+
+    return lib.to_wkb(
+        geometry,
+        np.bool_(hex),
+        np.intc(output_dimension),
+        np.intc(byte_order),
+        np.bool_(include_srid),
+        np.intc(flavor),
+        **kwargs,
+    )
+
+
+@requires_geos("3.10.0")
+def to_geojson(geometry, indent=None, **kwargs):
+    """Convert to the GeoJSON representation of a Geometry.
+
+    The GeoJSON format is defined in the `RFC 7946 <https://geojson.org/>`__.
+    NaN (not-a-number) coordinates will be written as 'null'.
+
+    The following are currently unsupported:
+
+    - Geometries of type LINEARRING: these are output as 'null'.
+    - Three-dimensional geometries: the third dimension is ignored.
+
+    Parameters
+    ----------
+    geometry : str, bytes or array_like
+        Geometry or geometries to convert to GeoJSON.
+    indent : int, optional
+        If indent is a non-negative integer, then GeoJSON will be formatted.
+        An indent level of 0 will only insert newlines. None (the default)
+        selects the most compact representation.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Point
+    >>> point = Point(1, 1)
+    >>> shapely.to_geojson(point)
+    '{"type":"Point","coordinates":[1.0,1.0]}'
+    >>> print(shapely.to_geojson(point, indent=2))
+    {
+      "type": "Point",
+      "coordinates": [
+          1.0,
+          1.0
+      ]
+    }
+
+    """
+    # GEOS Tickets:
+    # - handle linearrings: https://trac.osgeo.org/geos/ticket/1140
+    # - support 3D: https://trac.osgeo.org/geos/ticket/1141
+    if indent is None:
+        indent = -1
+    elif not np.isscalar(indent):
+        raise TypeError("indent only accepts scalar values")
+    elif indent < 0:
+        raise ValueError("indent cannot be negative")
+
+    return lib.to_geojson(geometry, np.intc(indent), **kwargs)
+
+
+def from_wkt(geometry, on_invalid="raise", **kwargs):
+    """Create geometries from the Well-Known Text (WKT) representation.
+
+    The Well-known Text format is defined in the `OGC Simple Features
+    Specification for SQL <https://www.opengeospatial.org/standards/sfs>`__.
+
+    Parameters
+    ----------
+    geometry : str or array_like
+        The WKT string(s) to convert.
+    on_invalid : {"raise", "warn", "ignore", "fix"}, default "raise"
+        Indicates what to do when an invalid WKT string is encountered. Note
+        that the validations involved are very basic, e.g. the minimum number of
+        points for the geometry type. For a thorough check, use
+        :func:`is_valid` after conversion to geometries. Valid options are:
+
+        - raise: an exception will be raised if any input geometry is invalid.
+        - warn: a warning will be raised and invalid WKT geometries will be
+          returned as ``None``.
+        - ignore: invalid geometries will be returned as ``None`` without a
+          warning.
+        - fix: an effort is made to fix invalid input geometries (currently just
+          unclosed rings). If this is not possible, they are returned as
+          ``None`` without a warning. Requires GEOS >= 3.11.
+
+          .. versionadded:: 2.1.0
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> shapely.from_wkt('POINT (0 0)')
+    <POINT (0 0)>
+
+    """
+    if not np.isscalar(on_invalid):
+        raise TypeError("on_invalid only accepts scalar values")
+
+    invalid_handler = np.uint8(DecodingErrorOptions.get_value(on_invalid))
+
+    return lib.from_wkt(geometry, invalid_handler, **kwargs)
+
+
+def from_wkb(geometry, on_invalid="raise", **kwargs):
+    r"""Create geometries from the Well-Known Binary (WKB) representation.
+
+    The Well-Known Binary format is defined in the `OGC Simple Features
+    Specification for SQL <https://www.opengeospatial.org/standards/sfs>`__.
+
+    Parameters
+    ----------
+    geometry : str or array_like
+        The WKB byte object(s) to convert.
+    on_invalid : {"raise", "warn", "ignore", "fix"}, default "raise"
+        Indicates what to do when an invalid WKB is encountered. Note that the
+        validations involved are very basic, e.g. the minimum number of points
+        for the geometry type. For a thorough check, use :func:`is_valid` after
+        conversion to geometries. Valid options are:
+
+        - raise: an exception will be raised if any input geometry is invalid.
+        - warn: a warning will be raised and invalid WKT geometries will be
+          returned as ``None``.
+        - ignore: invalid geometries will be returned as ``None`` without a
+          warning.
+        - fix: an effort is made to fix invalid input geometries (currently just
+          unclosed rings). If this is not possible, they are returned as
+          ``None`` without a warning. Requires GEOS >= 3.11.
+
+          .. versionadded:: 2.1.0
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> shapely.from_wkb(b'\x01\x01\x00\x00\x00\x00\x00\x00\x00\x00\x00\xf0?\x00\x00\x00\x00\x00\x00\xf0?')
+    <POINT (1 1)>
+
+    """  # noqa: E501
+    if not np.isscalar(on_invalid):
+        raise TypeError("on_invalid only accepts scalar values")
+
+    invalid_handler = np.uint8(DecodingErrorOptions.get_value(on_invalid))
+
+    # ensure the input has object dtype, to avoid numpy inferring it as a
+    # fixed-length string dtype (which removes trailing null bytes upon access
+    # of array elements)
+    geometry = np.asarray(geometry, dtype=object)
+    return lib.from_wkb(geometry, invalid_handler, **kwargs)
+
+
+@requires_geos("3.10.1")
+def from_geojson(geometry, on_invalid="raise", **kwargs):
+    """Create geometries from GeoJSON representations (strings).
+
+    If a GeoJSON is a FeatureCollection, it is read as a single geometry
+    (with type GEOMETRYCOLLECTION). This may be unpacked using
+    :meth:`shapely.get_parts`. Properties are not read.
+
+    The GeoJSON format is defined in `RFC 7946 <https://geojson.org/>`__.
+
+    The following are currently unsupported:
+
+    - Three-dimensional geometries: the third dimension is ignored.
+    - Geometries having 'null' in the coordinates.
+
+    Parameters
+    ----------
+    geometry : str, bytes or array_like
+        The GeoJSON string or byte object(s) to convert.
+    on_invalid : {"raise", "warn", "ignore"}, default "raise"
+        - raise: an exception will be raised if an input GeoJSON is invalid.
+        - warn: a warning will be raised and invalid input geometries will be
+          returned as ``None``.
+        - ignore: invalid input geometries will be returned as ``None`` without
+          a warning.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    See Also
+    --------
+    get_parts
+
+    Examples
+    --------
+    >>> import shapely
+    >>> shapely.from_geojson('{"type": "Point","coordinates": [1, 2]}')
+    <POINT (1 2)>
+
+    """
+    # GEOS Tickets:
+    # - support 3D: https://trac.osgeo.org/geos/ticket/1141
+    # - handle null coordinates: https://trac.osgeo.org/geos/ticket/1142
+    if not np.isscalar(on_invalid):
+        raise TypeError("on_invalid only accepts scalar values")
+
+    invalid_handler = np.uint8(DecodingErrorOptions.get_value(on_invalid))
+
+    # ensure the input has object dtype, to avoid numpy inferring it as a
+    # fixed-length string dtype (which removes trailing null bytes upon access
+    # of array elements)
+    geometry = np.asarray(geometry, dtype=object)
+
+    return lib.from_geojson(geometry, invalid_handler, **kwargs)

.venv/lib/python3.10/site-packages/shapely/linear.py

@@ -0,0 +1,236 @@
+"""Linear geometry functions."""
+
+from shapely import lib
+from shapely.decorators import deprecate_positional, multithreading_enabled
+from shapely.errors import UnsupportedGEOSVersionError
+
+__all__ = [
+    "line_interpolate_point",
+    "line_locate_point",
+    "line_merge",
+    "shared_paths",
+    "shortest_line",
+]
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   line_interpolate_point(line, distance, normalized=False, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'normalized' arg
+#    same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'normalized'
+#   line_interpolate_point(line, distance, *, normalized=False, **kwargs)
+
+
+@deprecate_positional(["normalized"], category=DeprecationWarning)
+@multithreading_enabled
+def line_interpolate_point(line, distance, normalized=False, **kwargs):
+    """Return a point interpolated at given distance on a line.
+
+    Parameters
+    ----------
+    line : Geometry or array_like
+        For multilinestrings or geometrycollections, the first geometry is taken
+        and the rest is ignored. This function raises a TypeError for non-linear
+        geometries. For empty linear geometries, empty points are returned.
+    distance : float or array_like
+        Negative values measure distance from the end of the line. Out-of-range
+        values will be clipped to the line endings.
+    normalized : bool, default False
+        If True, the distance is a fraction of the total
+        line length instead of the absolute distance.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line = LineString([(0, 2), (0, 10)])
+    >>> shapely.line_interpolate_point(line, 2)
+    <POINT (0 4)>
+    >>> shapely.line_interpolate_point(line, 100)
+    <POINT (0 10)>
+    >>> shapely.line_interpolate_point(line, -2)
+    <POINT (0 8)>
+    >>> shapely.line_interpolate_point(line, [0.25, -0.25], normalized=True).tolist()
+    [<POINT (0 4)>, <POINT (0 8)>]
+    >>> shapely.line_interpolate_point(LineString(), 1)
+    <POINT EMPTY>
+
+    """
+    if normalized:
+        return lib.line_interpolate_point_normalized(line, distance)
+    else:
+        return lib.line_interpolate_point(line, distance)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   line_locate_point(line, other, normalized=False, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'normalized' arg
+#    same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'normalized'
+#   line_locate_point(line, other, *, normalized=False, **kwargs)
+
+
+@deprecate_positional(["normalized"], category=DeprecationWarning)
+@multithreading_enabled
+def line_locate_point(line, other, normalized=False, **kwargs):
+    """Return the distance to the line origin of given point.
+
+    If given point does not intersect with the line, the point will first be
+    projected onto the line after which the distance is taken.
+
+    Parameters
+    ----------
+    line : Geometry or array_like
+        Line or lines to calculate the distance to.
+    other : Geometry or array_like
+        Point or points to calculate the distance from.
+    normalized : bool, default False
+        If True, the distance is a fraction of the total line length instead of
+        the absolute distance.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point
+    >>> line = LineString([(0, 2), (0, 10)])
+    >>> point = Point(4, 4)
+    >>> shapely.line_locate_point(line, point)
+    2.0
+    >>> shapely.line_locate_point(line, point, normalized=True)
+    0.25
+    >>> shapely.line_locate_point(line, Point(0, 18))
+    8.0
+    >>> shapely.line_locate_point(LineString(), point)
+    nan
+
+    """
+    if normalized:
+        return lib.line_locate_point_normalized(line, other)
+    else:
+        return lib.line_locate_point(line, other)
+
+
+@multithreading_enabled
+def line_merge(line, directed=False, **kwargs):
+    """Return (Multi)LineStrings formed by combining the lines in a MultiLineString.
+
+    Lines are joined together at their endpoints in case two lines are
+    intersecting. Lines are not joined when 3 or more lines are intersecting at
+    the endpoints. Line elements that cannot be joined are kept as is in the
+    resulting MultiLineString.
+
+    The direction of each merged LineString will be that of the majority of the
+    LineStrings from which it was derived. Except if ``directed=True`` is
+    specified, then the operation will not change the order of points within
+    lines and so only lines which can be joined with no change in direction
+    are merged.
+
+    Parameters
+    ----------
+    line : Geometry or array_like
+        Linear geometry or geometries to merge.
+    directed : bool, default False
+        Only combine lines if possible without changing point order.
+        Requires GEOS >= 3.11.0
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import MultiLineString
+    >>> shapely.line_merge(MultiLineString([[(0, 2), (0, 10)], [(0, 10), (5, 10)]]))
+    <LINESTRING (0 2, 0 10, 5 10)>
+    >>> shapely.line_merge(MultiLineString([[(0, 2), (0, 10)], [(0, 11), (5, 10)]]))
+    <MULTILINESTRING ((0 2, 0 10), (0 11, 5 10))>
+    >>> shapely.line_merge(MultiLineString())
+    <GEOMETRYCOLLECTION EMPTY>
+    >>> shapely.line_merge(MultiLineString([[(0, 0), (1, 0)], [(0, 0), (3, 0)]]))
+    <LINESTRING (1 0, 0 0, 3 0)>
+    >>> shapely.line_merge(MultiLineString([[(0, 0), (1, 0)], [(0, 0), (3, 0)]]), \
+directed=True)
+    <MULTILINESTRING ((0 0, 1 0), (0 0, 3 0))>
+
+    """
+    if directed:
+        if lib.geos_version < (3, 11, 0):
+            raise UnsupportedGEOSVersionError(
+                "'{}' requires at least GEOS {}.{}.{}.".format(
+                    "line_merge", *(3, 11, 0)
+                )
+            )
+        return lib.line_merge_directed(line, **kwargs)
+    return lib.line_merge(line, **kwargs)
+
+
+@multithreading_enabled
+def shared_paths(a, b, **kwargs):
+    """Return the shared paths between a and b.
+
+    Both geometries should be linestrings or arrays of linestrings.
+    A geometrycollection or array of geometrycollections is returned
+    with two elements in each geometrycollection. The first element is a
+    multilinestring containing shared paths with the same direction
+    for both inputs. The second element is a multilinestring containing
+    shared paths with the opposite direction for the two inputs.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Linestring or linestrings to compare.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line1 = LineString([(0, 0), (1, 0), (1, 1), (0, 1), (0, 0)])
+    >>> line2 = LineString([(1, 0), (2, 0), (2, 1), (1, 1), (1, 0)])
+    >>> shapely.shared_paths(line1, line2).wkt
+    'GEOMETRYCOLLECTION (MULTILINESTRING EMPTY, MULTILINESTRING ((1 0, 1 1)))'
+    >>> line3 = LineString([(1, 1), (0, 1)])
+    >>> shapely.shared_paths(line1, line3).wkt
+    'GEOMETRYCOLLECTION (MULTILINESTRING ((1 1, 0 1)), MULTILINESTRING EMPTY)'
+
+    """
+    return lib.shared_paths(a, b, **kwargs)
+
+
+@multithreading_enabled
+def shortest_line(a, b, **kwargs):
+    """Return the shortest line between two geometries.
+
+    The resulting line consists of two points, representing the nearest
+    points between the geometry pair. The line always starts in the first
+    geometry `a` and ends in the second geometry `b`. The endpoints of the
+    line will not necessarily be existing vertices of the input geometries
+    `a` and `b`, but can also be a point along a line segment.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to compare.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    See Also
+    --------
+    prepare : improve performance by preparing ``a`` (the first argument)
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line1 = LineString([(0, 0), (1, 0), (1, 1), (0, 1), (0, 0)])
+    >>> line2 = LineString([(0, 3), (3, 0), (5, 3)])
+    >>> shapely.shortest_line(line1, line2)
+    <LINESTRING (1 1, 1.5 1.5)>
+
+    """
+    return lib.shortest_line(a, b, **kwargs)

.venv/lib/python3.10/site-packages/shapely/measurement.py

@@ -0,0 +1,355 @@
+"""Methods for measuring (between) geometries."""
+
+import warnings
+
+import numpy as np
+
+from shapely import lib
+from shapely.decorators import multithreading_enabled
+
+__all__ = [
+    "area",
+    "bounds",
+    "distance",
+    "frechet_distance",
+    "hausdorff_distance",
+    "length",
+    "minimum_bounding_radius",
+    "minimum_clearance",
+    "total_bounds",
+]
+
+
+@multithreading_enabled
+def area(geometry, **kwargs):
+    """Compute the area of a (multi)polygon.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries for which to compute the area.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import MultiPolygon, Polygon
+    >>> polygon = Polygon([(0, 0), (0, 10), (10, 10), (10, 0), (0, 0)])
+    >>> shapely.area(polygon)
+    100.0
+    >>> polygon2 = Polygon([(10, 10), (10, 20), (20, 20), (20, 10), (10, 10)])
+    >>> shapely.area(MultiPolygon([polygon, polygon2]))
+    200.0
+    >>> shapely.area(Polygon())
+    0.0
+    >>> shapely.area(None)
+    nan
+
+    """
+    return lib.area(geometry, **kwargs)
+
+
+@multithreading_enabled
+def distance(a, b, **kwargs):
+    """Compute the Cartesian distance between two geometries.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to compute the distance between.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point, Polygon
+    >>> point = Point(0, 0)
+    >>> shapely.distance(Point(10, 0), point)
+    10.0
+    >>> shapely.distance(LineString([(1, 1), (1, -1)]), point)
+    1.0
+    >>> shapely.distance(Polygon([(3, 0), (5, 0), (5, 5), (3, 5), (3, 0)]), point)
+    3.0
+    >>> shapely.distance(Point(), point)
+    nan
+    >>> shapely.distance(None, point)
+    nan
+
+    """
+    return lib.distance(a, b, **kwargs)
+
+
+@multithreading_enabled
+def bounds(geometry, **kwargs):
+    """Compute the bounds (extent) of a geometry.
+
+    For each geometry these 4 numbers are returned: min x, min y, max x, max y.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries for which to compute the bounds.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point, Polygon
+    >>> shapely.bounds(Point(2, 3)).tolist()
+    [2.0, 3.0, 2.0, 3.0]
+    >>> shapely.bounds(LineString([(0, 0), (0, 2), (3, 2)])).tolist()
+    [0.0, 0.0, 3.0, 2.0]
+    >>> shapely.bounds(Polygon()).tolist()
+    [nan, nan, nan, nan]
+    >>> shapely.bounds(None).tolist()
+    [nan, nan, nan, nan]
+
+    """
+    return lib.bounds(geometry, **kwargs)
+
+
+def total_bounds(geometry, **kwargs):
+    """Compute the total bounds (extent) of the geometry.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries for which to compute the total bounds.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Returns
+    -------
+    numpy ndarray of [xmin, ymin, xmax, ymax]
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point, Polygon
+    >>> shapely.total_bounds(Point(2, 3)).tolist()
+    [2.0, 3.0, 2.0, 3.0]
+    >>> shapely.total_bounds([Point(2, 3), Point(4, 5)]).tolist()
+    [2.0, 3.0, 4.0, 5.0]
+    >>> shapely.total_bounds([
+    ...     LineString([(0, 1), (0, 2), (3, 2)]),
+    ...     LineString([(4, 4), (4, 6), (6, 7)])
+    ... ]).tolist()
+    [0.0, 1.0, 6.0, 7.0]
+    >>> shapely.total_bounds(Polygon()).tolist()
+    [nan, nan, nan, nan]
+    >>> shapely.total_bounds([Polygon(), Point(2, 3)]).tolist()
+    [2.0, 3.0, 2.0, 3.0]
+    >>> shapely.total_bounds(None).tolist()
+    [nan, nan, nan, nan]
+
+    """
+    b = bounds(geometry, **kwargs)
+    if b.ndim == 1:
+        return b
+
+    with warnings.catch_warnings():
+        # ignore 'All-NaN slice encountered' warnings
+        warnings.simplefilter("ignore", RuntimeWarning)
+        return np.array(
+            [
+                np.nanmin(b[..., 0]),
+                np.nanmin(b[..., 1]),
+                np.nanmax(b[..., 2]),
+                np.nanmax(b[..., 3]),
+            ]
+        )
+
+
+@multithreading_enabled
+def length(geometry, **kwargs):
+    """Compute the length of a (multi)linestring or polygon perimeter.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries for which to compute the length.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, MultiLineString, Polygon
+    >>> shapely.length(LineString([(0, 0), (0, 2), (3, 2)]))
+    5.0
+    >>> shapely.length(MultiLineString([
+    ...     LineString([(0, 0), (1, 0)]),
+    ...     LineString([(1, 0), (2, 0)])
+    ... ]))
+    2.0
+    >>> shapely.length(Polygon([(0, 0), (0, 10), (10, 10), (10, 0), (0, 0)]))
+    40.0
+    >>> shapely.length(LineString())
+    0.0
+    >>> shapely.length(None)
+    nan
+
+    """
+    return lib.length(geometry, **kwargs)
+
+
+@multithreading_enabled
+def hausdorff_distance(a, b, densify=None, **kwargs):
+    """Compute the discrete Hausdorff distance between two geometries.
+
+    The Hausdorff distance is a measure of similarity: it is the greatest
+    distance between any point in A and the closest point in B. The discrete
+    distance is an approximation of this metric: only vertices are considered.
+    The parameter 'densify' makes this approximation less coarse by splitting
+    the line segments between vertices before computing the distance.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to compute the distance between.
+    densify : float or array_like, optional
+        The value of densify is required to be between 0 and 1.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line1 = LineString([(130, 0), (0, 0), (0, 150)])
+    >>> line2 = LineString([(10, 10), (10, 150), (130, 10)])
+    >>> shapely.hausdorff_distance(line1, line2)
+    14.142135623730951
+    >>> shapely.hausdorff_distance(line1, line2, densify=0.5)
+    70.0
+    >>> shapely.hausdorff_distance(line1, LineString())
+    nan
+    >>> shapely.hausdorff_distance(line1, None)
+    nan
+
+    """
+    if densify is None:
+        return lib.hausdorff_distance(a, b, **kwargs)
+    else:
+        return lib.hausdorff_distance_densify(a, b, densify, **kwargs)
+
+
+@multithreading_enabled
+def frechet_distance(a, b, densify=None, **kwargs):
+    """Compute the discrete Fréchet distance between two geometries.
+
+    The Fréchet distance is a measure of similarity: it is the greatest
+    distance between any point in A and the closest point in B. The discrete
+    distance is an approximation of this metric: only vertices are considered.
+    The parameter 'densify' makes this approximation less coarse by splitting
+    the line segments between vertices before computing the distance.
+
+    Fréchet distance sweep continuously along their respective curves
+    and the direction of curves is significant. This makes it a better measure
+    of similarity than Hausdorff distance for curve or surface matching.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to compute the distance between.
+    densify : float or array_like, optional
+        The value of densify is required to be between 0 and 1.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line1 = LineString([(0, 0), (100, 0)])
+    >>> line2 = LineString([(0, 0), (50, 50), (100, 0)])
+    >>> shapely.frechet_distance(line1, line2)
+    70.71067811865476
+    >>> shapely.frechet_distance(line1, line2, densify=0.5)
+    50.0
+    >>> shapely.frechet_distance(line1, LineString())
+    nan
+    >>> shapely.frechet_distance(line1, None)
+    nan
+
+    """
+    if densify is None:
+        return lib.frechet_distance(a, b, **kwargs)
+    return lib.frechet_distance_densify(a, b, densify, **kwargs)
+
+
+@multithreading_enabled
+def minimum_clearance(geometry, **kwargs):
+    """Compute the Minimum Clearance distance.
+
+    A geometry's "minimum clearance" is the smallest distance by which
+    a vertex of the geometry could be moved to produce an invalid geometry.
+
+    If no minimum clearance exists for a geometry (for example, a single
+    point, or an empty geometry), infinity is returned.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries for which to compute the minimum clearance.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Polygon
+    >>> polygon = Polygon([(0, 0), (0, 10), (5, 6), (10, 10), (10, 0), (5, 4), (0, 0)])
+    >>> shapely.minimum_clearance(polygon)
+    2.0
+    >>> shapely.minimum_clearance(Polygon())
+    inf
+    >>> shapely.minimum_clearance(None)
+    nan
+
+    See Also
+    --------
+    minimum_clearance_line
+
+    """
+    return lib.minimum_clearance(geometry, **kwargs)
+
+
+@multithreading_enabled
+def minimum_bounding_radius(geometry, **kwargs):
+    """Compute the radius of the minimum bounding circle of an input geometry.
+
+    Parameters
+    ----------
+    geometry : Geometry or array_like
+        Geometry or geometries for which to compute the minimum bounding radius.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import GeometryCollection, LineString, MultiPoint, Point, Polygon
+    >>> shapely.minimum_bounding_radius(
+    ...     Polygon([(0, 5), (5, 10), (10, 5), (5, 0), (0, 5)])
+    ... )
+    5.0
+    >>> shapely.minimum_bounding_radius(LineString([(1, 1), (1, 10)]))
+    4.5
+    >>> shapely.minimum_bounding_radius(MultiPoint([(2, 2), (4, 2)]))
+    1.0
+    >>> shapely.minimum_bounding_radius(Point(0, 1))
+    0.0
+    >>> shapely.minimum_bounding_radius(GeometryCollection())
+    0.0
+
+    See Also
+    --------
+    minimum_bounding_circle
+
+    """
+    return lib.minimum_bounding_radius(geometry, **kwargs)

.venv/lib/python3.10/site-packages/shapely/ops.py

@@ -0,0 +1,720 @@
+"""Support for various GEOS geometry operations."""
+
+import shapely
+from shapely.algorithms.polylabel import polylabel  # noqa
+from shapely.errors import GeometryTypeError
+from shapely.geometry import (
+    GeometryCollection,
+    LineString,
+    MultiLineString,
+    MultiPoint,
+    Point,
+    Polygon,
+    shape,
+)
+from shapely.geometry.base import BaseGeometry
+from shapely.prepared import prep
+
+__all__ = [
+    "clip_by_rect",
+    "linemerge",
+    "nearest_points",
+    "operator",
+    "orient",
+    "polygonize",
+    "polygonize_full",
+    "shared_paths",
+    "snap",
+    "split",
+    "substring",
+    "transform",
+    "triangulate",
+    "unary_union",
+    "validate",
+    "voronoi_diagram",
+]
+
+
+class CollectionOperator:
+    def shapeup(self, ob):
+        if isinstance(ob, BaseGeometry):
+            return ob
+        else:
+            try:
+                return shape(ob)
+            except (ValueError, AttributeError):
+                return LineString(ob)
+
+    def polygonize(self, lines):
+        """Create polygons from a source of lines.
+
+        The source may be a MultiLineString, a sequence of LineString objects,
+        or a sequence of objects than can be adapted to LineStrings.
+        """
+        source = getattr(lines, "geoms", None) or lines
+        try:
+            source = iter(source)
+        except TypeError:
+            source = [source]
+        finally:
+            obs = [self.shapeup(line) for line in source]
+        collection = shapely.polygonize(obs)
+        return collection.geoms
+
+    def polygonize_full(self, lines):
+        """Create polygons from a source of lines.
+
+        The polygons and leftover geometries are returned as well.
+
+        The source may be a MultiLineString, a sequence of LineString objects,
+        or a sequence of objects than can be adapted to LineStrings.
+
+        Returns a tuple of objects: (polygons, cut edges, dangles, invalid ring
+        lines). Each are a geometry collection.
+
+        Dangles are edges which have one or both ends which are not incident on
+        another edge endpoint. Cut edges are connected at both ends but do not
+        form part of polygon. Invalid ring lines form rings which are invalid
+        (bowties, etc).
+        """
+        source = getattr(lines, "geoms", None) or lines
+        try:
+            source = iter(source)
+        except TypeError:
+            source = [source]
+        finally:
+            obs = [self.shapeup(line) for line in source]
+        return shapely.polygonize_full(obs)
+
+    def linemerge(self, lines, directed=False):
+        """Merge all connected lines from a source.
+
+        The source may be a MultiLineString, a sequence of LineString objects,
+        or a sequence of objects than can be adapted to LineStrings.  Returns a
+        LineString or MultiLineString when lines are not contiguous.
+        """
+        source = None
+        if getattr(lines, "geom_type", None) == "MultiLineString":
+            source = lines
+        elif hasattr(lines, "geoms"):
+            # other Multi geometries
+            source = MultiLineString([ls.coords for ls in lines.geoms])
+        elif hasattr(lines, "__iter__"):
+            try:
+                source = MultiLineString([ls.coords for ls in lines])
+            except AttributeError:
+                source = MultiLineString(lines)
+        if source is None:
+            raise ValueError(f"Cannot linemerge {lines}")
+        return shapely.line_merge(source, directed=directed)
+
+    def unary_union(self, geoms):
+        """Return the union of a sequence of geometries.
+
+        Usually used to convert a collection into the smallest set of polygons
+        that cover the same area.
+        """
+        return shapely.union_all(geoms, axis=None)
+
+
+operator = CollectionOperator()
+polygonize = operator.polygonize
+polygonize_full = operator.polygonize_full
+linemerge = operator.linemerge
+unary_union = operator.unary_union
+
+
+def triangulate(geom, tolerance=0.0, edges=False):
+    """Create the Delaunay triangulation and return a list of geometries.
+
+    The source may be any geometry type. All vertices of the geometry will be
+    used as the points of the triangulation.
+
+    From the GEOS documentation:
+    tolerance is the snapping tolerance used to improve the robustness of
+    the triangulation computation. A tolerance of 0.0 specifies that no
+    snapping will take place.
+
+    If edges is False, a list of Polygons (triangles) will be returned.
+    Otherwise the list of LineString edges is returned.
+
+    """
+    collection = shapely.delaunay_triangles(geom, tolerance=tolerance, only_edges=edges)
+    return list(collection.geoms)
+
+
+def voronoi_diagram(geom, envelope=None, tolerance=0.0, edges=False):
+    """Construct a Voronoi Diagram [1] from the given geometry.
+
+    Returns a list of geometries.
+
+    Parameters
+    ----------
+    geom: geometry
+        the input geometry whose vertices will be used to calculate
+        the final diagram.
+    envelope: geometry, None
+        clipping envelope for the returned diagram, automatically
+        determined if None. The diagram will be clipped to the larger
+        of this envelope or an envelope surrounding the sites.
+    tolerance: float, 0.0
+        sets the snapping tolerance used to improve the robustness
+        of the computation. A tolerance of 0.0 specifies that no
+        snapping will take place.
+    edges: bool, False
+        If False, return regions as polygons. Else, return only
+        edges e.g. LineStrings.
+
+    GEOS documentation can be found at [2]
+
+    Returns
+    -------
+    GeometryCollection
+        geometries representing the Voronoi regions.
+
+    Notes
+    -----
+    The tolerance `argument` can be finicky and is known to cause the
+    algorithm to fail in several cases. If you're using `tolerance`
+    and getting a failure, try removing it. The test cases in
+    tests/test_voronoi_diagram.py show more details.
+
+
+    References
+    ----------
+    [1] https://en.wikipedia.org/wiki/Voronoi_diagram
+    [2] https://geos.osgeo.org/doxygen/geos__c_8h_source.html  (line 730)
+
+    """
+    try:
+        result = shapely.voronoi_polygons(
+            geom, tolerance=tolerance, extend_to=envelope, only_edges=edges
+        )
+    except shapely.GEOSException as err:
+        errstr = "Could not create Voronoi Diagram with the specified inputs "
+        errstr += f"({err!s})."
+        if tolerance:
+            errstr += " Try running again with default tolerance value."
+        raise ValueError(errstr) from err
+
+    if result.geom_type != "GeometryCollection":
+        return GeometryCollection([result])
+    return result
+
+
+def validate(geom):
+    """Return True if the geometry is valid."""
+    return shapely.is_valid_reason(geom)
+
+
+def transform(func, geom):
+    """Apply `func` to all coordinates of `geom`.
+
+    Returns a new geometry of the same type from the transformed coordinates.
+
+    `func` maps x, y, and optionally z to output xp, yp, zp. The input
+    parameters may iterable types like lists or arrays or single values.
+    The output shall be of the same type. Scalars in, scalars out.
+    Lists in, lists out.
+
+    For example, here is an identity function applicable to both types
+    of input.
+
+      def id_func(x, y, z=None):
+          return tuple(filter(None, [x, y, z]))
+
+      g2 = transform(id_func, g1)
+
+    Using pyproj >= 2.1, this example will accurately project Shapely geometries:
+
+      import pyproj
+
+      wgs84 = pyproj.CRS('EPSG:4326')
+      utm = pyproj.CRS('EPSG:32618')
+
+      project = pyproj.Transformer.from_crs(wgs84, utm, always_xy=True).transform
+
+      g2 = transform(project, g1)
+
+    Note that the always_xy kwarg is required here as Shapely geometries only support
+    X,Y coordinate ordering.
+
+    Lambda expressions such as the one in
+
+      g2 = transform(lambda x, y, z=None: (x+1.0, y+1.0), g1)
+
+    also satisfy the requirements for `func`.
+    """
+    if geom.is_empty:
+        return geom
+    if geom.geom_type in ("Point", "LineString", "LinearRing", "Polygon"):
+        # First we try to apply func to x, y, z sequences. When func is
+        # optimized for sequences, this is the fastest, though zipping
+        # the results up to go back into the geometry constructors adds
+        # extra cost.
+        try:
+            if geom.geom_type in ("Point", "LineString", "LinearRing"):
+                return type(geom)(zip(*func(*zip(*geom.coords))))
+            elif geom.geom_type == "Polygon":
+                shell = type(geom.exterior)(zip(*func(*zip(*geom.exterior.coords))))
+                holes = [
+                    type(ring)(zip(*func(*zip(*ring.coords))))
+                    for ring in geom.interiors
+                ]
+                return type(geom)(shell, holes)
+
+        # A func that assumes x, y, z are single values will likely raise a
+        # TypeError, in which case we'll try again.
+        except TypeError:
+            if geom.geom_type in ("Point", "LineString", "LinearRing"):
+                return type(geom)([func(*c) for c in geom.coords])
+            elif geom.geom_type == "Polygon":
+                shell = type(geom.exterior)([func(*c) for c in geom.exterior.coords])
+                holes = [
+                    type(ring)([func(*c) for c in ring.coords])
+                    for ring in geom.interiors
+                ]
+                return type(geom)(shell, holes)
+
+    elif geom.geom_type.startswith("Multi") or geom.geom_type == "GeometryCollection":
+        return type(geom)([transform(func, part) for part in geom.geoms])
+    else:
+        raise GeometryTypeError(f"Type {geom.geom_type!r} not recognized")
+
+
+def nearest_points(g1, g2):
+    """Return the calculated nearest points in the input geometries.
+
+    The points are returned in the same order as the input geometries.
+    """
+    seq = shapely.shortest_line(g1, g2)
+    if seq is None:
+        if g1.is_empty:
+            raise ValueError("The first input geometry is empty")
+        else:
+            raise ValueError("The second input geometry is empty")
+
+    p1 = shapely.get_point(seq, 0)
+    p2 = shapely.get_point(seq, 1)
+    return (p1, p2)
+
+
+def snap(g1, g2, tolerance):
+    """Snaps an input geometry (g1) to reference (g2) geometry's vertices.
+
+    Parameters
+    ----------
+    g1 : geometry
+        The first geometry
+    g2 : geometry
+        The second geometry
+    tolerance : float
+        The snapping tolerance
+
+    Refer to :func:`shapely.snap` for full documentation.
+
+    """
+    return shapely.snap(g1, g2, tolerance)
+
+
+def shared_paths(g1, g2):
+    """Find paths shared between the two given lineal geometries.
+
+    Returns a GeometryCollection with two elements:
+     - First element is a MultiLineString containing shared paths with the
+       same direction for both inputs.
+     - Second element is a MultiLineString containing shared paths with the
+       opposite direction for the two inputs.
+
+    Parameters
+    ----------
+    g1 : geometry
+        The first geometry
+    g2 : geometry
+        The second geometry
+
+    """
+    if not isinstance(g1, LineString):
+        raise GeometryTypeError("First geometry must be a LineString")
+    if not isinstance(g2, LineString):
+        raise GeometryTypeError("Second geometry must be a LineString")
+    return shapely.shared_paths(g1, g2)
+
+
+class SplitOp:
+    @staticmethod
+    def _split_polygon_with_line(poly, splitter):
+        """Split a Polygon with a LineString."""
+        if not isinstance(poly, Polygon):
+            raise GeometryTypeError("First argument must be a Polygon")
+        if not isinstance(splitter, (LineString, MultiLineString)):
+            raise GeometryTypeError("Second argument must be a (Multi)LineString")
+
+        union = poly.boundary.union(splitter)
+
+        # greatly improves split performance for big geometries with many
+        # holes (the following contains checks) with minimal overhead
+        # for common cases
+        poly = prep(poly)
+
+        # some polygonized geometries may be holes, we do not want them
+        # that's why we test if the original polygon (poly) contains
+        # an inner point of polygonized geometry (pg)
+        return [
+            pg for pg in polygonize(union) if poly.contains(pg.representative_point())
+        ]
+
+    @staticmethod
+    def _split_line_with_line(line, splitter):
+        """Split a LineString with another (Multi)LineString or (Multi)Polygon."""
+        # if splitter is a polygon, pick it's boundary
+        if splitter.geom_type in ("Polygon", "MultiPolygon"):
+            splitter = splitter.boundary
+
+        if not isinstance(line, LineString):
+            raise GeometryTypeError("First argument must be a LineString")
+        if not isinstance(splitter, LineString) and not isinstance(
+            splitter, MultiLineString
+        ):
+            raise GeometryTypeError(
+                "Second argument must be either a LineString or a MultiLineString"
+            )
+
+        # |    s\l   | Interior | Boundary | Exterior |
+        # |----------|----------|----------|----------|
+        # | Interior |  0 or F  |    *     |    *     |   At least one of these two must be 0  # noqa: E501
+        # | Boundary |  0 or F  |    *     |    *     |   So either '0********' or '[0F]**0*****'  # noqa: E501
+        # | Exterior |    *     |    *     |    *     |   No overlapping interiors ('1********')  # noqa: E501
+        relation = splitter.relate(line)
+        if relation[0] == "1":
+            # The lines overlap at some segment (linear intersection of interiors)
+            raise ValueError("Input geometry segment overlaps with the splitter.")
+        elif relation[0] == "0" or relation[3] == "0":
+            # The splitter crosses or touches the line's interior
+            # --> return multilinestring from the split
+            return line.difference(splitter)
+        else:
+            # The splitter does not cross or touch the line's interior
+            # --> return collection with identity line
+            return [line]
+
+    @staticmethod
+    def _split_line_with_point(line, splitter):
+        """Split a LineString with a Point."""
+        if not isinstance(line, LineString):
+            raise GeometryTypeError("First argument must be a LineString")
+        if not isinstance(splitter, Point):
+            raise GeometryTypeError("Second argument must be a Point")
+
+        # check if point is in the interior of the line
+        if not line.relate_pattern(splitter, "0********"):
+            # point not on line interior --> return collection with single identity line
+            # (REASONING: Returning a list with the input line reference and creating a
+            # GeometryCollection at the general split function prevents unnecessary
+            # copying of linestrings in multipoint splitting function)
+            return [line]
+        elif line.coords[0] == splitter.coords[0]:
+            # if line is a closed ring the previous test doesn't behave as desired
+            return [line]
+
+        # point is on line, get the distance from the first point on line
+        distance_on_line = line.project(splitter)
+        coords = list(line.coords)
+        # split the line at the point and create two new lines
+        current_position = 0.0
+        for i in range(len(coords) - 1):
+            point1 = coords[i]
+            point2 = coords[i + 1]
+            dx = point1[0] - point2[0]
+            dy = point1[1] - point2[1]
+            segment_length = (dx**2 + dy**2) ** 0.5
+            current_position += segment_length
+            if distance_on_line == current_position:
+                # splitter is exactly on a vertex
+                return [LineString(coords[: i + 2]), LineString(coords[i + 1 :])]
+            elif distance_on_line < current_position:
+                # splitter is between two vertices
+                return [
+                    LineString(coords[: i + 1] + [splitter.coords[0]]),
+                    LineString([splitter.coords[0]] + coords[i + 1 :]),
+                ]
+        return [line]
+
+    @staticmethod
+    def _split_line_with_multipoint(line, splitter):
+        """Split a LineString with a MultiPoint."""
+        if not isinstance(line, LineString):
+            raise GeometryTypeError("First argument must be a LineString")
+        if not isinstance(splitter, MultiPoint):
+            raise GeometryTypeError("Second argument must be a MultiPoint")
+
+        chunks = [line]
+        for pt in splitter.geoms:
+            new_chunks = []
+            for chunk in filter(lambda x: not x.is_empty, chunks):
+                # add the newly split 2 lines or the same line if not split
+                new_chunks.extend(SplitOp._split_line_with_point(chunk, pt))
+            chunks = new_chunks
+
+        return chunks
+
+    @staticmethod
+    def split(geom, splitter):
+        """Split a geometry by another geometry and return a collection of geometries.
+
+        This function is the theoretical opposite of the union of
+        the split geometry parts. If the splitter does not split the geometry, a
+        collection with a single geometry equal to the input geometry is
+        returned.
+
+        The function supports:
+          - Splitting a (Multi)LineString by a (Multi)Point or (Multi)LineString
+            or (Multi)Polygon
+          - Splitting a (Multi)Polygon by a LineString
+
+        It may be convenient to snap the splitter with low tolerance to the
+        geometry. For example in the case of splitting a line by a point, the
+        point must be exactly on the line, for the line to be correctly split.
+        When splitting a line by a polygon, the boundary of the polygon is used
+        for the operation. When splitting a line by another line, a ValueError
+        is raised if the two overlap at some segment.
+
+        Parameters
+        ----------
+        geom : geometry
+            The geometry to be split
+        splitter : geometry
+            The geometry that will split the input geom
+
+        Examples
+        --------
+        >>> import shapely.ops
+        >>> from shapely import Point, LineString
+        >>> pt = Point((1, 1))
+        >>> line = LineString([(0,0), (2,2)])
+        >>> result = shapely.ops.split(line, pt)
+        >>> result.wkt
+        'GEOMETRYCOLLECTION (LINESTRING (0 0, 1 1), LINESTRING (1 1, 2 2))'
+
+        """
+        if geom.geom_type in ("MultiLineString", "MultiPolygon"):
+            return GeometryCollection(
+                [i for part in geom.geoms for i in SplitOp.split(part, splitter).geoms]
+            )
+
+        elif geom.geom_type == "LineString":
+            if splitter.geom_type in (
+                "LineString",
+                "MultiLineString",
+                "Polygon",
+                "MultiPolygon",
+            ):
+                split_func = SplitOp._split_line_with_line
+            elif splitter.geom_type == "Point":
+                split_func = SplitOp._split_line_with_point
+            elif splitter.geom_type == "MultiPoint":
+                split_func = SplitOp._split_line_with_multipoint
+            else:
+                raise GeometryTypeError(
+                    f"Splitting a LineString with a {splitter.geom_type} is "
+                    "not supported"
+                )
+
+        elif geom.geom_type == "Polygon":
+            if splitter.geom_type in ("LineString", "MultiLineString"):
+                split_func = SplitOp._split_polygon_with_line
+            else:
+                raise GeometryTypeError(
+                    f"Splitting a Polygon with a {splitter.geom_type} is not supported"
+                )
+
+        else:
+            raise GeometryTypeError(
+                f"Splitting {geom.geom_type} geometry is not supported"
+            )
+
+        return GeometryCollection(split_func(geom, splitter))
+
+
+split = SplitOp.split
+
+
+def substring(geom, start_dist, end_dist, normalized=False):
+    """Return a line segment between specified distances along a LineString.
+
+    Negative distance values are taken as measured in the reverse
+    direction from the end of the geometry. Out-of-range index
+    values are handled by clamping them to the valid range of values.
+
+    If the start distance equals the end distance, a Point is returned.
+
+    If the start distance is actually beyond the end distance, then the
+    reversed substring is returned such that the start distance is
+    at the first coordinate.
+
+    Parameters
+    ----------
+    geom : LineString
+        The geometry to get a substring of.
+    start_dist : float
+        The distance along `geom` of the start of the substring.
+    end_dist : float
+        The distance along `geom` of the end of the substring.
+    normalized : bool, False
+        Whether the distance parameters are interpreted as a
+        fraction of the geometry's length.
+
+    Returns
+    -------
+    Union[Point, LineString]
+        The substring between `start_dist` and `end_dist` or a Point
+        if they are at the same location.
+
+    Raises
+    ------
+    TypeError
+        If `geom` is not a LineString.
+
+    Examples
+    --------
+    >>> from shapely.geometry import LineString
+    >>> from shapely.ops import substring
+    >>> ls = LineString((i, 0) for i in range(6))
+    >>> ls.wkt
+    'LINESTRING (0 0, 1 0, 2 0, 3 0, 4 0, 5 0)'
+    >>> substring(ls, start_dist=1, end_dist=3).wkt
+    'LINESTRING (1 0, 2 0, 3 0)'
+    >>> substring(ls, start_dist=3, end_dist=1).wkt
+    'LINESTRING (3 0, 2 0, 1 0)'
+    >>> substring(ls, start_dist=1, end_dist=-3).wkt
+    'LINESTRING (1 0, 2 0)'
+    >>> substring(ls, start_dist=0.2, end_dist=-0.6, normalized=True).wkt
+    'LINESTRING (1 0, 2 0)'
+
+    Returning a `Point` when `start_dist` and `end_dist` are at the
+    same location.
+
+    >>> substring(ls, 2.5, -2.5).wkt
+    'POINT (2.5 0)'
+
+    """
+    if not isinstance(geom, LineString):
+        raise GeometryTypeError(
+            "Can only calculate a substring of LineString geometries. "
+            f"A {geom.geom_type} was provided."
+        )
+
+    # Filter out cases in which to return a point
+    if start_dist == end_dist:
+        return geom.interpolate(start_dist, normalized=normalized)
+    elif not normalized and start_dist >= geom.length and end_dist >= geom.length:
+        return geom.interpolate(geom.length, normalized=normalized)
+    elif not normalized and -start_dist >= geom.length and -end_dist >= geom.length:
+        return geom.interpolate(0, normalized=normalized)
+    elif normalized and start_dist >= 1 and end_dist >= 1:
+        return geom.interpolate(1, normalized=normalized)
+    elif normalized and -start_dist >= 1 and -end_dist >= 1:
+        return geom.interpolate(0, normalized=normalized)
+
+    if normalized:
+        start_dist *= geom.length
+        end_dist *= geom.length
+
+    # Filter out cases where distances meet at a middle point from opposite ends.
+    if start_dist < 0 < end_dist and abs(start_dist) + end_dist == geom.length:
+        return geom.interpolate(end_dist)
+    elif end_dist < 0 < start_dist and abs(end_dist) + start_dist == geom.length:
+        return geom.interpolate(start_dist)
+
+    start_point = geom.interpolate(start_dist)
+    end_point = geom.interpolate(end_dist)
+
+    if start_dist < 0:
+        start_dist = geom.length + start_dist  # Values may still be negative,
+    if end_dist < 0:  # but only in the out-of-range
+        end_dist = geom.length + end_dist  # sense, not the wrap-around sense.
+
+    reverse = start_dist > end_dist
+    if reverse:
+        start_dist, end_dist = end_dist, start_dist
+
+    start_dist = max(start_dist, 0)  # to avoid duplicating the first vertex
+
+    if reverse:
+        vertex_list = [tuple(*end_point.coords)]
+    else:
+        vertex_list = [tuple(*start_point.coords)]
+
+    coords = list(geom.coords)
+    current_distance = 0
+    for p1, p2 in zip(coords, coords[1:]):  # noqa
+        if start_dist < current_distance < end_dist:
+            vertex_list.append(p1)
+        elif current_distance >= end_dist:
+            break
+
+        current_distance += ((p2[0] - p1[0]) ** 2 + (p2[1] - p1[1]) ** 2) ** 0.5
+
+    if reverse:
+        vertex_list.append(tuple(*start_point.coords))
+        # reverse direction result
+        vertex_list = reversed(vertex_list)
+    else:
+        vertex_list.append(tuple(*end_point.coords))
+
+    return LineString(vertex_list)
+
+
+def clip_by_rect(geom, xmin, ymin, xmax, ymax):
+    """Return the portion of a geometry within a rectangle.
+
+    The geometry is clipped in a fast but possibly dirty way. The output is
+    not guaranteed to be valid. No exceptions will be raised for topological
+    errors.
+
+    Parameters
+    ----------
+    geom : geometry
+        The geometry to be clipped
+    xmin : float
+        Minimum x value of the rectangle
+    ymin : float
+        Minimum y value of the rectangle
+    xmax : float
+        Maximum x value of the rectangle
+    ymax : float
+        Maximum y value of the rectangle
+
+    Notes
+    -----
+    New in 1.7.
+
+    """
+    if geom.is_empty:
+        return geom
+    return shapely.clip_by_rect(geom, xmin, ymin, xmax, ymax)
+
+
+def orient(geom, sign=1.0):
+    """Return a properly oriented copy of the given geometry.
+
+    The signed area of the result will have the given sign. A sign of
+    1.0 means that the coordinates of the product's exterior rings will
+    be oriented counter-clockwise.
+
+    It is recommended to use :func:`shapely.orient_polygons` instead.
+
+    Parameters
+    ----------
+    geom : Geometry
+        The original geometry. May be a Polygon, MultiPolygon, or
+        GeometryCollection.
+    sign : float, optional.
+        The sign of the result's signed area.
+
+    Returns
+    -------
+    Geometry
+
+    """
+    return shapely.orient_polygons(geom, exterior_cw=sign < 0)

.venv/lib/python3.10/site-packages/shapely/plotting.py

@@ -0,0 +1,222 @@
+"""Plot single geometries using Matplotlib.
+
+Note: this module is experimental, and mainly targeting (interactive)
+exploration, debugging and illustration purposes.
+
+"""
+
+import numpy as np
+
+import shapely
+
+
+def _default_ax():
+    import matplotlib.pyplot as plt
+
+    ax = plt.gca()
+    ax.grid(True)
+    ax.set_aspect("equal")
+    return ax
+
+
+def _path_from_polygon(polygon):
+    from matplotlib.path import Path
+
+    from shapely.ops import orient
+
+    if isinstance(polygon, shapely.MultiPolygon):
+        return Path.make_compound_path(
+            *[_path_from_polygon(poly) for poly in polygon.geoms]
+        )
+    else:
+        polygon = orient(polygon)
+        return Path.make_compound_path(
+            Path(np.asarray(polygon.exterior.coords)[:, :2]),
+            *[Path(np.asarray(ring.coords)[:, :2]) for ring in polygon.interiors],
+        )
+
+
+def patch_from_polygon(polygon, **kwargs):
+    """Get a Matplotlib patch from a (Multi)Polygon.
+
+    Note: this function is experimental, and mainly targeting (interactive)
+    exploration, debugging and illustration purposes.
+
+    Parameters
+    ----------
+    polygon : shapely.Polygon or shapely.MultiPolygon
+        The polygon to convert to a Matplotlib Patch.
+    **kwargs
+        Additional keyword arguments passed to the matplotlib Patch.
+
+    Returns
+    -------
+    Matplotlib artist (PathPatch)
+
+    """
+    from matplotlib.patches import PathPatch
+
+    return PathPatch(_path_from_polygon(polygon), **kwargs)
+
+
+def plot_polygon(
+    polygon,
+    ax=None,
+    add_points=True,
+    color=None,
+    facecolor=None,
+    edgecolor=None,
+    linewidth=None,
+    **kwargs,
+):
+    """Plot a (Multi)Polygon.
+
+    Note: this function is experimental, and mainly targeting (interactive)
+    exploration, debugging and illustration purposes.
+
+    Parameters
+    ----------
+    polygon : shapely.Polygon or shapely.MultiPolygon
+        The polygon to plot.
+    ax : matplotlib Axes, default None
+        The axes on which to draw the plot. If not specified, will get the
+        current active axes or create a new figure.
+    add_points : bool, default True
+        If True, also plot the coordinates (vertices) as points.
+    color : matplotlib color specification
+        Color for both the polygon fill (face) and boundary (edge). By default,
+        the fill is using an alpha of 0.3. You can specify `facecolor` and
+        `edgecolor` separately for greater control.
+    facecolor : matplotlib color specification
+        Color for the polygon fill.
+    edgecolor : matplotlib color specification
+        Color for the polygon boundary.
+    linewidth : float
+        The line width for the polygon boundary.
+    **kwargs
+        Additional keyword arguments passed to the matplotlib Patch.
+
+    Returns
+    -------
+    Matplotlib artist (PathPatch), if `add_points` is false.
+    A tuple of Matplotlib artists (PathPatch, Line2D), if `add_points` is true.
+
+    """
+    from matplotlib import colors
+
+    if ax is None:
+        ax = _default_ax()
+
+    if color is None:
+        color = "C0"
+    color = colors.to_rgba(color)
+
+    if facecolor is None:
+        facecolor = list(color)
+        facecolor[-1] = 0.3
+        facecolor = tuple(facecolor)
+
+    if edgecolor is None:
+        edgecolor = color
+
+    patch = patch_from_polygon(
+        polygon, facecolor=facecolor, edgecolor=edgecolor, linewidth=linewidth, **kwargs
+    )
+    ax.add_patch(patch)
+    ax.autoscale_view()
+
+    if add_points:
+        line = plot_points(polygon, ax=ax, color=color)
+        return patch, line
+
+    return patch
+
+
+def plot_line(line, ax=None, add_points=True, color=None, linewidth=2, **kwargs):
+    """Plot a (Multi)LineString/LinearRing.
+
+    Note: this function is experimental, and mainly targeting (interactive)
+    exploration, debugging and illustration purposes.
+
+    Parameters
+    ----------
+    line : shapely.LineString or shapely.LinearRing
+        The line to plot.
+    ax : matplotlib Axes, default None
+        The axes on which to draw the plot. If not specified, will get the
+        current active axes or create a new figure.
+    add_points : bool, default True
+        If True, also plot the coordinates (vertices) as points.
+    color : matplotlib color specification
+        Color for the line (edgecolor under the hood) and points.
+    linewidth : float, default 2
+        The line width for the polygon boundary.
+    **kwargs
+        Additional keyword arguments passed to the matplotlib Patch.
+
+    Returns
+    -------
+    Matplotlib artist (PathPatch)
+
+    """
+    from matplotlib.patches import PathPatch
+    from matplotlib.path import Path
+
+    if ax is None:
+        ax = _default_ax()
+
+    if color is None:
+        color = "C0"
+
+    if isinstance(line, shapely.MultiLineString):
+        path = Path.make_compound_path(
+            *[Path(np.asarray(mline.coords)[:, :2]) for mline in line.geoms]
+        )
+    else:
+        path = Path(np.asarray(line.coords)[:, :2])
+
+    patch = PathPatch(
+        path, facecolor="none", edgecolor=color, linewidth=linewidth, **kwargs
+    )
+    ax.add_patch(patch)
+    ax.autoscale_view()
+
+    if add_points:
+        line = plot_points(line, ax=ax, color=color)
+        return patch, line
+
+    return patch
+
+
+def plot_points(geom, ax=None, color=None, marker="o", **kwargs):
+    """Plot a Point/MultiPoint or the vertices of any other geometry type.
+
+    Parameters
+    ----------
+    geom : shapely.Geometry
+        Any shapely Geometry object, from which all vertices are extracted
+        and plotted.
+    ax : matplotlib Axes, default None
+        The axes on which to draw the plot. If not specified, will get the
+        current active axes or create a new figure.
+    color : matplotlib color specification
+        Color for the filled points. You can use `markeredgecolor` and
+        `markerfacecolor` to have different edge and fill colors.
+    marker : str, default "o"
+        The matplotlib marker for the points.
+    **kwargs
+        Additional keyword arguments passed to matplotlib `plot` (Line2D).
+
+    Returns
+    -------
+    Matplotlib artist (Line2D)
+
+    """
+    if ax is None:
+        ax = _default_ax()
+
+    coords = shapely.get_coordinates(geom)
+    (line,) = ax.plot(
+        coords[:, 0], coords[:, 1], linestyle="", marker=marker, color=color, **kwargs
+    )
+    return line

.venv/lib/python3.10/site-packages/shapely/prepared.py

@@ -0,0 +1,74 @@
+"""Support for GEOS prepared geometry operations."""
+
+from pickle import PicklingError
+
+import shapely
+
+
+class PreparedGeometry:
+    """A geometry prepared for efficient comparison to a set of other geometries.
+
+    Examples
+    --------
+    >>> from shapely.prepared import prep
+    >>> from shapely.geometry import Point, Polygon
+    >>> triangle = Polygon([(0.0, 0.0), (1.0, 1.0), (1.0, -1.0)])
+    >>> p = prep(triangle)
+    >>> p.intersects(Point(0.5, 0.5))
+    True
+
+    """
+
+    def __init__(self, context):
+        """Prepare a geometry for efficient comparison to other geometries."""
+        if isinstance(context, PreparedGeometry):
+            self.context = context.context
+        else:
+            shapely.prepare(context)
+            self.context = context
+        self.prepared = True
+
+    def contains(self, other):
+        """Return True if the geometry contains the other, else False."""
+        return self.context.contains(other)
+
+    def contains_properly(self, other):
+        """Return True if the geometry properly contains the other, else False."""
+        return self.context.contains_properly(other)
+
+    def covers(self, other):
+        """Return True if the geometry covers the other, else False."""
+        return self.context.covers(other)
+
+    def crosses(self, other):
+        """Return True if the geometries cross, else False."""
+        return self.context.crosses(other)
+
+    def disjoint(self, other):
+        """Return True if geometries are disjoint, else False."""
+        return self.context.disjoint(other)
+
+    def intersects(self, other):
+        """Return True if geometries intersect, else False."""
+        return self.context.intersects(other)
+
+    def overlaps(self, other):
+        """Return True if geometries overlap, else False."""
+        return self.context.overlaps(other)
+
+    def touches(self, other):
+        """Return True if geometries touch, else False."""
+        return self.context.touches(other)
+
+    def within(self, other):
+        """Return True if geometry is within the other, else False."""
+        return self.context.within(other)
+
+    def __reduce__(self):
+        """Pickling is not supported."""
+        raise PicklingError("Prepared geometries cannot be pickled.")
+
+
+def prep(ob):
+    """Create and return a prepared geometric object."""
+    return PreparedGeometry(ob)

.venv/lib/python3.10/site-packages/shapely/set_operations.py

@@ -0,0 +1,784 @@
+"""Set-theoretic operations on geometry objects."""
+
+import warnings
+
+import numpy as np
+
+from shapely import Geometry, GeometryType, lib
+from shapely.decorators import (
+    deprecate_positional,
+    multithreading_enabled,
+    requires_geos,
+)
+
+__all__ = [
+    "coverage_union",
+    "coverage_union_all",
+    "difference",
+    "disjoint_subset_union",
+    "disjoint_subset_union_all",
+    "intersection",
+    "intersection_all",
+    "symmetric_difference",
+    "symmetric_difference_all",
+    "unary_union",
+    "union",
+    "union_all",
+]
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   difference(a, b, grid_size=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'grid_size' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'b'
+#   difference(a, b, *, grid_size=None, **kwargs)
+
+
+@deprecate_positional(["grid_size"], category=DeprecationWarning)
+@multithreading_enabled
+def difference(a, b, grid_size=None, **kwargs):
+    """Return the part of geometry A that does not intersect with geometry B.
+
+    If grid_size is nonzero, input coordinates will be snapped to a precision
+    grid of that size and resulting coordinates will be snapped to that same
+    grid.  If 0, this operation will use double precision coordinates.  If None,
+    the highest precision of the inputs will be used, which may be previously
+    set using set_precision.  Note: returned geometry does not have precision
+    set unless specified previously by set_precision.
+
+    Parameters
+    ----------
+    a : Geometry or array_like
+        Geometry or geometries to subtract b from.
+    b : Geometry or array_like
+        Geometry or geometries to subtract from a.
+    grid_size : float, optional
+        Precision grid size; will use the highest precision of the inputs by default.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``grid_size`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    set_precision
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line = LineString([(0, 0), (2, 2)])
+    >>> shapely.difference(line, LineString([(1, 1), (3, 3)]))
+    <LINESTRING (0 0, 1 1)>
+    >>> shapely.difference(line, LineString())
+    <LINESTRING (0 0, 2 2)>
+    >>> shapely.difference(line, None) is None
+    True
+    >>> box1 = shapely.box(0, 0, 2, 2)
+    >>> box2 = shapely.box(1, 1, 3, 3)
+    >>> shapely.difference(box1, box2).normalize()
+    <POLYGON ((0 0, 0 2, 1 2, 1 1, 2 1, 2 0, 0 0))>
+    >>> box1 = shapely.box(0.1, 0.2, 2.1, 2.1)
+    >>> shapely.difference(box1, box2, grid_size=1)
+    <POLYGON ((2 0, 0 0, 0 2, 1 2, 1 1, 2 1, 2 0))>
+
+    """
+    if grid_size is not None:
+        if not np.isscalar(grid_size):
+            raise ValueError("grid_size parameter only accepts scalar values")
+
+        return lib.difference_prec(a, b, grid_size, **kwargs)
+
+    return lib.difference(a, b, **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   intersection(a, b, grid_size=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'grid_size' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'b'
+#   intersection(a, b, *, grid_size=None, **kwargs)
+
+
+@deprecate_positional(["grid_size"], category=DeprecationWarning)
+@multithreading_enabled
+def intersection(a, b, grid_size=None, **kwargs):
+    """Return the geometry that is shared between input geometries.
+
+    If grid_size is nonzero, input coordinates will be snapped to a precision
+    grid of that size and resulting coordinates will be snapped to that same
+    grid.  If 0, this operation will use double precision coordinates.  If None,
+    the highest precision of the inputs will be used, which may be previously
+    set using set_precision.  Note: returned geometry does not have precision
+    set unless specified previously by set_precision.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to intersect with.
+    grid_size : float, optional
+        Precision grid size; will use the highest precision of the inputs by default.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``grid_size`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    intersection_all
+    set_precision
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line = LineString([(0, 0), (2, 2)])
+    >>> shapely.intersection(line, LineString([(1, 1), (3, 3)]))
+    <LINESTRING (1 1, 2 2)>
+    >>> box1 = shapely.box(0, 0, 2, 2)
+    >>> box2 = shapely.box(1, 1, 3, 3)
+    >>> shapely.intersection(box1, box2).normalize()
+    <POLYGON ((1 1, 1 2, 2 2, 2 1, 1 1))>
+    >>> box1 = shapely.box(0.1, 0.2, 2.1, 2.1)
+    >>> shapely.intersection(box1, box2, grid_size=1)
+    <POLYGON ((2 2, 2 1, 1 1, 1 2, 2 2))>
+
+    """
+    if grid_size is not None:
+        if not np.isscalar(grid_size):
+            raise ValueError("grid_size parameter only accepts scalar values")
+
+        return lib.intersection_prec(a, b, grid_size, **kwargs)
+
+    return lib.intersection(a, b, **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   intersection_all(geometries, axis=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'axis' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'geometries'
+#   intersection_all(geometries, *, axis=None, **kwargs)
+
+
+@deprecate_positional(["axis"], category=DeprecationWarning)
+@multithreading_enabled
+def intersection_all(geometries, axis=None, **kwargs):
+    """Return the intersection of multiple geometries.
+
+    This function ignores None values when other Geometry elements are present.
+    If all elements of the given axis are None, an empty GeometryCollection is
+    returned.
+
+    Parameters
+    ----------
+    geometries : array_like
+        Geometries to calculate the intersection of.
+    axis : int, optional
+        Axis along which the operation is performed. The default (None)
+        performs the operation over all axes, returning a scalar value.
+        Axis may be negative, in which case it counts from the last to the
+        first axis.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``axis`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    intersection
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line1 = LineString([(0, 0), (2, 2)])
+    >>> line2 = LineString([(1, 1), (3, 3)])
+    >>> shapely.intersection_all([line1, line2])
+    <LINESTRING (1 1, 2 2)>
+    >>> shapely.intersection_all([[line1, line2, None]], axis=1).tolist()
+    [<LINESTRING (1 1, 2 2)>]
+    >>> shapely.intersection_all([line1, None])
+    <LINESTRING (0 0, 2 2)>
+
+    """
+    geometries = np.asarray(geometries)
+    if axis is None:
+        geometries = geometries.ravel()
+    else:
+        geometries = np.rollaxis(geometries, axis=axis, start=geometries.ndim)
+
+    return lib.intersection_all(geometries, **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   symmetric_difference(a, b, grid_size=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'grid_size' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'b'
+#   symmetric_difference(a, b, *, grid_size=None, **kwargs)
+
+
+@deprecate_positional(["grid_size"], category=DeprecationWarning)
+@multithreading_enabled
+def symmetric_difference(a, b, grid_size=None, **kwargs):
+    """Return the geometry with the portions of input geometries that do not intersect.
+
+    If grid_size is nonzero, input coordinates will be snapped to a precision
+    grid of that size and resulting coordinates will be snapped to that same
+    grid.  If 0, this operation will use double precision coordinates.  If None,
+    the highest precision of the inputs will be used, which may be previously
+    set using set_precision.  Note: returned geometry does not have precision
+    set unless specified previously by set_precision.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to evaluate symmetric difference with.
+    grid_size : float, optional
+        Precision grid size; will use the highest precision of the inputs by default.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``grid_size`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    symmetric_difference_all
+    set_precision
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line = LineString([(0, 0), (2, 2)])
+    >>> shapely.symmetric_difference(line, LineString([(1, 1), (3, 3)]))
+    <MULTILINESTRING ((0 0, 1 1), (2 2, 3 3))>
+    >>> box1 = shapely.box(0, 0, 2, 2)
+    >>> box2 = shapely.box(1, 1, 3, 3)
+    >>> shapely.symmetric_difference(box1, box2).normalize()
+    <MULTIPOLYGON (((1 2, 1 3, 3 3, 3 1, 2 1, 2 2, 1 2)), ((0 0, 0 2, 1 2, 1 1, ...>
+    >>> box1 = shapely.box(0.1, 0.2, 2.1, 2.1)
+    >>> shapely.symmetric_difference(box1, box2, grid_size=1)
+    <MULTIPOLYGON (((2 0, 0 0, 0 2, 1 2, 1 1, 2 1, 2 0)), ((2 2, 1 2, 1 3, 3 3, ...>
+
+    """
+    if grid_size is not None:
+        if not np.isscalar(grid_size):
+            raise ValueError("grid_size parameter only accepts scalar values")
+
+        return lib.symmetric_difference_prec(a, b, grid_size, **kwargs)
+
+    return lib.symmetric_difference(a, b, **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   symmetric_difference_all(geometries, axis=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'axis' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'geometries'
+#   symmetric_difference_all(geometries, *, axis=None, **kwargs)
+
+
+@deprecate_positional(["axis"], category=DeprecationWarning)
+@multithreading_enabled
+def symmetric_difference_all(geometries, axis=None, **kwargs):
+    """Return the symmetric difference of multiple geometries.
+
+    This function ignores None values when other Geometry elements are present.
+    If all elements of the given axis are None an empty GeometryCollection is
+    returned.
+
+    .. deprecated:: 2.1.0
+
+        This function behaves incorrectly and will be removed in a future
+        version. See https://github.com/shapely/shapely/issues/2027 for more
+        details.
+
+    Parameters
+    ----------
+    geometries : array_like
+        Geometries to calculate the combined symmetric difference of.
+    axis : int, optional
+        Axis along which the operation is performed. The default (None)
+        performs the operation over all axes, returning a scalar value.
+        Axis may be negative, in which case it counts from the last to the
+        first axis.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``axis`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    symmetric_difference
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line1 = LineString([(0, 0), (2, 2)])
+    >>> line2 = LineString([(1, 1), (3, 3)])
+    >>> shapely.symmetric_difference_all([line1, line2])
+    <MULTILINESTRING ((0 0, 1 1), (2 2, 3 3))>
+    >>> shapely.symmetric_difference_all([[line1, line2, None]], axis=1).tolist()
+    [<MULTILINESTRING ((0 0, 1 1), (2 2, 3 3))>]
+    >>> shapely.symmetric_difference_all([line1, None])
+    <LINESTRING (0 0, 2 2)>
+    >>> shapely.symmetric_difference_all([None, None])
+    <GEOMETRYCOLLECTION EMPTY>
+
+    """
+    warnings.warn(
+        "The symmetric_difference_all function behaves incorrectly and will be "
+        "removed in a future version. "
+        "See https://github.com/shapely/shapely/issues/2027 for more details.",
+        DeprecationWarning,
+        stacklevel=2,
+    )
+    geometries = np.asarray(geometries)
+    if axis is None:
+        geometries = geometries.ravel()
+    else:
+        geometries = np.rollaxis(geometries, axis=axis, start=geometries.ndim)
+
+    return lib.symmetric_difference_all(geometries, **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   union(a, b, grid_size=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'grid_size' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'b'
+#   union(a, b, *, grid_size=None, **kwargs)
+
+
+@deprecate_positional(["grid_size"], category=DeprecationWarning)
+@multithreading_enabled
+def union(a, b, grid_size=None, **kwargs):
+    """Merge geometries into one.
+
+    If grid_size is nonzero, input coordinates will be snapped to a precision
+    grid of that size and resulting coordinates will be snapped to that same
+    grid.  If 0, this operation will use double precision coordinates.  If None,
+    the highest precision of the inputs will be used, which may be previously
+    set using set_precision.  Note: returned geometry does not have precision
+    set unless specified previously by set_precision.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to merge (union).
+    grid_size : float, optional
+        Precision grid size; will use the highest precision of the inputs by default.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``grid_size`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    union_all
+    set_precision
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString
+    >>> line = LineString([(0, 0), (2, 2)])
+    >>> shapely.union(line, LineString([(2, 2), (3, 3)]))
+    <MULTILINESTRING ((0 0, 2 2), (2 2, 3 3))>
+    >>> shapely.union(line, None) is None
+    True
+    >>> box1 = shapely.box(0, 0, 2, 2)
+    >>> box2 = shapely.box(1, 1, 3, 3)
+    >>> shapely.union(box1, box2).normalize()
+    <POLYGON ((0 0, 0 2, 1 2, 1 3, 3 3, 3 1, 2 1, 2 0, 0 0))>
+    >>> box1 = shapely.box(0.1, 0.2, 2.1, 2.1)
+    >>> shapely.union(box1, box2, grid_size=1)
+    <POLYGON ((2 0, 0 0, 0 2, 1 2, 1 3, 3 3, 3 1, 2 1, 2 0))>
+
+    """
+    if grid_size is not None:
+        if not np.isscalar(grid_size):
+            raise ValueError("grid_size parameter only accepts scalar values")
+
+        return lib.union_prec(a, b, grid_size, **kwargs)
+
+    return lib.union(a, b, **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   union_all(geometries, grid_size=None, axis=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'grid_size' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'geometries'
+#   union_all(geometries, *, grid_size=None, axis=None, **kwargs)
+
+
+@deprecate_positional(["grid_size", "axis"], category=DeprecationWarning)
+@multithreading_enabled
+def union_all(geometries, grid_size=None, axis=None, **kwargs):
+    """Return the union of multiple geometries.
+
+    This function ignores None values when other Geometry elements are present.
+    If all elements of the given axis are None an empty GeometryCollection is
+    returned.
+
+    If grid_size is nonzero, input coordinates will be snapped to a precision
+    grid of that size and resulting coordinates will be snapped to that same
+    grid.  If 0, this operation will use double precision coordinates.  If None,
+    the highest precision of the inputs will be used, which may be previously
+    set using set_precision.  Note: returned geometry does not have precision
+    set unless specified previously by set_precision.
+
+    `unary_union` is an alias of `union_all`.
+
+    Parameters
+    ----------
+    geometries : array_like
+        Geometries to merge/union.
+    grid_size : float, optional
+        Precision grid size; will use the highest precision of the inputs by default.
+    axis : int, optional
+        Axis along which the operation is performed. The default (None)
+        performs the operation over all axes, returning a scalar value.
+        Axis may be negative, in which case it counts from the last to the
+        first axis.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``grid_size`` or ``axis`` are
+        specified as positional arguments. In a future release, these will
+        need to be specified as keyword arguments.
+
+    See Also
+    --------
+    union
+    set_precision
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import LineString, Point
+    >>> line1 = LineString([(0, 0), (2, 2)])
+    >>> line2 = LineString([(2, 2), (3, 3)])
+    >>> shapely.union_all([line1, line2])
+    <MULTILINESTRING ((0 0, 2 2), (2 2, 3 3))>
+    >>> shapely.union_all([[line1, line2, None]], axis=1).tolist()
+    [<MULTILINESTRING ((0 0, 2 2), (2 2, 3 3))>]
+    >>> box1 = shapely.box(0, 0, 2, 2)
+    >>> box2 = shapely.box(1, 1, 3, 3)
+    >>> shapely.union_all([box1, box2]).normalize()
+    <POLYGON ((0 0, 0 2, 1 2, 1 3, 3 3, 3 1, 2 1, 2 0, 0 0))>
+    >>> box1 = shapely.box(0.1, 0.2, 2.1, 2.1)
+    >>> shapely.union_all([box1, box2], grid_size=1)
+    <POLYGON ((2 0, 0 0, 0 2, 1 2, 1 3, 3 3, 3 1, 2 1, 2 0))>
+    >>> shapely.union_all([None, Point(0, 1)])
+    <POINT (0 1)>
+    >>> shapely.union_all([None, None])
+    <GEOMETRYCOLLECTION EMPTY>
+    >>> shapely.union_all([])
+    <GEOMETRYCOLLECTION EMPTY>
+
+    """
+    # for union_all, GEOS provides an efficient route through first creating
+    # GeometryCollections
+    # first roll the aggregation axis backwards
+    geometries = np.asarray(geometries)
+    if axis is None:
+        geometries = geometries.ravel()
+    else:
+        geometries = np.rollaxis(geometries, axis=axis, start=geometries.ndim)
+
+    # create_collection acts on the inner axis
+    collections = lib.create_collection(
+        geometries, np.intc(GeometryType.GEOMETRYCOLLECTION)
+    )
+
+    if grid_size is not None:
+        if not np.isscalar(grid_size):
+            raise ValueError("grid_size parameter only accepts scalar values")
+
+        return lib.unary_union_prec(collections, grid_size, **kwargs)
+
+    return lib.unary_union(collections, **kwargs)
+
+
+unary_union = union_all
+
+
+@multithreading_enabled
+def coverage_union(a, b, **kwargs):
+    """Merge multiple polygons into one.
+
+    This is an optimized version of union which assumes the polygons to be
+    non-overlapping.
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to merge (union).
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    See Also
+    --------
+    coverage_union_all
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Polygon
+    >>> polygon_1 = Polygon([(0, 0), (0, 1), (1, 1), (1, 0), (0, 0)])
+    >>> polygon_2 = Polygon([(1, 0), (1, 1), (2, 1), (2, 0), (1, 0)])
+    >>> shapely.coverage_union(polygon_1, polygon_2).normalize()
+    <POLYGON ((0 0, 0 1, 1 1, 2 1, 2 0, 1 0, 0 0))>
+
+    Union with None returns same polygon
+
+    >>> shapely.coverage_union(polygon_1, None).normalize()
+    <POLYGON ((0 0, 0 1, 1 1, 1 0, 0 0))>
+
+    """
+    return coverage_union_all([a, b], **kwargs)
+
+
+# Note: future plan is to change this signature over a few releases:
+# shapely 2.0:
+#   coverage_union_all(geometries, axis=None, **kwargs)
+# shapely 2.1: shows deprecation warning about positional 'axis' arg
+#   same signature as 2.0
+# shapely 2.2(?): enforce keyword-only arguments after 'geometries'
+#   coverage_union_all(geometries, *, axis=None, **kwargs)
+
+
+@deprecate_positional(["axis"], category=DeprecationWarning)
+@multithreading_enabled
+def coverage_union_all(geometries, axis=None, **kwargs):
+    """Return the union of multiple polygons of a geometry collection.
+
+    This is an optimized version of union which assumes the polygons
+    to be non-overlapping.
+
+    This function ignores None values when other Geometry elements are present.
+    If all elements of the given axis are None, an empty GeometryCollection is
+    returned (before GEOS 3.12 this was an empty MultiPolygon).
+
+    Parameters
+    ----------
+    geometries : array_like
+        Geometries to merge/union.
+    axis : int, optional
+        Axis along which the operation is performed. The default (None)
+        performs the operation over all axes, returning a scalar value.
+        Axis may be negative, in which case it counts from the last to the
+        first axis.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    Notes
+    -----
+
+    .. deprecated:: 2.1.0
+        A deprecation warning is shown if ``axis`` is specified as a
+        positional argument. This will need to be specified as a keyword
+        argument in a future release.
+
+    See Also
+    --------
+    coverage_union
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Polygon
+    >>> polygon_1 = Polygon([(0, 0), (0, 1), (1, 1), (1, 0), (0, 0)])
+    >>> polygon_2 = Polygon([(1, 0), (1, 1), (2, 1), (2, 0), (1, 0)])
+    >>> shapely.coverage_union_all([polygon_1, polygon_2]).normalize()
+    <POLYGON ((0 0, 0 1, 1 1, 2 1, 2 0, 1 0, 0 0))>
+    >>> shapely.coverage_union_all([polygon_1, None]).normalize()
+    <POLYGON ((0 0, 0 1, 1 1, 1 0, 0 0))>
+    >>> shapely.coverage_union_all([None, None]).normalize()
+    <GEOMETRYCOLLECTION EMPTY>
+
+    """
+    # coverage union in GEOS works over GeometryCollections
+    # first roll the aggregation axis backwards
+    geometries = np.asarray(geometries)
+    if axis is None:
+        geometries = geometries.ravel()
+    else:
+        geometries = np.rollaxis(
+            np.asarray(geometries), axis=axis, start=geometries.ndim
+        )
+    # create_collection acts on the inner axis
+    collections = lib.create_collection(
+        geometries, np.intc(GeometryType.GEOMETRYCOLLECTION)
+    )
+    return lib.coverage_union(collections, **kwargs)
+
+
+@requires_geos("3.12.0")
+@multithreading_enabled
+def disjoint_subset_union(a, b, **kwargs):
+    """Merge multiple polygons into one using algorithm optimised for subsets.
+
+    This is an optimized version of union which assumes inputs can be
+    divided into subsets that do not intersect.
+
+    If there is only one such subset, performance can be expected to be worse than
+    :func:`union`. As such, it is recommeded to use ``disjoint_subset_union`` with
+    GeometryCollections rather than individual geometries.
+
+    .. versionadded:: 2.1.0
+
+    Parameters
+    ----------
+    a, b : Geometry or array_like
+        Geometry or geometries to merge (union).
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    See Also
+    --------
+    union
+    coverage_union
+    disjoint_subset_union_all
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Polygon
+    >>> polygon_1 = Polygon([(0, 0), (0, 1), (1, 1), (1, 0), (0, 0)])
+    >>> polygon_2 = Polygon([(1, 0), (1, 1), (2, 1), (2, 0), (1, 0)])
+    >>> shapely.disjoint_subset_union(polygon_1, polygon_2).normalize()
+    <POLYGON ((0 0, 0 1, 1 1, 2 1, 2 0, 1 0, 0 0))>
+
+    Union with None returns same polygon:
+
+    >>> shapely.disjoint_subset_union(polygon_1, None).normalize()
+    <POLYGON ((0 0, 0 1, 1 1, 1 0, 0 0))>
+    """
+    if (isinstance(a, Geometry) or a is None) and (
+        isinstance(b, Geometry) or b is None
+    ):
+        pass
+    elif isinstance(a, Geometry) or a is None:
+        a = np.full_like(b, a)
+    elif isinstance(b, Geometry) or b is None:
+        b = np.full_like(a, b)
+    elif len(a) != len(b):
+        raise ValueError("Arrays a and b must have the same length")
+    return disjoint_subset_union_all([a, b], axis=0, **kwargs)
+
+
+@requires_geos("3.12.0")
+@multithreading_enabled
+def disjoint_subset_union_all(geometries, *, axis=None, **kwargs):
+    """Return the union of multiple polygons.
+
+    This is an optimized version of union which assumes inputs can be divided into
+    subsets that do not intersect.
+
+    If there is only one such subset, performance can be expected to be worse than
+    :func:`union_all`.
+
+    This function ignores None values when other Geometry elements are present.
+    If all elements of the given axis are None, an empty GeometryCollection is
+    returned.
+
+    .. versionadded:: 2.1.0
+
+    Parameters
+    ----------
+    geometries : array_like
+        Geometries to union.
+    axis : int, optional
+        Axis along which the operation is performed. The default (None)
+        performs the operation over all axes, returning a scalar value.
+        Axis may be negative, in which case it counts from the last to the
+        first axis.
+    **kwargs
+        See :ref:`NumPy ufunc docs <ufuncs.kwargs>` for other keyword arguments.
+
+    See Also
+    --------
+    coverage_union_all
+    union_all
+    disjoint_subset_union
+
+    Examples
+    --------
+    >>> import shapely
+    >>> from shapely import Polygon
+    >>> polygon_1 = Polygon([(0, 0), (0, 1), (1, 1), (1, 0), (0, 0)])
+    >>> polygon_2 = Polygon([(1, 0), (1, 1), (2, 1), (2, 0), (1, 0)])
+    >>> shapely.disjoint_subset_union_all([polygon_1, polygon_2]).normalize()
+    <POLYGON ((0 0, 0 1, 1 1, 2 1, 2 0, 1 0, 0 0))>
+    >>> shapely.disjoint_subset_union_all([polygon_1, None]).normalize()
+    <POLYGON ((0 0, 0 1, 1 1, 1 0, 0 0))>
+    >>> shapely.disjoint_subset_union_all([None, None]).normalize()
+    <GEOMETRYCOLLECTION EMPTY>
+    """
+    geometries = np.asarray(geometries)
+    if axis is None:
+        geometries = geometries.ravel()
+    else:
+        geometries = np.rollaxis(
+            np.asarray(geometries), axis=axis, start=geometries.ndim
+        )
+    # create_collection acts on the inner axis
+    collections = lib.create_collection(
+        geometries, np.intc(GeometryType.GEOMETRYCOLLECTION)
+    )
+
+    return lib.disjoint_subset_union(collections, **kwargs)

.venv/lib/python3.10/site-packages/shapely/speedups.py

@@ -0,0 +1,39 @@
+"""Speedups for Shapely geometry operations.
+
+.. deprecated:: 2.0
+        Deprecated in Shapely 2.0, and will be removed in a future version.
+
+"""
+
+import warnings
+
+__all__ = ["available", "disable", "enable", "enabled"]
+
+
+available = True
+enabled = True
+
+
+_MSG = (
+    "This function has no longer any effect, and will be removed in a "
+    "future release. Starting with Shapely 2.0, equivalent speedups are "
+    "always available"
+)
+
+
+def enable():
+    """Will be removed in a future release and has no longer any effect.
+
+    Previously, this function enabled cython-based speedups. Starting with
+    Shapely 2.0, equivalent speedups are available in every installation.
+    """
+    warnings.warn(_MSG, FutureWarning, stacklevel=2)
+
+
+def disable():
+    """Will be removed in a future release and has no longer any effect.
+
+    Previously, this function enabled cython-based speedups. Starting with
+    Shapely 2.0, equivalent speedups are available in every installation.
+    """
+    warnings.warn(_MSG, FutureWarning, stacklevel=2)

.venv/lib/python3.10/site-packages/shapely/strtree.py

@@ -0,0 +1,550 @@
+"""STRtree spatial index for efficient spatial queries."""
+
+from collections.abc import Iterable
+from typing import Any
+
+import numpy as np
+
+from shapely import lib
+from shapely._enum import ParamEnum
+from shapely.decorators import UnsupportedGEOSVersionError
+from shapely.geometry.base import BaseGeometry
+from shapely.predicates import is_empty, is_missing
+
+__all__ = ["STRtree"]
+
+
+class BinaryPredicate(ParamEnum):
+    """The enumeration of GEOS binary predicates types."""
+
+    intersects = 1
+    within = 2
+    contains = 3
+    overlaps = 4
+    crosses = 5
+    touches = 6
+    covers = 7
+    covered_by = 8
+    contains_properly = 9
+
+
+class STRtree:
+    """A query-only R-tree spatial index.
+
+    It is created using the Sort-Tile-Recursive (STR) [1]_ algorithm.
+
+    The tree indexes the bounding boxes of each geometry.  The tree is
+    constructed directly at initialization and nodes cannot be added or
+    removed after it has been created.
+
+    All operations return indices of the input geometries.  These indices
+    can be used to index into anything associated with the input geometries,
+    including the input geometries themselves, or custom items stored in
+    another object of the same length as the geometries.
+
+    Bounding boxes limited to two dimensions and are axis-aligned (equivalent to
+    the ``bounds`` property of a geometry); any Z values present in geometries
+    are ignored for purposes of indexing within the tree.
+
+    Any mixture of geometry types may be stored in the tree.
+
+    Note: the tree is more efficient for querying when there are fewer
+    geometries that have overlapping bounding boxes and where there is greater
+    similarity between the outer boundary of a geometry and its bounding box.
+    For example, a MultiPolygon composed of widely-spaced individual Polygons
+    will have a large overall bounding box compared to the boundaries of its
+    individual Polygons, and the bounding box may also potentially overlap many
+    other geometries within the tree.  This means that the resulting tree may be
+    less efficient to query than a tree constructed from individual Polygons.
+
+    Parameters
+    ----------
+    geoms : sequence
+        A sequence of geometry objects.
+    node_capacity : int, default 10
+        The maximum number of child nodes per parent node in the tree.
+
+    References
+    ----------
+    .. [1] Leutenegger, Scott T.; Edgington, Jeffrey M.; Lopez, Mario A.
+       (February 1997). "STR: A Simple and Efficient Algorithm for
+       R-Tree Packing".
+       https://ia600900.us.archive.org/27/items/nasa_techdoc_19970016975/19970016975.pdf
+
+    """
+
+    def __init__(self, geoms: Iterable[BaseGeometry], node_capacity: int = 10):
+        """Create a new STRtree spatial index."""
+        # Keep references to geoms in a copied array so that this array is not
+        # modified while the tree depends on it remaining the same
+        self._geometries = np.array(geoms, dtype=np.object_, copy=True)
+
+        # initialize GEOS STRtree
+        self._tree = lib.STRtree(self.geometries, node_capacity)
+
+    def __len__(self):
+        """Return the number of geometries in the tree."""
+        return self._tree.count
+
+    def __reduce__(self):
+        """Pickle support."""
+        return (STRtree, (self.geometries,))
+
+    @property
+    def geometries(self):
+        """Geometries stored in the tree in the order used to construct the tree.
+
+        The order of this array corresponds to the tree indices returned by
+        other STRtree methods.
+
+        Do not attempt to modify items in the returned array.
+
+        Returns
+        -------
+        ndarray of Geometry objects
+
+        """
+        return self._geometries
+
+    def query(self, geometry, predicate=None, distance=None):
+        """Get the index combinations of all possibly intersecting geometries.
+
+        Returns the integer indices of all combinations of each input geometry
+        and tree geometries where the bounding box of each input geometry
+        intersects the bounding box of a tree geometry.
+
+        If the input geometry is a scalar, this returns an array of shape (n, ) with
+        the indices of the matching tree geometries.  If the input geometry is an
+        array_like, this returns an array with shape (2,n) where the subarrays
+        correspond to the indices of the input geometries and indices of the
+        tree geometries associated with each.  To generate an array of pairs of
+        input geometry index and tree geometry index, simply transpose the
+        result.
+
+        If a predicate is provided, the tree geometries are first queried based
+        on the bounding box of the input geometry and then are further filtered
+        to those that meet the predicate when comparing the input geometry to
+        the tree geometry:
+        predicate(geometry, tree_geometry)
+
+        The 'dwithin' predicate requires GEOS >= 3.10.
+
+        Bounding boxes are limited to two dimensions and are axis-aligned
+        (equivalent to the ``bounds`` property of a geometry); any Z values
+        present in input geometries are ignored when querying the tree.
+
+        Any input geometry that is None or empty will never match geometries in
+        the tree.
+
+        Parameters
+        ----------
+        geometry : Geometry or array_like
+            Input geometries to query the tree and filter results using the
+            optional predicate.
+        predicate : {None, 'intersects', 'within', 'contains', 'overlaps', 'crosses',\
+'touches', 'covers', 'covered_by', 'contains_properly', 'dwithin'}, optional
+            The predicate to use for testing geometries from the tree
+            that are within the input geometry's bounding box.
+        distance : number or array_like, optional
+            Distances around each input geometry within which to query the tree
+            for the 'dwithin' predicate.  If array_like, shape must be
+            broadcastable to shape of geometry.  Required if predicate='dwithin'.
+
+        Returns
+        -------
+        ndarray with shape (n,) if geometry is a scalar
+            Contains tree geometry indices.
+
+        OR
+
+        ndarray with shape (2, n) if geometry is an array_like
+            The first subarray contains input geometry indices.
+            The second subarray contains tree geometry indices.
+
+        Examples
+        --------
+        >>> from shapely import box, Point, STRtree
+        >>> import numpy as np
+        >>> points = [Point(0, 0), Point(1, 1), Point(2,2), Point(3, 3)]
+        >>> tree = STRtree(points)
+
+        Query the tree using a scalar geometry:
+
+        >>> indices = tree.query(box(0, 0, 1, 1))
+        >>> indices.tolist()
+        [0, 1]
+
+        Query using an array of geometries:
+
+        >>> boxes = np.array([box(0, 0, 1, 1), box(2, 2, 3, 3)])
+        >>> arr_indices = tree.query(boxes)
+        >>> arr_indices.tolist()
+        [[0, 0, 1, 1], [0, 1, 2, 3]]
+
+        Or transpose to get all pairs of input and tree indices:
+
+        >>> arr_indices.T.tolist()
+        [[0, 0], [0, 1], [1, 2], [1, 3]]
+
+        Retrieve the tree geometries by results of query:
+
+        >>> tree.geometries.take(indices).tolist()
+        [<POINT (0 0)>, <POINT (1 1)>]
+
+        Retrieve all pairs of input and tree geometries:
+
+        >>> np.array([boxes.take(arr_indices[0]),\
+tree.geometries.take(arr_indices[1])]).T.tolist()
+        [[<POLYGON ((1 0, 1 1, 0 1, 0 0, 1 0))>, <POINT (0 0)>],
+         [<POLYGON ((1 0, 1 1, 0 1, 0 0, 1 0))>, <POINT (1 1)>],
+         [<POLYGON ((3 2, 3 3, 2 3, 2 2, 3 2))>, <POINT (2 2)>],
+         [<POLYGON ((3 2, 3 3, 2 3, 2 2, 3 2))>, <POINT (3 3)>]]
+
+        Query using a predicate:
+
+        >>> tree = STRtree([box(0, 0, 0.5, 0.5), box(0.5, 0.5, 1, 1), box(1, 1, 2, 2)])
+        >>> tree.query(box(0, 0, 1, 1), predicate="contains").tolist()
+        [0, 1]
+        >>> tree.query(Point(0.75, 0.75), predicate="dwithin", distance=0.5).tolist()
+        [0, 1, 2]
+
+        >>> tree.query(boxes, predicate="contains").tolist()
+        [[0, 0], [0, 1]]
+        >>> tree.query(boxes, predicate="dwithin", distance=0.5).tolist()
+        [[0, 0, 0, 1], [0, 1, 2, 2]]
+
+        Retrieve custom items associated with tree geometries (records can
+        be in whatever data structure so long as geometries and custom data
+        can be extracted into arrays of the same length and order):
+
+        >>> records = [
+        ...     {"geometry": Point(0, 0), "value": "A"},
+        ...     {"geometry": Point(2, 2), "value": "B"}
+        ... ]
+        >>> tree = STRtree([record["geometry"] for record in records])
+        >>> items = np.array([record["value"] for record in records])
+        >>> items.take(tree.query(box(0, 0, 1, 1))).tolist()
+        ['A']
+
+
+        Notes
+        -----
+        In the context of a spatial join, input geometries are the "left"
+        geometries that determine the order of the results, and tree geometries
+        are "right" geometries that are joined against the left geometries. This
+        effectively performs an inner join, where only those combinations of
+        geometries that can be joined based on overlapping bounding boxes or
+        optional predicate are returned.
+
+        """
+        geometry = np.asarray(geometry)
+        is_scalar = False
+        if geometry.ndim == 0:
+            geometry = np.expand_dims(geometry, 0)
+            is_scalar = True
+
+        if predicate is None:
+            indices = self._tree.query(geometry, 0)
+            return indices[1] if is_scalar else indices
+
+        # Requires GEOS >= 3.10
+        elif predicate == "dwithin":
+            if lib.geos_version < (3, 10, 0):
+                raise UnsupportedGEOSVersionError(
+                    "dwithin predicate requires GEOS >= 3.10"
+                )
+            if distance is None:
+                raise ValueError(
+                    "distance parameter must be provided for dwithin predicate"
+                )
+            distance = np.asarray(distance, dtype="float64")
+            if distance.ndim > 1:
+                raise ValueError("Distance array should be one dimensional")
+
+            try:
+                distance = np.broadcast_to(distance, geometry.shape)
+            except ValueError:
+                raise ValueError("Could not broadcast distance to match geometry")
+
+            indices = self._tree.dwithin(geometry, distance)
+            return indices[1] if is_scalar else indices
+
+        predicate = BinaryPredicate.get_value(predicate)
+        indices = self._tree.query(geometry, predicate)
+        return indices[1] if is_scalar else indices
+
+    def nearest(self, geometry) -> Any | None:
+        """Return the index of the nearest geometry in the tree.
+
+        This is determined for each input geometry based on distance within
+        two-dimensional Cartesian space.
+
+        This distance will be 0 when input geometries intersect tree geometries.
+
+        If there are multiple equidistant or intersected geometries in the tree,
+        only a single result is returned for each input geometry, based on the
+        order that tree geometries are visited; this order may be
+        nondeterministic.
+
+        If any input geometry is None or empty, an error is raised.  Any Z
+        values present in input geometries are ignored when finding nearest
+        tree geometries.
+
+        Parameters
+        ----------
+        geometry : Geometry or array_like
+            Input geometries to query the tree.
+
+        Returns
+        -------
+        scalar or ndarray
+            Indices of geometries in tree. Return value will have the same shape
+            as the input.
+
+            None is returned if this index is empty. This may change in
+            version 2.0.
+
+        See Also
+        --------
+        query_nearest: returns all equidistant geometries, exclusive geometries, \
+and optional distances
+
+        Examples
+        --------
+        >>> from shapely import Point, STRtree
+        >>> tree = STRtree([Point(i, i) for i in range(10)])
+
+        Query the tree for nearest using a scalar geometry:
+
+        >>> index = tree.nearest(Point(2.2, 2.2))
+        >>> index
+        2
+        >>> tree.geometries.take(index)
+        <POINT (2 2)>
+
+        Query the tree for nearest using an array of geometries:
+
+        >>> indices = tree.nearest([Point(2.2, 2.2), Point(4.4, 4.4)])
+        >>> indices.tolist()
+        [2, 4]
+        >>> tree.geometries.take(indices).tolist()
+        [<POINT (2 2)>, <POINT (4 4)>]
+
+        Nearest only return one object if there are multiple equidistant results:
+
+        >>> tree = STRtree ([Point(0, 0), Point(0, 0)])
+        >>> tree.nearest(Point(0, 0))
+        0
+
+        """
+        if self._tree.count == 0:
+            return None
+
+        geometry_arr = np.asarray(geometry, dtype=object)
+        if is_missing(geometry_arr).any() or is_empty(geometry_arr).any():
+            raise ValueError(
+                "Cannot determine nearest geometry for empty geometry or "
+                "missing value (None)."
+            )
+        # _tree.nearest returns ndarray with shape (2, 1) -> index in input
+        # geometries and index into tree geometries
+        indices = self._tree.nearest(np.atleast_1d(geometry_arr))[1]
+
+        if geometry_arr.ndim == 0:
+            return indices[0]
+        else:
+            return indices
+
+    def query_nearest(
+        self,
+        geometry,
+        max_distance=None,
+        return_distance=False,
+        exclusive=False,
+        all_matches=True,
+    ):
+        """Return the index of the nearest geometries in the tree.
+
+        This is determined for each input geometry based on distance within
+        two-dimensional Cartesian space.
+
+        This distance will be 0 when input geometries intersect tree geometries.
+
+        If there are multiple equidistant or intersected geometries in tree and
+        `all_matches` is True (the default), all matching tree geometries are
+        returned; otherwise only the first matching tree geometry is returned.
+        Tree indices are returned in the order they are visited for each input
+        geometry and may not be in ascending index order; no meaningful order is
+        implied.
+
+        The max_distance used to search for nearest items in the tree may have a
+        significant impact on performance by reducing the number of input
+        geometries that are evaluated for nearest items in the tree.  Only those
+        input geometries with at least one tree geometry within +/- max_distance
+        beyond their envelope will be evaluated.  However, using a large
+        max_distance may have a negative performance impact because many tree
+        geometries will be queried for each input geometry.
+
+        The distance, if returned, will be 0 for any intersected geometries in
+        the tree.
+
+        Any geometry that is None or empty in the input geometries is omitted
+        from the output.  Any Z values present in input geometries are ignored
+        when finding nearest tree geometries.
+
+        Parameters
+        ----------
+        geometry : Geometry or array_like
+            Input geometries to query the tree.
+        max_distance : float, optional
+            Maximum distance within which to query for nearest items in tree.
+            Must be greater than 0.
+        return_distance : bool, default False
+            If True, will return distances in addition to indices.
+        exclusive : bool, default False
+            If True, the nearest tree geometries that are equal to the input
+            geometry will not be returned.
+        all_matches : bool, default True
+            If True, all equidistant and intersected geometries will be returned
+            for each input geometry.
+            If False, only the first nearest geometry will be returned.
+
+        Returns
+        -------
+        tree indices or tuple of (tree indices, distances) if geometry is a scalar
+            indices is an ndarray of shape (n, ) and distances (if present) an
+            ndarray of shape (n, )
+
+        OR
+
+        indices or tuple of (indices, distances)
+            indices is an ndarray of shape (2,n) and distances (if present) an
+            ndarray of shape (n).
+            The first subarray of indices contains input geometry indices.
+            The second subarray of indices contains tree geometry indices.
+
+        See Also
+        --------
+        nearest: returns singular nearest geometry for each input
+
+        Examples
+        --------
+        >>> import numpy as np
+        >>> from shapely import box, Point, STRtree
+        >>> points = [Point(0, 0), Point(1, 1), Point(2,2), Point(3, 3)]
+        >>> tree = STRtree(points)
+
+        Find the nearest tree geometries to a scalar geometry:
+
+        >>> indices = tree.query_nearest(Point(0.25, 0.25))
+        >>> indices.tolist()
+        [0]
+
+        Retrieve the tree geometries by results of query:
+
+        >>> tree.geometries.take(indices).tolist()
+        [<POINT (0 0)>]
+
+        Find the nearest tree geometries to an array of geometries:
+
+        >>> query_points = np.array([Point(2.25, 2.25), Point(1, 1)])
+        >>> arr_indices = tree.query_nearest(query_points)
+        >>> arr_indices.tolist()
+        [[0, 1], [2, 1]]
+
+        Or transpose to get all pairs of input and tree indices:
+
+        >>> arr_indices.T.tolist()
+        [[0, 2], [1, 1]]
+
+        Retrieve all pairs of input and tree geometries:
+
+        >>> list(zip(query_points.take(arr_indices[0]), tree.geometries.take(arr_indices[1])))
+        [(<POINT (2.25 2.25)>, <POINT (2 2)>), (<POINT (1 1)>, <POINT (1 1)>)]
+
+        All intersecting geometries in the tree are returned by default:
+
+        >>> tree.query_nearest(box(1,1,3,3)).tolist()
+        [1, 2, 3]
+
+        Set all_matches to False to to return a single match per input geometry:
+
+        >>> tree.query_nearest(box(1,1,3,3), all_matches=False).tolist()
+        [1]
+
+        Return the distance to each nearest tree geometry:
+
+        >>> index, distance = tree.query_nearest(Point(0.5, 0.5), return_distance=True)
+        >>> index.tolist()
+        [0, 1]
+        >>> distance.round(4).tolist()
+        [0.7071, 0.7071]
+
+        Return the distance for each input and nearest tree geometry for an array
+        of geometries:
+
+        >>> indices, distance = tree.query_nearest([Point(0.5, 0.5), Point(1, 1)], return_distance=True)
+        >>> indices.tolist()
+        [[0, 0, 1], [0, 1, 1]]
+        >>> distance.round(4).tolist()
+        [0.7071, 0.7071, 0.0]
+
+        Retrieve custom items associated with tree geometries (records can
+        be in whatever data structure so long as geometries and custom data
+        can be extracted into arrays of the same length and order):
+
+        >>> records = [
+        ...     {"geometry": Point(0, 0), "value": "A"},
+        ...     {"geometry": Point(2, 2), "value": "B"}
+        ... ]
+        >>> tree = STRtree([record["geometry"] for record in records])
+        >>> items = np.array([record["value"] for record in records])
+        >>> items.take(tree.query_nearest(Point(0.5, 0.5))).tolist()
+        ['A']
+
+        """  # noqa: E501
+        geometry = np.asarray(geometry, dtype=object)
+        is_scalar = False
+        if geometry.ndim == 0:
+            geometry = np.expand_dims(geometry, 0)
+            is_scalar = True
+
+        if max_distance is not None:
+            if not np.isscalar(max_distance):
+                raise ValueError("max_distance parameter only accepts scalar values")
+
+            if max_distance <= 0:
+                raise ValueError("max_distance must be greater than 0")
+
+        # a distance of 0 means no max_distance is used
+        max_distance = max_distance or 0
+
+        if not np.isscalar(exclusive):
+            raise ValueError("exclusive parameter only accepts scalar values")
+
+        if exclusive not in {True, False}:
+            raise ValueError("exclusive parameter must be boolean")
+
+        if not np.isscalar(all_matches):
+            raise ValueError("all_matches parameter only accepts scalar values")
+
+        if all_matches not in {True, False}:
+            raise ValueError("all_matches parameter must be boolean")
+
+        results = self._tree.query_nearest(
+            geometry, max_distance, exclusive, all_matches
+        )
+
+        # output indices are shape (n, )
+        if is_scalar:
+            if not return_distance:
+                return results[0][1]
+
+            else:
+                return (results[0][1], results[1])
+
+        # output indices are shape (2, n)
+        if not return_distance:
+            return results[0]
+
+        return results

.venv/lib/python3.10/site-packages/shapely/testing.py

@@ -0,0 +1,208 @@
+"""Utilities for testing with shapely geometries."""
+
+from functools import partial
+
+import numpy as np
+
+import shapely
+
+__all__ = ["assert_geometries_equal"]
+
+
+def _equals_exact_with_ndim(x, y, tolerance):
+    dimension_equals = shapely.get_coordinate_dimension(
+        x
+    ) == shapely.get_coordinate_dimension(y)
+    with np.errstate(invalid="ignore"):
+        # Suppress 'invalid value encountered in equals_exact' with nan coordinates
+        geometry_equals = shapely.equals_exact(x, y, tolerance=tolerance)
+    return dimension_equals & geometry_equals
+
+
+def _replace_nan(arr):
+    return np.where(np.isnan(arr), 0.0, arr)
+
+
+def _assert_nan_coords_same(x, y, tolerance, err_msg, verbose):
+    x, y = np.broadcast_arrays(x, y)
+    x_coords = shapely.get_coordinates(x, include_z=True)
+    y_coords = shapely.get_coordinates(y, include_z=True)
+
+    # Check the shapes (condition is copied from numpy test_array_equal)
+    if x_coords.shape != y_coords.shape:
+        return False
+
+    # Check NaN positional equality
+    x_id = np.isnan(x_coords)
+    y_id = np.isnan(y_coords)
+    if not (x_id == y_id).all():
+        msg = build_err_msg(
+            [x, y],
+            err_msg + "\nx and y nan coordinate location mismatch:",
+            verbose=verbose,
+        )
+        raise AssertionError(msg)
+
+    # If this passed, replace NaN with a number to be able to use equals_exact
+    x_no_nan = shapely.transform(x, _replace_nan, include_z=True)
+    y_no_nan = shapely.transform(y, _replace_nan, include_z=True)
+
+    return _equals_exact_with_ndim(x_no_nan, y_no_nan, tolerance=tolerance)
+
+
+def _assert_none_same(x, y, err_msg, verbose):
+    x_id = shapely.is_missing(x)
+    y_id = shapely.is_missing(y)
+
+    if not (x_id == y_id).all():
+        msg = build_err_msg(
+            [x, y],
+            err_msg + "\nx and y None location mismatch:",
+            verbose=verbose,
+        )
+        raise AssertionError(msg)
+
+    # If there is a scalar, then here we know the array has the same
+    # flag as it everywhere, so we should return the scalar flag.
+    if x.ndim == 0:
+        return bool(x_id)
+    elif y.ndim == 0:
+        return bool(y_id)
+    else:
+        return y_id
+
+
+def assert_geometries_equal(
+    x,
+    y,
+    tolerance=1e-7,
+    equal_none=True,
+    equal_nan=True,
+    normalize=False,
+    err_msg="",
+    verbose=True,
+):
+    """Raise an AssertionError if two geometry array_like objects are not equal.
+
+    Given two array_like objects, check that the shape is equal and all elements
+    of these objects are equal. An exception is raised at shape mismatch or
+    conflicting values. In contrast to the standard usage in shapely, no
+    assertion is raised if both objects have NaNs/Nones in the same positions.
+
+    Parameters
+    ----------
+    x, y : Geometry or array_like
+        Geometry or geometries to compare.
+    tolerance: float, default 1e-7
+        The tolerance to use when comparing geometries.
+    equal_none : bool, default True
+        Whether to consider None elements equal to other None elements.
+    equal_nan : bool, default True
+        Whether to consider nan coordinates as equal to other nan coordinates.
+    normalize : bool, default False
+        Whether to normalize geometries prior to comparison.
+    err_msg : str, optional
+        The error message to be printed in case of failure.
+    verbose : bool, optional
+        If True, the conflicting values are appended to the error message.
+
+    """
+    __tracebackhide__ = True  # Hide traceback for py.test
+    if normalize:
+        x = shapely.normalize(x)
+        y = shapely.normalize(y)
+    x = np.asarray(x)
+    y = np.asarray(y)
+
+    is_scalar = x.ndim == 0 or y.ndim == 0
+
+    # Check the shapes (condition is copied from numpy test_array_equal)
+    if not (is_scalar or x.shape == y.shape):
+        msg = build_err_msg(
+            [x, y],
+            err_msg + f"\n(shapes {x.shape}, {y.shape} mismatch)",
+            verbose=verbose,
+        )
+        raise AssertionError(msg)
+
+    flagged = False
+    if equal_none:
+        flagged = _assert_none_same(x, y, err_msg, verbose)
+
+    if not np.isscalar(flagged):
+        x, y = x[~flagged], y[~flagged]
+        # Only do the comparison if actual values are left
+        if x.size == 0:
+            return
+    elif flagged:
+        # no sense doing comparison if everything is flagged.
+        return
+
+    is_equal = _equals_exact_with_ndim(x, y, tolerance=tolerance)
+    if is_scalar and not np.isscalar(is_equal):
+        is_equal = bool(is_equal[0])
+
+    if np.all(is_equal):
+        return
+    elif not equal_nan:
+        msg = build_err_msg(
+            [x, y],
+            err_msg + f"\nNot equal to tolerance {tolerance:g}",
+            verbose=verbose,
+        )
+        raise AssertionError(msg)
+
+    # Optionally refine failing elements if NaN should be considered equal
+    if not np.isscalar(is_equal):
+        x, y = x[~is_equal], y[~is_equal]
+        # Only do the NaN check if actual values are left
+        if x.size == 0:
+            return
+    elif is_equal:
+        # no sense in checking for NaN if everything is equal.
+        return
+
+    is_equal = _assert_nan_coords_same(x, y, tolerance, err_msg, verbose)
+    if not np.all(is_equal):
+        msg = build_err_msg(
+            [x, y],
+            err_msg + f"\nNot equal to tolerance {tolerance:g}",
+            verbose=verbose,
+        )
+        raise AssertionError(msg)
+
+
+## BELOW A COPY FROM numpy.testing._private.utils (numpy version 1.20.2)
+
+
+def build_err_msg(
+    arrays,
+    err_msg,
+    header="Geometries are not equal:",
+    verbose=True,
+    names=("x", "y"),
+    precision=8,
+):
+    msg = ["\n" + header]
+    if err_msg:
+        if err_msg.find("\n") == -1 and len(err_msg) < 79 - len(header):
+            msg = [msg[0] + " " + err_msg]
+        else:
+            msg.append(err_msg)
+    if verbose:
+        for i, a in enumerate(arrays):
+            if isinstance(a, np.ndarray):
+                # precision argument is only needed if the objects are ndarrays
+                r_func = partial(np.array_repr, precision=precision)
+            else:
+                r_func = repr
+
+            try:
+                r = r_func(a)
+            except Exception as exc:
+                r = f"[repr failed for <{type(a).__name__}>: {exc}]"
+            if r.count("\n") > 3:
+                r = "\n".join(r.splitlines()[:3])
+                r += "..."
+            msg.append(f" {names[i]}: {r}")
+    return "\n".join(msg)