# Unexpected behaviour for a CF conforming datasets

**URL:** <https://discourse.holoviz.org/t/unexpected-behaviour-for-a-cf-conforming-datasets/5851>\
**Category:** hvPlot\
**Tags:** holoviews\
**Created:** [August 11, 2023, 10:15pm UTC](https://discourse.holoviz.org/t/unexpected-behaviour-for-a-cf-conforming-datasets/5851 "2023-08-11T22:15:44Z")\
**Posts on this page:** 1\
**Page:** 1

<div class="post-metadata">

**Author:** ![yt87](https://discourse.holoviz.org/user_avatar/discourse.holoviz.org/yt87/32/5021_2.png) [@yt87](https://discourse.holoviz.org/u/yt87)\
**Post date:** [August 11, 2023, 10:15pm UTC](https://discourse.holoviz.org/t/unexpected-behaviour-for-a-cf-conforming-datasets/5851/1 "2023-08-11T22:15:44Z")

</div>

I am trying to plot a CF conforming dataset. A simplified definition:

```auto
import os
from pathlib import Path

import cartopy.crs as ccrs
import cf_xarray
import hvplot.xarray
import numpy as np
import pyproj

import xarray as xr

proj_cf = {
    "semi_major_axis": 6378137.0,
    "semi_minor_axis": 6356752.314245179,
    "inverse_flattening": 298.257223563,
    "reference_ellipsoid_name": "WGS 84",
    "longitude_of_prime_meridian": 0.0,
    "prime_meridian_name": "Greenwich",
    "geographic_crs_name": "unnamed",
    "horizontal_datum_name": "Unknown based on WGS84 ellipsoid",
    "grid_mapping_name": "rotated_latitude_longitude",
    "grid_north_pole_latitude": 23.5,
    "grid_north_pole_longitude": 15.0,
    "north_pole_grid_longitude": 0.0,
}

crs = pyproj.CRS.from_cf(proj_cf)
proj_cf["crs_wkt"] = crs.to_wkt()

xlim = (-2, 2)
ylim = (-2, 2)
eps = 1e-5
rlons = np.arange(xlim[0] + 1 / 2, xlim[1] + eps)
rlats = np.arange(ylim[0] + 1 / 2, ylim[1] + eps)

ds = xr.Dataset(
    {
        "dummy": (
            ("rlats", "rlons"),
            np.ones((rlats.size, rlons.size)),
            {"grid_mapping": "rotated_latitude_longitude"},
        )
    },
    coords={
        "rlon": (
            ("rlon"),
            rlons,
            {
                "standard_name": "grid_longitude",
                "units": "degrees",
                "long_name": "longitude in rotated pole grid",
                "axis": "X",
            },
        ),
        "rlat": (
            ("rlat"),
            rlats,
            {
                "standard_name": "grid_latitude",
                "units": "degrees",
                "long_name": "latitude in rotated pole grid",
                "axis": "Y",
            },
        ),
        "rotated_latitude_longitude": ((), 0, proj_cf), # CF convention
    },
)

proj_attrs = ds.coords["rotated_latitude_longitude"].attrs
crs_cartopy = ccrs.RotatedPole(
    pole_longitude=proj_attrs["grid_north_pole_longitude"],
    pole_latitude=proj_attrs["grid_north_pole_latitude"],
    central_rotated_longitude=proj_attrs["north_pole_grid_longitude"],
    globe=ccrs.Globe(ellipse="WGS84"),
)

```

The first attempt to plot it:

```auto
ds["dummy"].hvplot.quadmesh(projection=crs_cartopy)

```

Getting a traceback

```auto
---------------------------------------------------------------------------
DataError Traceback (most recent call last)
Cell In[2], line 2
      1 'trying ds'
----> 2 ds["dummy"].hvplot.quadmesh(projection=crs_cartopy)

File ~/mambaforge/envs/icec/lib/python3.10/site-packages/hvplot/plotting/core.py:2067, in hvPlot.quadmesh(self, x, y, z, colorbar, **kwds)
   2018 def quadmesh(self, x=None, y=None, z=None, colorbar=True, **kwds):
   2019 """
   2020 QuadMesh plot
   2021 
   (...)
   2065 - HoloViews: https://holoviews.org/reference/elements/bokeh/QuadMesh.html
   2066 """
-> 2067 return self(x, y, z=z, kind="quadmesh", colorbar=colorbar, **kwds)

File ~/mambaforge/envs/icec/lib/python3.10/site-packages/hvplot/plotting/core.py:92, in hvPlotBase. __call__ (self, x, y, kind, **kwds)
     89 plot = self._get_converter(x, y, kind, **kwds)(kind, x, y)
     90 return pn.panel(plot, **panel_dict)
---> 92 return self._get_converter(x, y, kind, **kwds)(kind, x, y)

File ~/mambaforge/envs/icec/lib/python3.10/site-packages/hvplot/converter.py:1241, in HoloViewsConverter. __call__ (self, kind, x, y)
   1239 else:
   1240 name = data.name or self.label or self.value_label
-> 1241 dataset = Dataset(data, self.indexes, name)
   1242 else:
   1243 try:

File ~/mambaforge/envs/icec/lib/python3.10/site-packages/holoviews/core/data/ __init__.py:331, in Dataset. __init__ (self, data, kdims, vdims, **kwargs)
    328 kdims, vdims = kwargs.get('kdims'), kwargs.get('vdims')
    330 validate_vdims = kwargs.pop('_validate_vdims', True)
--> 331 initialized = Interface.initialize(type(self), data, kdims, vdims,
    332 datatype=kwargs.get('datatype'))
    333 (data, self.interface, dims, extra_kws) = initialized
    334 super(). __init__ (data, **dict(kwargs,**dict(dims, **extra_kws)))

File ~/mambaforge/envs/icec/lib/python3.10/site-packages/holoviews/core/data/interface.py:253, in Interface.initialize(cls, eltype, data, kdims, vdims, datatype)
    251 continue
    252 try:
--> 253 (data, dims, extra_kws) = interface.init(eltype, data, kdims, vdims)
    254 break
    255 except DataError:

File ~/mambaforge/envs/icec/lib/python3.10/site-packages/holoviews/core/data/xarray.py:207, in XArrayInterface.init(cls, eltype, data, kdims, vdims)
    205 raise TypeError('Data must be be an xarray Dataset type.')
    206 elif not_found:
--> 207 raise DataError("xarray Dataset must define coordinates "
    208 "for all defined kdims, %s coordinates not found."
    209 % not_found, cls)
    211 for vdim in vdims:
    212 if packed:

DataError: xarray Dataset must define coordinates for all defined kdims, [Dimension('rlons'), Dimension('rlats')] coordinates not found.

XArrayInterface expects gridded data, for more information on supported datatypes see http://holoviews.org/user_guide/Gridded_Datasets.html

```

So, it looks that holoviews can not get projection from the dataset directly.  
Second attempt: adding `crs=crs_cartopy` to the arguments.

```auto
ds["dummy"].hvplot.quadmesh(crs=crs_cartopy, projection=crs_cartopy)

```

Surprisingly, the traceback is the same.

Third attempt: removing projection definition in the dataset

```auto
ds2 = ds.drop_vars('rotated_latitude_longitude')
ds2['dummy'].attrs = {}
ds2["dummy"].hvplot.quadmesh(projection=crs_cartopy)

```

It works, however the plot is not right. As expected, since source projection is not defined.

The final attempt produces correct output:

```auto
ds2["dummy"].hvplot.quadmesh(crs=crs_cartopy, projection=crs_cartopy)

```

That is, I have to specify both `crs` and `projection` arguments as cartopy projections, but also have to drop reference to the projection definition within the dataset. This looks to me as a bug, specifying `crs` as an argument should override search for projection in the dataset. Should I open an issue on GitHub?