Set cartopy crs automatically

[larsbuntemeyer]

Woud be nice to have a cf function that generates CRS from grid_mappings that cartopy can understand for plotting. Would be nice to simply work without having to know which mapping is present, e.g., no need to pick the required cartopy CRS by hand.

Continue the dicussions started here:

• 2D bounds - simple version #370
• add cf.grid_mapping_names #391 (comment)

Related issues:

• from_cf constructor inherited from pyproj is broken. SciTools/cartopy#2099

[larsbuntemeyer]

For the cartopy case, this seems to have been an issue for some time now.

• Create Cartopy projection from pyproj.Proj SciTools/cartopy#813

I am unsure, why this seems to be complicated (also it should not?)...

A possible workaround before going upstream: we could use a couple of mappings that simply map cf grid_mapping_name to the appropriate cartopy CRS (including **kwargs).

[trexfeathers]

Those of us developing Iris are equally frustrated around this topic:

• Adopt pyproj.CRS SciTools/iris#4643
• Generate ccrs.Globe with datum specifiers other than the proj4 string datum ID SciTools/cartopy#2039

Apologies for my incomplete understanding of this area, I believe the two links I have provided are relevant, although there may be others.

[larsbuntemeyer]

@dcherian I have some code that i now use regularly to do transformations based on grid mappings using pyproj:

https://github.com/euro-cordex/py-cordex/blob/90e18d1bf020b3b869d9d2ffcc7039a63a92f176/cordex/transform.py#L77-L175

I also use that to compute boundaries by transforming the linear x, y boundaries according to the grid mapping info. I was wondering if that could have a place here in a kind of mapping class, e.g., as an option in _guess_bounds... also pinging @aulemahal here...

[dcherian]

I think the transformations should be in rioxarray already?

cc @snowman2

[larsbuntemeyer]

Ahh, nice! Yes, you mentioned that in the roadmap. I'll check it out.

[snowman2]

https://pyproj4.github.io/pyproj/stable/build_crs_cf.html

[snowman2]

https://corteva.github.io/rioxarray/stable/getting_started/crs_management.html

[dcherian]

FWIW I'd be happy to add .cf.crs -> pyproj.CRS and .cf.crs_cartopy if someone has heuristics for cartopy_crs lying around.

[larsbuntemeyer]

I started something here for ccrs some time ago (it's the naive approach). Would be happy to move something to cf_xarray!

[aulemahal]

Pyproj already has almost everything here https://github.com/pyproj4/pyproj/blob/69fb46f2c007f71a3534a7c935b6f8813a45dfa6/pyproj/crs/_cf1x8.py

When I last looked at this, fixing SciTools/cartopy#2099 really looked like the most elegant way forward... Also, if a naive and explicit mapping from CF to Cartopy is the solution, maybe it would be most useful in cartopy directly ?

[dcherian]

I think this is mostly solved with pyproj.CRS.from_cf

[aulemahal]

Is it solved ? Testing with the latest pyproj and cartopy still seems to have the issue, no ? Is there a new way to do this ?

import pyproj
import cartopy.crs as ccrs
import matplotlib.pyplot as plt
import xarray as xr

ds = xr.open_dataset('https://pavics.ouranos.ca/twitcher/ows/proxy/thredds/dodsC/datasets/simulations/RCM-CMIP6/CORDEX/NAM-12/1hr/NAM-12_CNRM-ESM2-1_historical_r1i1p1f2_OURANOS_CRCM5_v1-r1_1hr_195001010000-201412312300.ncml')

# pyproj only 
crs = pyproj.CRS.from_cf(ds.crs.attrs)
fig, ax = plt.subplots(figsize=(12, 9), subplot_kw={'projection': ccrs.PlateCarree()})
ds.tas.isel(time=0).plot(ax=ax, transform=crs)

# Or cartopy ?
crs = ccrs.CRS.from_cf(ds.crs.attrs)
ds.tas.isel(time=0).plot(ax=ax, transform=RP)

Both raise AttributeError: 'DerivedGeographicCRS' object has no attribute 'contains_branch_seperately'

By which I mean that even if a pyproj CRS is obtainable through the cf attrs, it seems the main common use of such an object is still not possible.