Add CRS accessor for cartopy

[aulemahal]

This fixes #395.

• [x] Add a crs accessor to both the Dataset and DataArray accessors.
  • For Datasets, we can only return something if only one grid mapping is defined.
  • For DataArrays, the grid mapping must be a coordinate
  • For both, the absence of a grid mapping and the presence of a "longitude" entry in the dataset is assumed to mean that the grid mapping is "latitude_longitude". CF allows an absent grid mapping for this CRS.
• [ ] Include this in the plotting function (?)
• [x] Add tests

[dcherian]

Also, +1 to automatic use in the plotting function, that was my initial intent!

[aulemahal]

Ugh... So it's not as easy as I thought. ccrs.Projection(pyproj.CRS.from_cf(CF)) is not equivalent to the native cartopy projection because it seems to be missing information about the limits of the projection domain. If I try to set this PR's cartopy_crs as the "projection" of a matplotlib figure and also use it as the transform of the data, I get this error:

---------------------------------------------------------------------------
NotImplementedError                       Traceback (most recent call last)
Cell In[23], line 1
----> 1 ds.Temperature_surface.isel(time=0).cf.plot()

File ~/Projets/cf-xarray/cf_xarray/accessor.py:1154, in __call__(self, *args, **kwargs)
   1145 def __call__(self, *args, **kwargs):
   1146     """
   1147     Allows .plot()
   1148     """
   1149     plot = _getattr(
   1150         obj=self._obj,
   1151         attr="plot",
   1152         accessor=self.accessor,
   1153         key_mappers=dict.fromkeys(self._keys, (_single(_get_all),)),
-> 1154     )
   1155     return self._plot_decorator(plot)(*args, **kwargs)

File <string>:40, in _plot_wrapper(*args, **kwargs)

File ~/Projets/cf-xarray/cf_xarray/accessor.py:763, in _getattr.<locals>.wrapper(*args, **kwargs)
    759 posargs, arguments = accessor._process_signature(
    760     func, args, kwargs, key_mappers
    761 )
    762 final_func = extra_decorator(func) if extra_decorator else func
--> 763 result = final_func(*posargs, **arguments)
    764 if wrap_classes and isinstance(result, _WRAPPED_CLASSES):
    765     result = _CFWrappedClass(result, accessor)

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/xarray/plot/accessor.py:48, in DataArrayPlotAccessor.__call__(self, **kwargs)
     46 @functools.wraps(dataarray_plot.plot, assigned=("__doc__", "__annotations__"))
     47 def __call__(self, **kwargs) -> Any:
---> 48     return dataarray_plot.plot(self._da, **kwargs)

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/xarray/plot/dataarray_plot.py:310, in plot(darray, row, col, col_wrap, ax, hue, subplot_kws, **kwargs)
    306     plotfunc = hist
    308 kwargs["ax"] = ax
--> 310 return plotfunc(darray, **kwargs)

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/xarray/plot/dataarray_plot.py:1604, in _plot2d.<locals>.newplotfunc(***failed resolving arguments***)
   1600 if "imshow" == plotfunc.__name__ and isinstance(aspect, str):
   1601     # forbid usage of mpl strings
   1602     raise ValueError("plt.imshow's `aspect` kwarg is not available in xarray")
-> 1604 ax = get_axis(figsize, size, aspect, ax, **subplot_kws)
   1606 primitive = plotfunc(
   1607     xplt,
   1608     yplt,
   (...)   1615     **kwargs,
   1616 )
   1618 # Label the plot with metadata

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/xarray/plot/utils.py:497, in get_axis(figsize, size, aspect, ax, **subplot_kws)
    494     raise ValueError("cannot use subplot_kws with existing ax")
    496 if ax is None:
--> 497     ax = _maybe_gca(**subplot_kws)
    499 return ax

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/xarray/plot/utils.py:513, in _maybe_gca(**subplot_kws)
    509 if f.axes:
    510     # can not pass kwargs to active axes
    511     return plt.gca()
--> 513 return plt.axes(**subplot_kws)

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/matplotlib/pyplot.py:1353, in axes(arg, **kwargs)
   1351 if arg is None:
   1352     if pos is None:
-> 1353         return fig.add_subplot(**kwargs)
   1354     else:
   1355         return fig.add_axes(pos, **kwargs)

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/matplotlib/figure.py:768, in FigureBase.add_subplot(self, *args, **kwargs)
    766         args = tuple(map(int, str(args[0])))
    767     projection_class, pkw = self._process_projection_requirements(**kwargs)
--> 768     ax = projection_class(self, *args, **pkw)
    769     key = (projection_class, pkw)
    770 return self._add_axes_internal(ax, key)

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/cartopy/mpl/geoaxes.py:407, in GeoAxes.__init__(self, *args, **kwargs)
    403     raise ValueError("A GeoAxes can only be created with a "
    404                      "projection of type cartopy.crs.Projection")
    405 self.projection = projection
--> 407 super().__init__(*args, **kwargs)
    408 self.img_factories = []
    409 self._done_img_factory = False

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/matplotlib/axes/_base.py:732, in _AxesBase.__init__(self, fig, facecolor, frameon, sharex, sharey, label, xscale, yscale, box_aspect, forward_navigation_events, *args, **kwargs)
    729 self.set_axisbelow(mpl.rcParams['axes.axisbelow'])
    731 self._rasterization_zorder = None
--> 732 self.clear()
    734 # funcs used to format x and y - fall back on major formatters
    735 self.fmt_xdata = None

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/cartopy/mpl/geoaxes.py:587, in GeoAxes.clear(self)
    585 """Clear the current Axes and add boundary lines."""
    586 result = super().clear()
--> 587 self.__clear()
    588 return result

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/cartopy/mpl/geoaxes.py:575, in GeoAxes.__clear(self)
    572 self._tight = True
    573 self.set_aspect('equal')
--> 575 self._boundary()
    577 # XXX consider a margin - but only when the map is not global...
    578 # self._xmargin = 0.15
    579 # self._ymargin = 0.15
    581 self.dataLim.intervalx = self.projection.x_limits

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/cartopy/mpl/geoaxes.py:1538, in GeoAxes._boundary(self)
   1531 def _boundary(self):
   1532     """
   1533     Add the map's boundary to this GeoAxes.
   1534 
   1535     The :data:`.patch` and :data:`.spines['geo']` are updated to match.
   1536 
   1537     """
-> 1538     path, = cpatch.geos_to_path(self.projection.boundary)
   1540     # Get the outline path in terms of self.transData
   1541     proj_to_data = self.projection._as_mpl_transform(self) - self.transData

File ~/miniforge3/envs/cfxr-dev/lib/python3.13/site-packages/cartopy/crs.py:701, in Projection.boundary(self)
    698 @property
    699 def boundary(self):
    700     if self.bounds is None:
--> 701         raise NotImplementedError
    702     x0, x1, y0, y1 = self.bounds
    703     return sgeom.LineString([(x0, y0), (x0, y1), (x1, y1), (x1, y0),
    704                              (x0, y0)])

NotImplementedError: 

If instead I use "PlateCarree" as the projection for the figure, my rotated pole test data plots like this:

Here is the same plot using cartopy's RotatedPole projection instead of the "trick" implemented in this PR:

Thus, it seems my hopes are voided, this trick is not sufficient. Maybe it's only a problem for Rotated Pole, but that's my bread and butter. I have some code to re-implement from_cf directly in cartopy, I will try to push that instead. We'll see.

[dcherian]

Can we set the bounds on the cartopy crs since we have all the info?