Merge pull request #2455 from rcomer/path-shapely

Introduce `path_to_shapely` and `shapely_to_path`

.github/workflows/ci-testing.yml

@@ -41,7 +41,7 @@ jobs:
           (github.event_name == 'push' || github.event_name == 'pull_request')
         id: minimum-packages
         run: |
-          pip install cython==0.29.28 matplotlib==3.6 numpy==1.23 owslib==0.27 pyproj==3.3.1 scipy==1.9 shapely==1.8 pyshp==2.3.1
+          pip install cython==0.29.28 matplotlib==3.6 numpy==1.23 owslib==0.27 pyproj==3.3.1 scipy==1.9 shapely==2.0 pyshp==2.3.1
 
       - name: Coverage packages
         id: coverage

INSTALL

@@ -66,7 +66,7 @@ Further information about the required dependencies can be found here:
     Python package for 2D plotting. Python package required for any
     graphical capabilities.
 
-**Shapely** 1.8 or later (https://github.com/shapely/shapely)
+**Shapely** 2.0 or later (https://github.com/shapely/shapely)
     Python package for the manipulation and analysis of planar geometric objects.
 
 **pyshp** 2.3 or later (https://pypi.python.org/pypi/pyshp)

docs/source/reference/matplotlib.rst

@@ -64,7 +64,7 @@ Artist extensions
     SlippyImageArtist
 
 Patch
-~~~~~~~~~~~~~~~~~~~~~
+~~~~~
 
 .. automodule:: cartopy.mpl.patch
 
@@ -74,3 +74,14 @@ Patch
     geos_to_path
     path_segments
     path_to_geos
+
+Path
+~~~~
+
+.. automodule:: cartopy.mpl.path
+
+.. autosummary::
+    :toctree: generated/
+
+    path_to_shapely
+    shapely_to_path

docs/source/whatsnew/index.rst

@@ -14,6 +14,7 @@ Versions
 .. toctree::
    :maxdepth: 2
 
+   v0.25
    v0.24
    v0.23
    v0.22

docs/source/whatsnew/v0.25.rst

@@ -0,0 +1,30 @@
+Version 0.25 (Date TBD)
+=======================
+
+The new minimum supported versions of dependencies that have been updated are:
+
+* Shapely 2.0
+
+
+Features
+--------
+
+* Ruth Comer introduced `~cartopy.mpl.path.shapely_to_path` and
+  `~cartopy.mpl.path.path_to_shapely` which map a single Shapely geometry or
+  collection to a single Matplotlib path and *vice versa*. (:pull:`2455`)
+
+
+Deprecations and Removals
+-------------------------
+
+* `~cartopy.mpl.patch.path_to_geos` and `~cartopy.mpl.patch.geos_to_path` are
+  deprecated.  Use `~cartopy.mpl.path.path_to_shapely` and
+  `~cartopy.mpl.path.shapely_to_path` instead.
+
+* `~cartopy.mpl.patch.path_segments` is deprecated without replacement.  The
+  implementation is simply
+
+  .. code-block:: python
+
+    pth = path.cleaned(**kwargs)
+    return pth.vertices[:-1, :], pth.codes[:-1]

environment.yml

@@ -10,7 +10,7 @@ channels:
 dependencies:
   - cython>=0.29.28
   - numpy>=1.23
-  - shapely>=1.8
+  - shapely>=2.0
   - pyshp>=2.3
   - pyproj>=3.3.1
   - packaging>=21

lib/cartopy/crs.py

@@ -667,6 +667,7 @@ class Projection(CRS, metaclass=ABCMeta):
         'MultiPoint': '_project_multipoint',
         'MultiLineString': '_project_multiline',
         'MultiPolygon': '_project_multipolygon',
+        'GeometryCollection': '_project_geometry_collection'
     }
     # Whether or not this projection can handle wrapped coordinates
     _wrappable = False
@@ -835,7 +836,8 @@ def _project_line_string(self, geometry, src_crs):
     def _project_linear_ring(self, linear_ring, src_crs):
         """
         Project the given LinearRing from the src_crs into this CRS and
-        returns a list of LinearRings and a single MultiLineString.
+        returns a GeometryCollection containing zero or more LinearRings and
+        a single MultiLineString.
 
         """
         debug = False
@@ -915,7 +917,7 @@ def _project_linear_ring(self, linear_ring, src_crs):
         if rings:
             multi_line_string = sgeom.MultiLineString(line_strings)
 
-        return rings, multi_line_string
+        return sgeom.GeometryCollection([*rings, multi_line_string])
 
     def _project_multipoint(self, geometry, src_crs):
         geoms = []
@@ -939,6 +941,11 @@ def _project_multipolygon(self, geometry, src_crs):
                 geoms.extend(r.geoms)
         return sgeom.MultiPolygon(geoms)
 
+    def _project_geometry_collection(self, geometry, src_crs):
+        return sgeom.GeometryCollection(
+            [self.project_geometry(geom, src_crs) for geom in geometry.geoms])
+
+
     def _project_polygon(self, polygon, src_crs):
         """
         Return the projected polygon(s) derived from the given polygon.
@@ -957,7 +964,8 @@ def _project_polygon(self, polygon, src_crs):
         rings = []
         multi_lines = []
         for src_ring in [polygon.exterior] + list(polygon.interiors):
-            p_rings, p_mline = self._project_linear_ring(src_ring, src_crs)
+            geom_collection = self._project_linear_ring(src_ring, src_crs)
+            *p_rings, p_mline = geom_collection.geoms
             if p_rings:
                 rings.extend(p_rings)
             if len(p_mline.geoms) > 0:

lib/cartopy/mpl/feature_artist.py

@@ -14,12 +14,11 @@
 
 import matplotlib.artist
 import matplotlib.collections
-import matplotlib.path as mpath
 import numpy as np
 
 import cartopy.feature as cfeature
 from cartopy.mpl import _MPL_38
-import cartopy.mpl.patch as cpatch
+import cartopy.mpl.path as cpath
 
 
 class _GeomKey:
@@ -217,11 +216,7 @@ def draw(self, renderer):
                 else:
                     projected_geom = geom
 
-                geom_paths = cpatch.geos_to_path(projected_geom)
-
-                # The transform may have split the geometry into two paths, we only want
-                # one compound path.
-                geom_path = mpath.Path.make_compound_path(*geom_paths)
+                geom_path = cpath.shapely_to_path(projected_geom)
                 mapping[key] = geom_path
 
             if self._styler is None:

lib/cartopy/mpl/geoaxes.py

@@ -42,7 +42,7 @@
 import cartopy.mpl.contour
 import cartopy.mpl.feature_artist as feature_artist
 import cartopy.mpl.geocollection
-import cartopy.mpl.patch as cpatch
+import cartopy.mpl.path as cpath
 from cartopy.mpl.slippy_image_artist import SlippyImageArtist
 
 
@@ -170,26 +170,11 @@ def transform_path_non_affine(self, src_path):
         if src_path.vertices.shape == (1, 2):
             return mpath.Path(self.transform(src_path.vertices))
 
-        transformed_geoms = []
-        geoms = cpatch.path_to_geos(src_path)
+        geom = cpath.path_to_shapely(src_path)
+        transformed_geom = self.target_projection.project_geometry(
+            geom, self.source_projection)
 
-        for geom in geoms:
-            proj_geom = self.target_projection.project_geometry(
-                geom, self.source_projection)
-            transformed_geoms.append(proj_geom)
-
-        if not transformed_geoms:
-            result = mpath.Path(np.empty([0, 2]))
-        else:
-            paths = cpatch.geos_to_path(transformed_geoms)
-            if not paths:
-                return mpath.Path(np.empty([0, 2]))
-            points, codes = list(zip(*[cpatch.path_segments(path,
-                                                            curves=False,
-                                                            simplify=False)
-                                       for path in paths]))
-            result = mpath.Path(np.concatenate(points, 0),
-                                np.concatenate(codes))
+        result = cpath.shapely_to_path(transformed_geom)
 
         # store the result in the cache for future performance boosts
         key = (self.source_projection, self.target_projection)
@@ -228,14 +213,14 @@ def set_transform(self, transform):
         super().set_transform(self._trans_wrap)
 
     def set_boundary(self, path, transform):
-        self._original_path = cpatch._ensure_path_closed(path)
+        self._original_path = cpath._ensure_path_closed(path)
         self.set_transform(transform)
         self.stale = True
 
     def _adjust_location(self):
         if self.stale:
             self.set_path(
-                cpatch._ensure_path_closed(
+                cpath._ensure_path_closed(
                     self._original_path.clip_to_bbox(self.axes.viewLim)))
             # Some places in matplotlib's transform stack cache the actual
             # path so we trigger an update by invalidating the transform.
@@ -255,13 +240,13 @@ def __init__(self, axes, **kwargs):
 
     def set_boundary(self, path, transform):
         # Make sure path is closed (required by "Path.clip_to_bbox")
-        self._original_path = cpatch._ensure_path_closed(path)
+        self._original_path = cpath._ensure_path_closed(path)
         self.set_transform(transform)
         self.stale = True
 
     def _adjust_location(self):
         if self.stale:
-            self._path = cpatch._ensure_path_closed(
+            self._path = cpath._ensure_path_closed(
                 self._original_path.clip_to_bbox(self.axes.viewLim)
                 )
 
@@ -1535,7 +1520,7 @@ def _boundary(self):
         The :data:`.patch` and :data:`.spines['geo']` are updated to match.
 
         """
-        path, = cpatch.geos_to_path(self.projection.boundary)
+        path = cpath.shapely_to_path(self.projection.boundary)
 
         # Get the outline path in terms of self.transData
         proj_to_data = self.projection._as_mpl_transform(self) - self.transData

lib/cartopy/mpl/patch.py

@@ -12,54 +12,30 @@
 
 """
 
+import warnings
+
 from matplotlib.path import Path
 import numpy as np
 import shapely.geometry as sgeom
 
+import cartopy.mpl.path as cpath
 
-def _ensure_path_closed(path):
-    """
-    Method to ensure that a path contains only closed sub-paths.
-
-    Parameters
-    ----------
-    path
-        A :class:`matplotlib.path.Path` instance.
-
-    Returns
-    -------
-    path
-        A :class:`matplotlib.path.Path` instance with only closed polygons.
-
-    """
-    # Split path into potential sub-paths and close all polygons
-    # (explicitly disable path simplification applied in to_polygons)
-    should_simplify = path.should_simplify
-    try:
-        path.should_simplify = False
-        polygons = path.to_polygons()
-    finally:
-        path.should_simplify = should_simplify
-
-    codes, vertices = [], []
-    for poly in polygons:
-        vertices.extend([poly[0], *poly])
-        codes.extend([Path.MOVETO, *[Path.LINETO]*(len(poly) - 1), Path.CLOSEPOLY])
-
-    return Path(vertices, codes)
 
 def geos_to_path(shape):
     """
     Create a list of :class:`matplotlib.path.Path` objects that describe
     a shape.
 
+    .. deprecated:: 0.25
+       Use `cartopy.mpl.path.shapely_to_path` instead.
+
     Parameters
     ----------
     shape
         A list, tuple or single instance of any of the following
         types: :class:`shapely.geometry.point.Point`,
         :class:`shapely.geometry.linestring.LineString`,
-        :class:`shapely.geometry.linestring.LinearRing`,
+        :class:`shapely.geometry.polygon.LinearRing`,
         :class:`shapely.geometry.polygon.Polygon`,
         :class:`shapely.geometry.multipoint.MultiPoint`,
         :class:`shapely.geometry.multipolygon.MultiPolygon`,
@@ -73,6 +49,9 @@ def geos_to_path(shape):
         A list of :class:`matplotlib.path.Path` objects.
 
     """
+    warnings.warn("geos_to_path is deprecated and will be removed in a future release."
+                  "  Use cartopy.mpl.path.shapely_to_path instead.",
+                  DeprecationWarning, stacklevel=2)
     if isinstance(shape, (list, tuple)):
         paths = []
         for shp in shape:
@@ -115,6 +94,8 @@ def path_segments(path, **kwargs):
     Create an array of vertices and a corresponding array of codes from a
     :class:`matplotlib.path.Path`.
 
+    .. deprecated:: 0.25
+
     Parameters
     ----------
     path
@@ -123,7 +104,7 @@ def path_segments(path, **kwargs):
     Other Parameters
     ----------------
     kwargs
-        See :func:`matplotlib.path.iter_segments` for details of the keyword
+        See `matplotlib.path.Path.iter_segments` for details of the keyword
         arguments.
 
     Returns
@@ -135,15 +116,20 @@ def path_segments(path, **kwargs):
         codes and their meanings.
 
     """
-    pth = path.cleaned(**kwargs)
-    return pth.vertices[:-1, :], pth.codes[:-1]
+    warnings.warn(
+        "path_segments is deprecated and will be removed in a future release.",
+        DeprecationWarning, stacklevel=2)
+    return cpath._path_segments(path, **kwargs)
 
 
 def path_to_geos(path, force_ccw=False):
     """
     Create a list of Shapely geometric objects from a
     :class:`matplotlib.path.Path`.
 
+    .. deprecated:: 0.25
+       Use `cartopy.mpl.path.path_to_shapely` instead.
+
     Parameters
     ----------
     path
@@ -163,8 +149,11 @@ def path_to_geos(path, force_ccw=False):
         :class:`shapely.geometry.multilinestring.MultiLineString`.
 
     """
+    warnings.warn("path_to_geos is deprecated and will be removed in a future release."
+                  "  Use cartopy.mpl.path.path_to_shapely instead.",
+                  DeprecationWarning, stacklevel=2)
     # Convert path into numpy array of vertices (and associated codes)
-    path_verts, path_codes = path_segments(path, curves=False)
+    path_verts, path_codes = cpath._path_segments(path, curves=False)
 
     # Split into subarrays such that each subarray consists of connected
     # line segments based on the start of each one being marked by a