adcroft
commented
5 years ago @rabernat makes very good points and I agree with everything said.
The concept of a space-filling and contiguous mesh is fundamental to many finite-volume codes and I find it baffling that a more efficient description has not already replaced the current description.
The names of the mesh coordinates do not need to be fixed. They can be named in an attribute, e.g. :mesh, attached to the variable. Stylizing on the :coordinates attribute from 5.2. Two-Dimensional Latitude, Longitude, Coordinate Variables an example could look something like this:
dimensions:
i = 128 ;
j = 64 ;
k = 18 ;
im = 129 ;
jm = 65 ;
km = 19 ;
variables:
float thetao(k,j,i) ;
thetao:long_name = "seawater_potential_temperature" ;
thetao:units = "C" ;
thetao:coordinates = "lon lat z" ;
thetao:mesh = "meshlon meshlat zi" ;
float lon(j,i) ;
lon:long_name = "longitude" ;
lon:units = "degrees_east" ;
float lat(j,i) ;
lat:long_name = "latitude" ;
lat:units = "degrees_north" ;
float z(k) ;
z:long_name = "elevation" ;
z:units = "m" ;
z:positive = "up" ;
float meshlon(jm,im) ;
meshlon:long_name = "longitude" ;
meshlon:units = "degrees_east" ;
float meshlat(jm,im) ;
meshlat:long_name = "latitude" ;
meshlat:units = "degrees_north" ;
float zi(km) ;
zi:long_name = "elevation" ;
zi:units = "m" ;
zi:positive = "up" ;Several notes:
- I threw away the unused 1D coordinate variables for brevity.
- I'm assuming the
coordinatesattribute can be used for more than just the horizontal coordinates longitude-latitude. I think that's legitimate but didn't find an example in a quick scan. :mesh_coordinatesmight be a better name for the attribute but I like shorter names myself.- There should probably be a global attribute specifying the version of the standard being used in the file.
At the risk of derailing the thread I'll mention that we can arrive at an even more succinct description if we discard the idea of data coordinates which I think is unnecessary when providing the finite volume mesh information:
dimensions:
i = 128 ;
j = 64 ;
k = 18 ;
im = 129 ;
jm = 65 ;
km = 19 ;
variables:
float thetao(k,j,i) ;
thetao:long_name = "seawater_potential_temperature" ;
thetao:units = "C" ;
thetao:mesh = "meshlon meshlat zi" ;
float meshlon(jm,im) ;
meshlon:long_name = "longitude" ;
meshlon:units = "degrees_east" ;
float meshlat(jm,im) ;
meshlat:long_name = "latitude" ;
meshlat:units = "degrees_north" ;
float zi(km) ;
zi:long_name = "elevation" ;
zi:units = "m" ;
zi:positive = "up" ;I believe the latter example to be a more appropriate description of the data. However, I realize that the latter is not backward compatible since there is an expectation that there are coordinates and will cause trouble if they are missing. So I will be pleased enough if we can at least adopt an attribute that identifies mesh coordinates as in the earlier example.
rabernat
feggleton
JonathanGregory
JimBiardCics
davidhassell
neumannd
taylor13
Chilipp
ChrisBarker-NOAA
zklaus
I work every day with ocean models that use orthogonal curvilinear coordinates (MITgcm, MOM, POP, ROMS, NEMO, etc. etc.). This is an example tripolar grid from CESM:
The grid cells in such models are contiguous quads, with four points specifying the lat / lon vertex locations of each cell. CF conventions tell me (Section 7.1: Cell Boundaries) that I should use a boundary_variable.
This convention is general enough to accommodate potentially overlapping or non-contiguous quads, essentially
n x mtotally unrelated four-sided shapes.My main point: It's inefficient to store structured grid geometry this way.
I don't want to have to check this, I want the conventions to tell me. In our latest global high-resolution ocean models, I have a mesh that is of size
n=12960, m=17280, 223 million cells. I am interested in streamlining my analysis and visualization workflow as much as possible, which means minimizing the required memory and computational steps.Instead of specifying a boundary variable, I propose to introduce the concept of a mesh variable, with the following conventions:
A mesh variable will have the same number of dimensions as its associated coordinate or auxiliary coordinate variable, but with one extra element in each dimension.
In the case where the horizontal grid is described by two-dimensional auxiliary coordinate variables in latitude
lat(n,m)and longitudelon(n,m), and the associated cells are four-sided and contiguous, then the mesh variables are given in the formlatmesh(n+1, m+1)andlonmesh(n+1, m+1).It would not be hard to generate such data, since this is how most GCMs keep track of their own coordinate grids internally (e.g. MITgcm). This convention also aligns well with how most visualization software plots such data, e.g. matplotlib's pcolormesh function. So adding something like this to the CF conventions would streamline the path from model output to plotting, eliminating the potentially error-prone step of encoding, and then decoding, the "boundary variable" type coordinates. For the dataset I described above, the difference is about 3 GB of memory.
I don't feel strongly about what it's called. Maybe "mesh variable" is not the right choice. But I feel something like this is sorely needed.
cc @adcroft & @StephenGriffies, with whom this topic has come up repeatedly.