How to efficiently use DataArrays with Cartopy's add_cyclic_point utility?

[ktyle]

Hi, I'm not sure as to the best way to take a global NetCDF dataset (e.g. 2.5 x 2.5 deg GFS) whose longitude range does not "cycle" back to 0/360 and pass the variable of interest and the longitude array to Cartopy's add_cyclic_point() utility. I am first selecting the variable of interest and the longitude and converting them to arrays, and then passing them along to add_cyclic_point, but not sure then how to take the tuple that is returned by that utility and "reintegrating" it into xarray's plot.pcolormesh utility. Do I need to create a new Xarray DataArray from the output of add_cyclic_point?

[shoyer]

Do I need to create a new Xarray DataArray from the output of add_cyclic_point?

Yes, that would be my suggested approach.

This would be a good use case for xarray.apply (#964).

[shoyer]

Note that xarray.apply doesn't solve this yet, but hopefully it will soon. I'll update this issue when it does.

[lvankampenhout]

Any update on this issue? It would be great if add_cyclic_point could be applied to all variables automatically.

Just for other peoples reference, I now have this workaround, creating erai_jja_cy a 'cyclic' version of erai_jja:

dd, ll = add_cyclic_point(erai_jja.values, erai_jja.lon)
erai_jja_cy = xr.DataArray(dd, coords={'lat':erai_jja.lat, 'lon':ll}, dims=('lat','lon'))

[shoyer]

Something like this might work for you:

import xarray
import numpy as np

def add_cyclic_point(xarray_obj, dim, period=None):
  if period is None:
    period = xarray_obj.sizes[dim] / xarray_obj.coords[dim][:2].diff(dim).item()
  first_point = xarray_obj.isel({dim: slice(1)})
  first_point.coords[dim] = first_point.coords[dim]+period
  return xarray.concat([xarray_obj, first_point], dim=dim)

Example usage:

>>> ds = xarray.Dataset({'foo': ('x', np.arange(10))}, {'x': np.arange(10)})
>>> add_cyclic_point(ds, dim='x', period=10)
<xarray.Dataset>
Dimensions:  (x: 11)
Coordinates:
  * x        (x) int64 0 1 2 3 4 5 6 7 8 9 10
Data variables:
    foo      (x) int64 0 1 2 3 4 5 6 7 8 9 0

It would also make sense to add a more generic .pad() method, which this function would be a special case of (if mode='wrap').

[lvankampenhout]

Stephan, thanks a lot for your code snippet from December, this is an elegant solution to the problem. One minor correction though, because I found that it fails to infer the period if none is given. The divide should be a multiplication I believe, i.e.

import xarray
import numpy as np

def add_cyclic_point(xarray_obj, dim, period=None):
  if period is None:
    period = xarray_obj.sizes[dim] * xarray_obj.coords[dim][:2].diff(dim).item()
  first_point = xarray_obj.isel({dim: slice(1)})
  first_point.coords[dim] = first_point.coords[dim]+period
  return xarray.concat([xarray_obj, first_point], dim=dim)