check precision of vertex order (#361)

* changed precision of vertex order check comparison

* [pre-commit.ci] auto fixes from pre-commit.com hooks

for more information, see https://pre-commit.ci

* added rotated pole datasets with inexact bounds

* renamed test datasets

* added test on bounds with precision issues

* clean up

* change check coordinate

Co-authored-by: github action <[email protected]>
Co-authored-by: pre-commit-ci[bot] <66853113+pre-commit-ci[bot]@users.noreply.github.com>
Co-authored-by: Deepak Cherian <[email protected]>

Author: 4 people

cf_xarray/datasets.py

@@ -221,6 +221,105 @@ def _create_mollw_ds():
 
 mollwds = _create_mollw_ds()
 
+
+def _create_inexact_bounds():
+    # Dataset that creates rotated pole curvilinear coordinates with CF bounds in
+    # counterclockwise order that have precision issues.
+    lon = np.array(
+        [
+            [64.21746363939087, 64.42305921561967, 64.62774455060337],
+            [64.38488380781622, 64.5906340869802, 64.79546745916944],
+            [64.55349991538067, 64.75940009330206, 64.96437666717893],
+        ]
+    )
+
+    lat = np.array(
+        [
+            [66.63852914622773, 66.57692364510149, 66.51510790697783],
+            [66.72632651121215, 66.66454929137457, 66.60256232699675],
+            [66.81394503895442, 66.75199541041104, 66.68983654206977],
+        ]
+    )
+
+    lon_bounds = np.array(
+        [
+            [
+                [64.03109895962405, 64.23707042152387, 64.44213290724275],
+                [64.19784426308645, 64.40397614495, 64.60919235424163],
+                [64.36578231818473, 64.57206989228206, 64.77743506492257],
+            ],
+            [
+                [64.23707042152385, 64.44213290724275, 64.64629053434521],
+                [64.40397614494998, 64.60919235424163, 64.81349710155956],
+                [64.57206989228204, 64.77743506492257, 64.98188214131258],
+            ],
+            [
+                [64.40397614494998, 64.60919235424161, 64.81349710155956],
+                [64.57206989228204, 64.77743506492257, 64.98188214131258],
+                [64.74136272095073, 64.94687197648351, 65.15145647131322],
+            ],
+            [
+                [64.19784426308645, 64.40397614495, 64.60919235424161],
+                [64.36578231818473, 64.57206989228206, 64.77743506492257],
+                [64.53492430330854, 64.74136272095075, 64.94687197648351],
+            ],
+        ]
+    )
+
+    lat_bounds = np.array(
+        [
+            [
+                [66.62524451217388, 66.56383049714253, 66.50220514181308],
+                [66.71321640554079, 66.65163077691872, 66.58983429068188],
+                [66.801010821918, 66.73925288286632, 66.67728458111449],
+            ],
+            [
+                [66.56383049714253, 66.50220514181308, 66.44037016643905],
+                [66.65163077691872, 66.58983429068188, 66.52782867442114],
+                [66.73925288286632, 66.67728458111449, 66.61510765153199],
+            ],
+            [
+                [66.6516307769187, 66.58983429068185, 66.52782867442114],
+                [66.73925288286632, 66.67728458111449, 66.61510765153199],
+                [66.82669478850752, 66.76455398820585, 66.702205074604],
+            ],
+            [
+                [66.71321640554079, 66.6516307769187, 66.58983429068185],
+                [66.801010821918, 66.73925288286632, 66.67728458111449],
+                [66.88862573342278, 66.82669478850752, 66.76455398820585],
+            ],
+        ]
+    )
+
+    rotated = xr.Dataset(
+        coords=dict(
+            lon=xr.DataArray(
+                lon,
+                dims=("x", "y"),
+                attrs={"units": "degrees_east", "bounds": "lon_bounds"},
+            ),
+            lat=xr.DataArray(
+                lat,
+                dims=("x", "y"),
+                attrs={"units": "degrees_north", "bounds": "lat_bounds"},
+            ),
+        ),
+        data_vars=dict(
+            lon_bounds=xr.DataArray(
+                lon_bounds, dims=("bounds", "x", "y"), attrs={"units": "degrees_east"}
+            ),
+            lat_bounds=xr.DataArray(
+                lat_bounds, dims=("bounds", "x", "y"), attrs={"units": "degrees_north"}
+            ),
+        ),
+    )
+
+    return rotated
+
+
+rotds = _create_inexact_bounds()
+
+
 forecast = xr.decode_cf(
     xr.Dataset.from_dict(
         {

cf_xarray/helpers.py

@@ -95,7 +95,9 @@ def _bounds_helper(values, n_core_dims, nbounds, order):
             vertex_vals = np.block([[bot_left, bot_right], [top_left, top_right]])
         if order is None:  # We verify if the ccw version works.
             calc_bnds = np.moveaxis(vertices_to_bounds(vertex_vals).values, 0, -1)
-            order = "counterclockwise" if np.all(calc_bnds == values) else "clockwise"
+            order = (
+                "counterclockwise" if np.allclose(calc_bnds, values) else "clockwise"
+            )
         if order == "clockwise":
             bot_left = values[..., :, :, 0]
             top_left = values[..., -1:, :, 1]

cf_xarray/tests/test_helpers.py

@@ -3,7 +3,7 @@
 
 import cf_xarray as cfxr  # noqa
 
-from ..datasets import airds, mollwds
+from ..datasets import airds, mollwds, rotds
 
 try:
     from dask.array import Array as DaskArray
@@ -27,6 +27,13 @@ def test_bounds_to_vertices():
     assert_equal(mollwds.lat_vertices, lat_ccw)
     assert_equal(lat_no, lat_ccw)
 
+    # 2D case with precision issues, check if CF- order is "detected" correctly
+    lon_ccw = cfxr.bounds_to_vertices(
+        rotds.lon_bounds, bounds_dim="bounds", order="counterclockwise"
+    )
+    lon_no = cfxr.bounds_to_vertices(rotds.lon_bounds, bounds_dim="bounds", order=None)
+    assert_equal(lon_no, lon_ccw)
+
     # Transposing the array changes the bounds direction
     ds = mollwds.transpose("bounds", "y", "x", "y_vertices", "x_vertices")
     lon_c = cfxr.bounds_to_vertices(