Handling Duplicate Node Indices

[philipc2]

Some of our testing grids include duplicates of certain nodes, meaning that there existing multiple nodes that share the same latitude and longitude location, but each have their own index within our connectivity and data variables.

Example

Consider the following simple grid (

node_lon = (-180, 0, -180, 180)
node_lat = (-90, 0, -90, 90)
data_on_nodes = [1, 2, 1, 3]

Indices 0 and 1 contain the same latitude, longitude, and data value. This causes issues when constructing certain connectivity's (such as node_face, since some nodes will only belong to one face.

Affected Grids

The following grid files in our testing suite each have duplicate nodes.
geos-cs/c1/test-c12.native.nc4
ugrid/geoflow-small/grid.nc
esmf/ne30/ne30pg3.grid.nc

All of these are some form of a cube-sphere (or variation)

Proposed Solution

Instead of deleting the duplicates (which would require a complete re-ordering of connectivity and to modify the size of all node-centered data variables), a simple masking and merging of node indices at the connectivity level should do the trick.

Example

Let's say node indicies [0, 5, 10, 23] all share the same data value and lat/lon coordinates. We can itterate over the face_node, edge_node and node_node connectivity arrays and set the values of [5, 10, 23] to point to 0. This would allow us to leave the data arrays untouched and do a single pass over our connectivity's arrays with a direct replacement, without needing to reshape them in any way.

[philipc2]

Need to consider cases when points may lie on the pole or antemeridian (i.e. (-180, 90) and (180, 90))

[philipc2]

Consider raising an error when trying run certain functionality that would fail.

[rajeeja]

Adding to the list:

Add a tolerance with default, but user can specify if needed.

[hongyuchen1030]

Hi Team,

I wanted to add a couple of notes to help us move forward:

1. The error tolerance is already set here with a detailed description.
2. The handling of the pole and antimeridian cases is already incorporated in my coordinate transformation function, along with thorough documentation. For the pole point, it is uniformly set to [0, 0.5 * pi] or [0, 0.5 * -pi]. You can find this here.

For the antimeridian case, I initially set the range between [0, 2*pi) and clearly documented that all our algorithms are based on this range. It seems that the current range has been changed to (-pi, pi). Updating the API documentation for _xyz_to_lonlat_rad and _xyz_to_lonlat_deg to reflect this change would be beneficial.

Unifying the coordinate representation and providing detailed documentation will ensure the robustness of our project. I believe these enhancements will greatly benefit our project.

Thank you for all your efforts and contributions!