NCAS-CMS / cf-python

A CF-compliant Earth Science data analysis library
http://ncas-cms.github.io/cf-python
MIT License
126 stars 19 forks source link

Return (identifier for) construct added by `compute_vertical_coordinates` #802

Open sadielbartholomew opened 3 months ago

sadielbartholomew commented 3 months ago

When using compute_vertical_coordinates to compute the non-parametric vertical coordinates from the coordinate reference and parameters on the field, the method returns (when not in-place) the field that was operated on, but this does not include information about the specific coordinate created, or indeed whether a coordinate construct was successfully created, since running the method on fields without the necessary information returns the same field (see snippet 2 below).

I'm requesting a means to access either the coordinate itself that gets created, via the object itself or an identifier e.g. the unique key for it, either as:

  1. an alternative return output, e.g. configured via a flag keyword to change from the field object return, or;
  2. an additional output e.g. as a 2-tuple with the resultant field.

By returning None (or False, but None seems more appropriate) for cases where no such coordinate could be created, this would also enable a null operation to be identified directly.

Examples

Snippet 1: successful operation

>>> f  # a field with a coordinate reference with all necessary info, so the method should and does work
<CF Field: ncvar%T(time(61), atmosphere_hybrid_sigma_pressure_coordinate(50), ncdim%south_north(179), ncdim%west_east(139))>
>>> f.coordinate_reference()
<CF CoordinateReference: standard_name:atmosphere_hybrid_sigma_pressure_coordinate>
>>> f.coordinate_reference().dump()
Coordinate Reference: standard_name:atmosphere_hybrid_sigma_pressure_coordinate
    Coordinate conversion:computed_standard_name = air_pressure
    Coordinate conversion:standard_name = atmosphere_hybrid_sigma_pressure_coordinate
    Coordinate conversion:a = domainancillary0
    Coordinate conversion:b = domainancillary1
    Coordinate conversion:p0 = domainancillary2
    Coordinate conversion:ps = domainancillary3
    Coordinate: dimensioncoordinate1
>>> # Now compute from the coordinate reference
>>> f.compute_vertical_coordinates(inplace=True)
>>> g = f.compute_vertical_coordinates()
>>> # It is not clear, without manual inspection, what coordinates were added, if any,
>>> # without comparing g and f or noticing the 'air_pressure' aux. coord.
>>> g.equals(f)
False
>>> print(f)
Field: ncvar%T (ncvar%T)
------------------------
Data            : ncvar%T(time(61), atmosphere_hybrid_sigma_pressure_coordinate(50), ncdim%south_north(179), ncdim%west_east(139))
Dimension coords: time(61) = [2023-07-10 12:00:00, ..., 2023-07-13 00:00:00] proleptic_gregorian
                : atmosphere_hybrid_sigma_pressure_coordinate(50) = [0.9969073534011841, ..., 0.0019467439269647002]
                : forecast_reference_time(1) = [2023-07-10 12:00:00] proleptic_gregorian
Auxiliary coords: ncvar%XLAT(time(61), ncdim%south_north(179), ncdim%west_east(139)) = [[[36.98377990722656, ..., 69.75613403320312]]] degrees_north
                : ncvar%XLONG(time(61), ncdim%south_north(179), ncdim%west_east(139)) = [[[-22.838226318359375, ..., 20.62158203125]]] degrees_east
                : ncvar%XTIME(time(61)) = [0.0, ..., 3600.0]
                : air_pressure(atmosphere_hybrid_sigma_pressure_coordinate(50), time(61), ncdim%south_north(179), ncdim%west_east(139)) = [[[[1021.5154515355825, ..., 51.87520217895508]]]] 100 Pa
Coord references: standard_name:atmosphere_hybrid_sigma_pressure_coordinate
Domain ancils   : ncvar%A(atmosphere_hybrid_sigma_pressure_coordinate(50)) = [0.00019427388906478882, ..., 0.05119684338569641]
                : ncvar%B(atmosphere_hybrid_sigma_pressure_coordinate(50)) = [0.9968656897544861, ..., 0.0]
                : ncvar%p0() = 1013.25 hPa
                : ncvar%ps(time(61), ncdim%south_north(179), ncdim%west_east(139)) = [[[1024.52978515625, ..., 972.4573364257812]]] hPa

In this case, I would like to have the air_pressure auxiliary coordinate, or its key, returned by the method, since this is the one that is computed and created on the field.

Snippet 2: null operation

No coordinates can be created here, so the field is unchanged, but we can't tell that without comparing fields:

>>> # Example of a null operation
>>> c = cf.example_field(0)
>>> cf.log_level(-1)
>>> c1 = c.compute_vertical_coordinates()
>>> c1.equals(c)
True

In this case, I would like: