Axes and Coordinates

One powerful feature of cf_xarray is the ability to refer to named dimensions by standard axis or coordinate names in Dataset or DataArray methods.

For example, one can call ds.cf.mean("latitude") instead of ds.mean("lat")

from cf_xarray.datasets import airds

# identical to airds.mean("lat")
airds.cf.mean("latitude")
<xarray.Dataset> Size: 2kB
Dimensions:  (time: 4, lon: 50)
Coordinates:
  * lon      (lon) float32 200B 200.0 202.5 205.0 207.5 ... 317.5 320.0 322.5
  * time     (time) datetime64[ns] 32B 2013-01-01 ... 2013-01-01T18:00:00
Data variables:
    air      (time, lon) float64 2kB 279.4 279.7 279.7 ... 276.6 277.8 279.0

Tip

Most xarray methods are wrapped by cf-xarray. Simply access them as DataArray.cf.method(dim="latitude") for example and try it! If something does not work, please raise an issue.

Coordinate Criteria

How does this work? cf_xarray has an internal table of criteria (mostly copied from MetPy) that lets it identify specific coordinate names "latitude", "longitude", "vertical", "time".

Tip

See Custom Criteria to find out how to define your own custom criteria.

This table lists these internal criteria

latitude

longitude

time

vertical

_CoordinateAxisType

Lat

Lon

Time

axis

T

cartesian_axis

T

grads_dim

t

long_name

latitude

longitude

time

air_pressure, altitude, depth, geopotential_height, height, height_above_geopotential_datum, height_above_mean_sea_level, height_above_reference_ellipsoid

positive

down, up

standard_name

latitude

longitude

time

air_pressure, altitude, depth, geopotential_height, height, height_above_geopotential_datum, height_above_mean_sea_level, height_above_reference_ellipsoid

units

degreeN, degree_N, degree_north, degreesN, degrees_N, degrees_north

degreeE, degree_E, degree_east, degreesE, degrees_E, degrees_east

Any DataArray that has standard_name: "latitude" or _CoordinateAxisType: "Lat" or "units": "degrees_north" in its attrs will be identified as the "latitude" variable by cf-xarray. Similarly for other coordinate names.

Axis Names

Similar criteria exist for the concept of “axes”.

T

X

Y

Z

_CoordinateAxisType

Time

GeoX

GeoY

GeoZ, Height, Pressure

axis

T

X

Y

Z

cartesian_axis

T

X

Y

Z

grads_dim

t

x

y

z

long_name

time

cell index along first dimension, grid_longitude, projection_x_angular_coordinate, projection_x_coordinate

cell index along second dimension, grid_latitude, projection_y_angular_coordinate, projection_y_coordinate

atmosphere_hybrid_height_coordinate, atmosphere_hybrid_sigma_pressure_coordinate, atmosphere_ln_pressure_coordinate, atmosphere_sigma_coordinate, atmosphere_sleve_coordinate, model_level_number, ocean_double_sigma_coordinate, ocean_s_coordinate, ocean_s_coordinate_g1, ocean_s_coordinate_g2, ocean_sigma_coordinate, ocean_sigma_z_coordinate

positive

standard_name

time

grid_longitude, projection_x_angular_coordinate, projection_x_coordinate

grid_latitude, projection_y_angular_coordinate, projection_y_coordinate

atmosphere_hybrid_height_coordinate, atmosphere_hybrid_sigma_pressure_coordinate, atmosphere_ln_pressure_coordinate, atmosphere_sigma_coordinate, atmosphere_sleve_coordinate, model_level_number, ocean_double_sigma_coordinate, ocean_s_coordinate, ocean_s_coordinate_g1, ocean_s_coordinate_g2, ocean_sigma_coordinate, ocean_sigma_z_coordinate

units

.axes and .coordinates properties

Alternatively use the special properties DataArray.cf.axes or DataArray.cf.coordinates to access the variable names. These properties return dictionaries that map “CF names” to a list of variable names. Note that a list is always returned even if only one variable name matches the name "latitude" (for example).

airds.cf.axes
{'X': ['lon'], 'Y': ['lat'], 'T': ['time']}
airds.cf.coordinates
{'longitude': ['lon'], 'latitude': ['lat'], 'time': ['time']}

Axes or Coordinate?

TODO describe latitude vs Y; longitude vs X; vertical vs Z

Checking presence of axis or coordinate

Note that a given “CF name” is only present if there is at least one variable that can be identified with that name. The airds dataset has no "vertical" coordinate or "Z" axis, so those keys are not present. So to check whether a "vertical" coordinate or "Z" axis is present, one can

"Z" in airds.cf.axes
False
"vertical" in airds.cf.coordinates
False

Or one can check the dataset as a whole:

"Z" in airds.cf
False

Using the repr

It is always useful to check the variables identified by cf-xarray using the repr

airds.cf
Coordinates:
             CF Axes: * X: ['lon']
                      * Y: ['lat']
                      * T: ['time']
                        Z: n/a

      CF Coordinates: * longitude: ['lon']
                      * latitude: ['lat']
                      * time: ['time']
                        vertical: n/a

       Cell Measures:   area: ['cell_area']
                        volume: n/a

      Standard Names: * latitude: ['lat']
                      * longitude: ['lon']
                      * time: ['time']

              Bounds:   n/a

       Grid Mappings:   n/a

Data Variables:
       Cell Measures:   area, volume: n/a

      Standard Names:   air_temperature: ['air']

              Bounds:   n/a

       Grid Mappings:   n/a