Grid Mappings#
See
cf_xarray
understands the concept of coordinate projections using the grid_mapping attribute convention. For example, the dataset might contain two sets of coordinates:
native coordinates in which the data is defined, e.g., regular 1D coordinates
projected coordinates which probably denote some “real” coordinates in latitude and longitude
Due to the projection, those real coordinates are probably 2D data variables. The grid_mapping
attribute of a data variable makes a connection to another data variable defining the coordinate reference system (CRS) of those native coordinates. It should enable you to project the native coordinates into any other CRS, including the real 2D latitude and longitude coordinates. This is often useful for plotting, e.g., you can tell cartopy how to correctly plot coastlines for the CRS your data is defined in.
Extracting grid mapping info#
Dataset#
To access grid_mapping
attributes, consider this example:
from cf_xarray.datasets import rotds
rotds
<xarray.Dataset> Dimensions: (rlon: 3, rlat: 3, bounds: 4) Coordinates: * rlon (rlon) float64 17.93 18.05 18.16 * rlat (rlat) float64 21.61 21.73 21.84 lon (rlon, rlat) float64 64.22 64.42 64.63 ... 64.55 64.76 64.96 lat (rlon, rlat) float64 66.64 66.58 66.52 ... 66.81 66.75 66.69 Dimensions without coordinates: bounds Data variables: lon_bounds (bounds, rlon, rlat) float64 64.03 64.24 64.44 ... 64.74 64.95 lat_bounds (bounds, rlon, rlat) float64 66.63 66.56 66.5 ... 66.83 66.76 rotated_pole int32 0 temp (rlat, rlon) float64 0.6974 0.2508 0.7846 ... 0.9248 0.9372
The related grid mappings can be discovered using Dataset.cf.grid_mapping_names
which maps a
“grid mapping name” to the
appropriate variable name:
rotds.cf.grid_mapping_names
{'rotated_latitude_longitude': ['rotated_pole']}
Access the grid_mapping
variable as
rotds.cf["grid_mapping"]
<xarray.DataArray 'rotated_pole' ()> array(0, dtype=int32) Attributes: grid_mapping_name: rotated_latitude_longitude grid_north_pole_latitude: 39.25 grid_north_pole_longitude: -162.0
DataArrays#
Grid mapping variables are propagated when extracting DataArrays:
da = rotds.cf["temp"]
da
<xarray.DataArray 'temp' (rlat: 3, rlon: 3)> array([[0.69744995, 0.25078857, 0.78463255], [0.16888358, 0.87045378, 0.51160585], [0.9344603 , 0.92478874, 0.93716378]]) Coordinates: * rlon (rlon) float64 17.93 18.05 18.16 * rlat (rlat) float64 21.61 21.73 21.84 rotated_pole int32 0 Attributes: standard_name: air_temperature grid_mapping: rotated_pole
To find the grid mapping name use the singular DataArray.cf.grid_mapping_name
da.cf.grid_mapping_name
'rotated_latitude_longitude'
And to get the grid mapping variable
da.cf["grid_mapping"]
<xarray.DataArray 'rotated_pole' ()> array(0, dtype=int32) Attributes: grid_mapping_name: rotated_latitude_longitude grid_north_pole_latitude: 39.25 grid_north_pole_longitude: -162.0
Use grid_mapping
in projections#
The grid mapping information use very useful in projections, e.g., for plotting. pyproj understands CF conventions right away, e.g.
from pyproj import CRS
CRS.from_cf(rotds.cf["grid_mapping"].attrs)
gives you more details on the projection:
<Derived Geographic 2D CRS: {"$schema": "https://proj.org/schemas/v0.2/projjso ...>
Name: undefined
Axis Info [ellipsoidal]:
- lon[east]: Longitude (degree)
- lat[north]: Latitude (degree)
Area of Use:
- undefined
Coordinate Operation:
- name: Pole rotation (netCDF CF convention)
- method: Pole rotation (netCDF CF convention)
Datum: World Geodetic System 1984
- Ellipsoid: WGS 84
- Prime Meridian: Greenwich
For use in cartopy, there is some more overhead due to this issue. So you should select the right cartopy CRS and just feed in the grid mapping info:
from cartopy import crs as ccrs
grid_mapping = rotds.cf["grid_mapping"]
pole_latitude = grid_mapping.grid_north_pole_latitude
pole_longitude = grid_mapping.grid_north_pole_longitude
ccrs.RotatedPole(pole_longitude, pole_latitude)