Repository navigation
Return (identifier for) construct added by compute_vertical_coordinates #802
Description
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:
- an alternative return output, e.g. configured via a flag keyword to change from the field object return, or;
- 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]]] hPaIn 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)
TrueIn this case, I would like:
- to have the method return
Noneto indicate that no new coordinate construct was added; - ideally more detaails to be shared about the inability from lack of relevant metadata, etc. to calculate the vertical coordinates, since even at
-1log level I am not told anything (see above).
Activity
sadielbartholomew commentedon Jan 30, 2025
Going to work on this next week, as it is has emerged as very useful for the VISION work. We have discussed a workaround for the VISION toolkit until then, but for this new enhancement the plan is:
keythat makes the return more flexible as follows:key=Truekey=False, same behaviour as present/beforeinplace=TrueNoneplus warn via logging that these kwargs are incompatibleNoneinplace=False, computed new coord. successfully(output_field, key_of_new_coord)output_fieldinplace=False, did/could not compute new coord.(output_field, None)and document thatNonefor thekeyposition indicates no coord. createdoutput_field@davidhassell before I code it up just to check:
key. My proposal iskeyhas aFalsedefault simply because thekeyoption provides more information which can be requested at will via=True, as useful as the extra information may be. What do you think?davidhassell commentedon Jan 31, 2025
Hi Sadie,
Nearly, as I saw it:
key=Truekey=False, same behaviour as present/beforeinplace=Trueraise ValueError[1]Noneinplace=False, computed new coord. successfully(output_field, key_of_new_coord)output_fieldinplace=False, no new computed coord. [2](output_field, None)and document thatNonefor thekeyposition indicates no coord. createdoutput_field[1] I think it's better to fail here, because there is ambiguity about what the user actually expected
[2] There is no distinction between "did not" and "could not", here. If there is a coordinate to compute, but that computation fails (e.g. in some circumstances, due to a missing domain ancillary) then an exception will always be raised.
Yes,
key=Falseshould be the default for the reason you give, and to preserve backwards compatibility.sadielbartholomew commentedon Jan 31, 2025
Lovely, thanks @davidhassell - happy to upgrade to a failure for [1] and I understand RE [2] and glad we agree on the sensible default. PR to go up sometime next week.