handle monotonicity of cftime objects (xarray-contrib#582)

Co-authored-by: Deepak Cherian <[email protected]>

Author: larsbuntemeyer

cf_xarray/helpers.py

@@ -1,5 +1,6 @@
 from __future__ import annotations
 
+import datetime
 from collections.abc import Hashable, Sequence
 
 import numpy as np
@@ -230,6 +231,10 @@ def _get_core_dim_orders(core_dim_coords: dict[str, np.ndarray]) -> dict[str, st
             # Cast to float64 for safe comparison
             diffs_float = diffs.astype("float64")
             nonzero_diffs = diffs_float[diffs_float != 0]
+        elif isinstance(diffs[0], datetime.timedelta):
+            # For datetime timedelta, we use the total_seconds method
+            diffs_float = np.array([x.total_seconds() for x in diffs])
+            nonzero_diffs = diffs_float[diffs_float != 0]
         else:
             zero = 0
             nonzero_diffs = diffs[diffs != zero]
@@ -360,9 +365,17 @@ def _is_bounds_monotonic(bounds: np.ndarray) -> bool:
     # Cannot cast ufunc 'greater' input 0 from dtype('<m8[ns]') to dtype('<m8')
     # with casting rule 'same_kind' To avoid this, always cast to float64 before
     # np.diff.
-    arr_numeric = bounds.astype("float64").flatten()
-    diffs = np.diff(arr_numeric)
-    nonzero_diffs = diffs[diffs != 0]
+
+    diffs = np.diff(bounds.flatten())
+
+    if isinstance(diffs[0], datetime.timedelta):
+        # For datetime timedelta, we use the total_seconds method
+        diffs_float = np.array([x.total_seconds() for x in diffs])
+
+    else:
+        diffs_float = diffs.astype("float64")
+
+    nonzero_diffs = diffs_float[diffs_float != 0]
 
     # All values are equal, treat as monotonic
     if nonzero_diffs.size == 0:

cf_xarray/tests/test_helpers.py

@@ -1,10 +1,12 @@
+import numpy as np
 import xarray as xr
 from numpy.testing import assert_array_equal
 from xarray.testing import assert_equal
 
 import cf_xarray as cfxr  # noqa
 
 from ..datasets import airds, mollwds, rotds
+from . import requires_cftime
 
 try:
     from dask.array import Array as DaskArray
@@ -121,3 +123,36 @@ def test_vertices_to_bounds() -> None:
     # 2D case
     lon_b = cfxr.vertices_to_bounds(mollwds.lon_vertices, out_dims=("bounds", "x", "y"))
     assert_array_equal(mollwds.lon_bounds, lon_b)
+
+
+@requires_cftime
+def test_bounds_to_vertices_cftime() -> None:
+    import cftime
+
+    # Create cftime objects for monthly bounds
+    periods = 3
+    # start = cftime.DatetimeGregorian(2000, 1, 1)
+    edges = [cftime.DatetimeGregorian(2000, m, 1) for m in range(1, periods + 2)]
+
+    # Bounds as [start, end) for each month
+    bnds = np.array([[edges[i], edges[i + 1]] for i in range(periods)])
+    mid = np.array([edges[i] + (edges[i + 1] - edges[i]) / 2 for i in range(periods)])
+
+    # Sample data
+    values = xr.DataArray(
+        np.arange(periods, dtype=float), dims=("time",), coords={"time": mid}
+    )
+
+    # Build dataset with CF-style bounds
+    ds = xr.Dataset(
+        {"foo": values},
+        coords={
+            "time": ("time", mid, {"bounds": "time_bounds"}),
+            "time_bounds": (("time", "bounds"), bnds),
+            "bounds": ("bounds", [0, 1]),
+        },
+    )
+
+    time_c = cfxr.bounds_to_vertices(ds.time_bounds, "bounds")
+    time_b = cfxr.vertices_to_bounds(time_c, out_dims=("bounds", "time"))
+    assert_array_equal(ds.time_bounds, time_b)