pydata / xarray

N-D labeled arrays and datasets in Python
https://xarray.dev
Apache License 2.0
3.62k stars 1.08k forks source link

Plotting inconsistencies with Cartopy #3169

Closed andreas-h closed 2 years ago

andreas-h commented 5 years ago

I'm on xarray 0.11.3, and currently do not have the possibility to test on master. However, the commit history suggests that this issue still exists. Please accept my apologies if I missed something.

I'm following the plotting documentation at https://xarray.pydata.org/en/stable/plotting.html#maps

When plotting two maps in the same call (exactly as written in the documentation),

p = air.isel(time=[0, 1]).plot(transform=ccrs.PlateCarree(), col='time',
                               subplot_kws={'projection': ccrs.Orthographic(-80, 35)})

for ax in p.axes.flat:
    ax.coastlines()
    ax.gridlines()

everything works fine.

However, when I only want to plot one map,

p = air.isel(time=[0]).plot(transform=ccrs.PlateCarree(),
                            subplot_kws={'projection': ccrs.Orthographic(-80, 35)})

for ax in p.axes.flat:
    ax.coastlines()
    ax.gridlines()

I get the following error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-63-2c5a2a533878> in <module>
      1 p = air.isel(time=[0]).plot(transform=ccrs.PlateCarree(),
----> 2                             subplot_kws={'projection': ccrs.Orthographic(-80, 35)})
      3 
      4 for ax in p.axes.flat:
      5     ax.coastlines()

/home/jupyterhub/conda/envs/prod/lib/python3.6/site-packages/xarray/plot/plot.py in __call__(self, **kwargs)
    550 
    551     def __call__(self, **kwargs):
--> 552         return plot(self._da, **kwargs)
    553 
    554     @functools.wraps(hist)

/home/jupyterhub/conda/envs/prod/lib/python3.6/site-packages/xarray/plot/plot.py in plot(darray, row, col, col_wrap, ax, hue, rtol, subplot_kws, **kwargs)
    185     kwargs['ax'] = ax
    186 
--> 187     return plotfunc(darray, **kwargs)
    188 
    189 

/home/jupyterhub/conda/envs/prod/lib/python3.6/site-packages/xarray/plot/plot.py in newplotfunc(darray, x, y, figsize, size, aspect, ax, row, col, col_wrap, xincrease, yincrease, add_colorbar, add_labels, vmin, vmax, cmap, center, robust, extend, levels, infer_intervals, colors, subplot_kws, cbar_ax, cbar_kwargs, xscale, yscale, xticks, yticks, xlim, ylim, norm, **kwargs)
    852                              vmax=cmap_params['vmax'],
    853                              norm=cmap_params['norm'],
--> 854                              **kwargs)
    855 
    856         # Label the plot with metadata

/home/jupyterhub/conda/envs/prod/lib/python3.6/site-packages/xarray/plot/plot.py in pcolormesh(x, y, z, ax, infer_intervals, **kwargs)
   1106             y = _infer_interval_breaks(y, axis=0)
   1107 
-> 1108     primitive = ax.pcolormesh(x, y, z, **kwargs)
   1109 
   1110     # by default, pcolormesh picks "round" values for bounds

/home/jupyterhub/conda/envs/prod/lib/python3.6/site-packages/matplotlib/__init__.py in inner(ax, data, *args, **kwargs)
   1808                         "the Matplotlib list!)" % (label_namer, func.__name__),
   1809                         RuntimeWarning, stacklevel=2)
-> 1810             return func(ax, *args, **kwargs)
   1811 
   1812         inner.__doc__ = _add_data_doc(inner.__doc__,

/home/jupyterhub/conda/envs/prod/lib/python3.6/site-packages/matplotlib/axes/_axes.py in pcolormesh(self, alpha, norm, cmap, vmin, vmax, shading, antialiased, *args, **kwargs)
   6011         if (not isinstance(t, mtransforms.Transform) and
   6012             hasattr(t, '_as_mpl_transform')):
-> 6013             t = t._as_mpl_transform(self.axes)
   6014 
   6015         if t and any(t.contains_branch_seperately(self.transData)):

lib/cartopy/_crs.pyx in cartopy._crs.CRS._as_mpl_transform()

ValueError: Axes should be an instance of GeoAxes, got <class 'matplotlib.axes._subplots.AxesSubplot'>

It is not clear (from the user perspective) why one should work and the other should not.

Of course I can create the axes instance before using xarray's plot function, but it would be nicer if I wouldn't have to.

shoyer commented 5 years ago

I agree, this is unfortunate. I suspect the source of the discrepancy is that we use a totally different codepath if either col or row is supplied, creating our own subplot grid.

ahuang11 commented 4 years ago

I think this is fixed in the latest master.

import xarray as xr
import numpy as np
import cartopy.crs as ccrs

da = xr.tutorial.open_dataset('air_temperature')['air']

p = da.isel(time=[0, 1]).plot(
    transform=ccrs.PlateCarree(), col='time',
    subplot_kws={'projection': ccrs.Orthographic(-80, 35)}
)

for ax in p.axes.flat:
    ax.coastlines()
    ax.gridlines()
p = da.isel(time=0).plot(
    transform=ccrs.PlateCarree(),
    subplot_kws={'projection': ccrs.Orthographic(-80, 35)}
)

download download (1)

However if you try iterating over the axes, it will crash because p is quadmesh type not a facetgrid type so I was wondering whether it'd be good to have an keyword like force_facetgrid_type=True which uses the FacetGrid class to plot a single map so that the user can keep the same code for single axes and facet axes?

p = da.isel(time=0).plot(
    transform=ccrs.PlateCarree(),
    subplot_kws={'projection': ccrs.Orthographic(-80, 35)}
)

for ax in p.axes.flat:
    ax.coastlines()
    ax.gridlines()
mathause commented 4 years ago

As a workaround you can use:

p = da.isel(time=[0]).plot.pcolormesh(
    col="time",
    transform=ccrs.PlateCarree(),
    subplot_kws={'projection': ccrs.Orthographic(-80, 35)}
)

Note the time=[0] - using time=0 squezes the time dimension. Also note the .plot.pcolormesh. Using just .plot also calls daarray.squeeze().

https://github.com/pydata/xarray/blob/8fab5a2449d8368251f96fc2b9d1eaa3040894e6/xarray/plot/plot.py#L163

stale[bot] commented 2 years ago

In order to maintain a list of currently relevant issues, we mark issues as stale after a period of inactivity

If this issue remains relevant, please comment here or remove the stale label; otherwise it will be marked as closed automatically