NOAA-ORR-ERD / gridded

A single API for accessing / working with gridded model results on multiple grid types
https://noaa-orr-erd.github.io/gridded/index.html
The Unlicense
64 stars 14 forks source link

Add method to build the face edge connectivity #60

Closed Chilipp closed 3 years ago

Chilipp commented 3 years ago

as I'll need this for the dual mesh, I thought I start with implementing the face_edge_connectivity. I use the pandas MultiIndex here as it is very efficient (4 seconds for a highres mesh with 1.2 million faces). It add's another dependency, though, but I think pandas is that common nowadays (even within the netCDF-community due to xarray) that this is not a problem. I made pandas an optional dependency (i.e. I only added it to the test requirements file).

What do you think about this approach @ChrisBarker-NOAA and @jay-hennen?

Chilipp commented 3 years ago

One thought: I can't remember the exact spec of the UGRID standard, but is there anything in there about clockwise vs anti-clockwise? or the relationship between the node number and the edge number?

Yep, you're right. The UGRID conventions on 2D triangular mesh topology say:

The attribute face_node_connectivity points to an index variable identifying for every face (here consistently triangle) the indices of its three corner nodes. The corner nodes should be specified in anticlockwise (also referred to as counterclockwise) direction as viewed from above (consistent with the CF-convention for bounds of p-sided cells).

[...]

Optionally the topology may have the following attributes: face_edge_connectivity pointing to an index variable identifying for every face (here consistently triangle) the indices of its three edges. The edges should be specified in anticlockwise direction as viewed from above.

I am calculating the face_edge_connectivity from the face_node_connectivity and maintain the order. Hence, everything should be sorted correctly (if the face_node_connectivity is correct)

Chilipp commented 3 years ago

I'd rather have a pure-numpy solution, but I don't have time to figure that out :-)

sure :sweat_smile: But pure numpy is a bit tricky if you want to avoid looping over the full array in python. one could use the scipys cKDTree if you prefer scipy over pandas. This would also make it twice as fast (i.e. 2 seconds for the 1.2 million triangles instead of 4). I think scipy is anyway already a requirement, is it?

Chilipp commented 3 years ago

alright, switched to scipy and remove all pandas stuff with https://github.com/NOAA-ORR-ERD/gridded/pull/60/commits/30cd5f224584d52bffac5802b54d8643753f636c

Chilipp commented 3 years ago

the CI is failing due to some problems with netCDF4 (i.e. unrelated to this PR). Feel free to merge it if you are happy with it or please let me know if there is more to it :)

Chilipp commented 3 years ago

Hey @ChrisBarker-NOAA! I did not get your comments, sorry. Don't know why they have not been submitted.

Anyway: concerning the -999. This is only in the test file. I am using a masked array there as this is how netcdf4 does handle it. As you said, for a mixed topology of quads and triangles, for instance, the UGRID conventions say to use a _Fillvalue. Hence, netcdf4 turns the faces variable in the read_netcdf.py into a masked array. But I don't know if the other constructor methods act the same way.

ChrisBarker-NOAA commented 3 years ago

yeah, I think we need to put this in an issue to make sure it's consistent.