People often complain about how slow Cartopy is lately since we moved over to pyproj and not the direct C bindings we were using. If a user knows what they are doing and wants something fast or for a small domain with no wrapped coordinates, we can help them out and give a much faster option, so perhaps we should do that. Below is a quick proposal for how this could work.
Feedback on this would be most welcome! This is just a thought at this point, but perhaps one that others have better ideas for.
Proposal
Add a fast-path to CRS.project_geometry() that immediately transforms all of the points in the geometry to the destination projection. This is what basemap does and it is much faster than handling interpolations/bisections between points to make smoother lines. It is equivalent to the fast_transform argument that we added to the contour() functions.
There are a few ways we could implement this and I'm not sure which would be best / most useful.
- A global configuration attribute somewhere:
cartopy.interpolate_paths = False.
Benefit: easy to change behavior
Downside: can't individually control which paths/projections are affected by this.
- Add more
fast_transform keyword decorators to plotting functions.
Benefit: we already have this written
Downside: would require some more work on which x/y arguments to transform for each plotting function, and would also probably require some work on the add_patch() etc... areas of the code that aren't directly plotting functions.
- Add a new keyword / property to CRS and subclasses that could control this on a per-CRS basis.
Benefit: user could decide which projection and artists to interpolate or not (i.e. you could have both
pc_slow = PlateCarree(interpolate_paths=True) andpc_fast = PlateCarree(interpolate_paths=False)`)
Downside: Would require more work on keeping CRS's consistent and subclassing of them would add another somewhat arbitrary keyword that isn't CRS-specific, bur rather Cartopy-specific.
Issues with this approach
There are some downsides to this associated with long paths and wrapped coordinates. For example, the "logo" example:
current (with interpolated paths): https://scitools.org.uk/cartopy/docs/latest/gallery/miscellanea/logo.html
proposed (no interpolation along paths):
Code to reproduce
Example of some minimal code that could be adapted/expanded upon.
diff --git a/lib/cartopy/crs.py b/lib/cartopy/crs.py
index 30c4170e..a77dbb4f 100644
--- a/lib/cartopy/crs.py
+++ b/lib/cartopy/crs.py
@@ -19,6 +19,7 @@ import math
import warnings
import numpy as np
+import shapely.ops
import shapely.geometry as sgeom
from pyproj import Transformer
from pyproj.exceptions import ProjError
@@ -123,7 +124,7 @@ class CRS(_CRS):
#: Whether this projection can handle ellipses.
_handles_ellipses = True
- def __init__(self, proj4_params, globe=None):
+ def __init__(self, proj4_params, globe=None, interpolate_paths=True):
"""
Parameters
----------
@@ -137,7 +138,11 @@ class CRS(_CRS):
globe: :class:`~cartopy.crs.Globe` instance, optional
If omitted, the default Globe instance will be created.
See :class:`~cartopy.crs.Globe` for details.
-
+ interpolate_paths: bool, default True
+ Whether to interpolate paths in this CRS. Setting this
+ to False can significantly improve performance, but
+ paths will not wrap and may look jagged if there
+ are not enough points in the initial geometry.
"""
# for compatibility with pyproj.CRS and rasterio.crs.CRS
try:
@@ -187,6 +192,7 @@ class CRS(_CRS):
init_items.append(f'+{k}')
self.proj4_init = ' '.join(init_items) + ' +no_defs'
super().__init__(self.proj4_init)
+ self.interpolate_paths = interpolate_paths
def __eq__(self, other):
if isinstance(other, CRS) and self.proj4_init == other.proj4_init:
@@ -801,6 +807,11 @@ class Projection(CRS, metaclass=ABCMeta):
elif not isinstance(src_crs, CRS):
raise TypeError('Source CRS must be an instance of CRS'
' or one of its subclasses, or None.')
+
+ if not getattr(src_crs, 'interpolate_paths', True):
+ pyproj_trans = _get_transformer_from_crs(src_crs, self).transform
+ return shapely.ops.transform(pyproj_trans, geometry)
+
geom_type = geometry.geom_type
method_name = self._method_map.get(geom_type)
if not method_name:
Type: Question
1 Jul 06, 2022
1 Jan 23, 2022
2.4k Jan 07, 2023
4k Jan 04, 2023
0 Sep 12, 2021
146 Nov 12, 2022
138 Dec 06, 2022
2 Dec 13, 2021
1 Jan 10, 2022
990 Jan 05, 2023
15 Dec 28, 2022
820 Jan 06, 2023
7 Oct 27, 2021
9 Mar 18, 2022
75 Dec 15, 2022
7 Apr 05, 2022
366 Oct 27, 2022
338 Nov 25, 2022
1.3k Dec 13, 2022
0 May 27, 2022