Center vs Edge Vertical Coordinates

[bena-nasa]

Question

I've been trying to parse the CF documentation to find an answer to this unsuccessfully. I'm developer of NASA GEOS model which runs on hybrid sigma pressure levels in the atmosphere. I'm writing some new code and I'm trying to do this "right", i.e. follow the CF conventions. Fields in the model can be defined on the full vertical level or on the edges (Edge pressure being the canonical example). So akin staggering in say the horizontal grid. Practically this means that say we have 72 levels, the NetCDF file would have a dimension that must be 72 and another that is 73. Now from my reading of the conventions, 7.1.4 is the example to follow, it says for the centers I should define a coordinate variable that is of size 72 (they use eta in the example), with a formula terms and bounds (eta_bnds) in the example. The file then has ak and bk of size 72 and ak_bnds and bk_bnds of size (72,2)

So any data variable on the centers will depend on the eta dimension that is size 72. Great, but now I also have a data variable that is defined on the edges of the levels in the file as well so it needs other to depend on a dimension of size 73. It still represents the same vertical levels, but the edges. I can't have the data variable depend on the say the eta_bnds this results in duplication of data that is unacceptable from a storage perspective.

So I can't for the life of me figure then what his file should look like to conform to what I'm reading in the CF documentation. Should the coordinate variable for the edges, still be defined with formula terms and just point to different variables in the file for ak and bk that are of size 73 but still point to the same surface pressure?

I just don't see anything about this center vs edge problem other than inferring if the coordinate variable has bounds, well it must be the "center" of the cell. Surely other groups must deal with this? And while having these extra bounds coordinates are fine for metadata describing the grid, the storage penalty the actual data is too big.

[taylor13]

If a variable is carried on 72 levels, you define a coordinate variable, say "levcenters(72)". If another variable is carried on 73 levels, then you define another coordinate variable, say "levedges(72)".

You may also indicate that the coordinates for centers represent grid cells of certain depths. You do that by defining the "bounds" associated with levcenters.

CF does not indicate the relationship between levcenters and levedges.

[bena-nasa]

I realize that I can define the bounds associated with the full levels/centers as in section 7.1.4 of the conventions where there is an example for hybrid sigma pressure coordinates. I need another coordinate variable for the edges that is one bigger than the centers, so that data variables that are defined on the edges will have a dimension they can use, and I can identify that they have a vertical coordinate. But, once again, they are both describing the same hybrid sigma levels at the end of the day. Just different "staggerings".

This doesn't any different than if say you had many a horizontal grid that had variables defined at different Arakawa staggers (A/C/D), it's all the same grid. I was a bit surprised to see nothing about his in the standard, unless I missed it. I realize you define cell, bounds but you can't use those for the data variables, too much storage duplication...

Something seems missing then from the standard.

But I guess the answer for now, is just define a coordinate variable for the edges, it can have it's own formula terms than that point to different ak/bk variables that just duplicate what's in the ak_bnds,bk_bnds for the centers at the end of the day

[ChrisBarker-NOAA]

This doesn't any different than if say you had many a horizontal grid that had variables defined at different Arakawa staggers (A/C/D), it's all the same grid. I was a bit surprised to see nothing about his in the standard, unless I missed it.

No, you didn't -- it is missing from the standard.

I realize you define cell, bounds but you can't use those for the data variables, too much storage duplication...

And even then, I'm not sure you can define what you want anyway -- defining the cells is overly redundant (I'd like to see a change there), but do-able, but then you need cell_methods to define what teh realtionship is between teh variable and the cell, and I not sure there's a way to express, e.g. flux through a cell side.

For the Arakawa staggers, there is the SGRID proposal:

https://sgrid.github.io/sgrid/

https://github.com/sgrid/sgrid

It was developed some years ago, and is in use at least a little bit (Deltares, I think uses it, and at least some ROMS implementations use it.) I think it should be considered for inclusion in CF, but haven't found the time or energy to propose that formally (nor has, apparently, anyone else involved with the project).

I'm not sure if it addresses your needs for the vertical, but I"d be interested to see what you think about it.

At the center of this is that CF is trying to be as general as possible, and not have more than one way to express the same thing. Which is a laudable goal. However, a model grid (or mesh) is not just a collection of individual cells, and I think a lot is lost if you lose the information contained in the grid itself, as a single "thing".

[JonathanGregory]

Dear @bena-nasa

As Karl and Chris have said, the staggered vertical grids are described independently in CF. The vertical grid of dimension 72 has a vertical coordinate (72) and layer bounds (72,2), in which all values apart from the top and bottom appear twice, if the layers are contiguous. If these 73 distinct values are the levels in another vertical grid, it has a dimension of size 73, a vertical coordinate variable (73) and layer bounds (73,2) if needed, containing 74 distinct values, of which 72 are the levels of the first grid. This certainly involves duplication of information, but it is more general and easy to use. Since the coordinate variables are only 1D, they are usually much smaller than the data variables, which are 2D, 3D or 4D, so the duplication isn't practically a problem in most cases. Is it a problem in your application?

You won't be surprised to know that this sort of question has been discussed many times in CF! Quite recently we talked about the corresponding question for horizontal coordinate variables. Below is an extract from one of my comments in that discussion.

Best wishes

Jonathan

---

The CF convention that bounds for 1D coordinate variables are dimensioned (n,2) means that in the common case of contiguous cells every bound is duplicated. Each of the (n,m,4) bounds of auxiliary lat-lon coordinate variables on a rectangular x-y grid has four copies. This wasn't an accident, but a deliberate decision at the start of CF. It was debated several times in early years. It is explicitly explained and acknowledged in Sect 7.1:

Note that the boundary variable for a set of N contiguous intervals is an array of shape (N,2). Although in this case there will be a duplication of the boundary coordinates between adjacent intervals, this representation has the advantage that it is general enough to handle, without modification, non-contiguous intervals, as well as intervals on an axis using the unlimited dimension.

Also it allows overlapping cells, for which we have use-cases. There is no CF way to indicate that the cells are contiguous. To find out if they are, you have to examine the bounds. We could provide a convention to indicate it, but that would redundant, and therefore liable to become inconsistent. We have a principle of not doing that in Sect 1.2.

To avoid potential inconsistency within the metadata, the conventions should minimise redundancy.

Although the (N,2) convention is wasteful of space in the common cases, it makes the document and software simpler if we have only one convention. I think there would need to be a very strong reason to add any alternative convention, for the reasons summarised Sect 1.2:

Because all previous versions must generally continue to be supported in software for the sake of archived datasets, and in order to limit the complexity of the conventions, there is a strong preference against introducing any new capability to the conventions when there is already some method that can adequately serve the same purpose (even if a different method would arguably be better than the existing one).

The wastefulness is not great for 2D grids made of 1D coordinates, because the bounds take size of order N, while the data is N^2.