From 58de6d202b1ca3c42d3d78a9984b8780fe84e374 Mon Sep 17 00:00:00 2001 From: Spencer Clark Date: Thu, 19 Feb 2026 11:54:10 -0500 Subject: [PATCH] Finalize change in behavior of `decode_timedelta=None` (#11173) --- doc/whats-new.rst | 11 +++++- xarray/coding/times.py | 20 +---------- xarray/conventions.py | 8 ++--- xarray/tests/test_coding_times.py | 58 +++++++++++-------------------- xarray/tests/test_conventions.py | 16 ++------- 5 files changed, 38 insertions(+), 75 deletions(-) diff --git a/doc/whats-new.rst b/doc/whats-new.rst index 1bcd2016324..f4d9819588e 100644 --- a/doc/whats-new.rst +++ b/doc/whats-new.rst @@ -17,7 +17,16 @@ New Features Breaking Changes ~~~~~~~~~~~~~~~~ - +- Xarray will no longer by default decode a variable into a + :py:class:`np.timedelta64` dtype based on the presence of a timedelta-like + ``"units"`` attribute alone. Instead it will rely on the presence of a + :py:class:`np.timedelta64` dtype attribute, which is now xarray's default way + of encoding :py:class:`np.timedelta64` values. The old decoding behavior can + be restored by specifying ``decode_timedelta=True`` or + ``decode_timedelta=CFTimedeltaCoder(decode_via_units=True)`` in + :py:meth:`open_dataset`. This finalizes the deprecation cycle initiated in + xarray version 2025.01.2 (:pull:`11173`). By `Spencer Clark + `_. Deprecations ~~~~~~~~~~~~ diff --git a/xarray/coding/times.py b/xarray/coding/times.py index d45e8e4dc9f..b60acbcf407 100644 --- a/xarray/coding/times.py +++ b/xarray/coding/times.py @@ -1492,13 +1492,12 @@ class CFTimedeltaCoder(VariableCoder): def __init__( self, time_unit: PDDatetimeUnitOptions | None = None, - decode_via_units: bool = True, + decode_via_units: bool = False, decode_via_dtype: bool = True, ) -> None: self.time_unit = time_unit self.decode_via_units = decode_via_units self.decode_via_dtype = decode_via_dtype - self._emit_decode_timedelta_future_warning = False def encode(self, variable: Variable, name: T_Name = None) -> Variable: if np.issubdtype(variable.dtype, np.timedelta64): @@ -1540,23 +1539,6 @@ def decode(self, variable: Variable, name: T_Name = None) -> Variable: else: time_unit = self.time_unit else: - if self._emit_decode_timedelta_future_warning: - var_string = f"the variable {name!r}" if name else "" - emit_user_level_warning( - "In a future version, xarray will not decode " - f"{var_string} into a timedelta64 dtype based on the " - "presence of a timedelta-like 'units' attribute by " - "default. Instead it will rely on the presence of a " - "timedelta64 'dtype' attribute, which is now xarray's " - "default way of encoding timedelta64 values.\n" - "To continue decoding into a timedelta64 dtype, either " - "set `decode_timedelta=True` when opening this " - "dataset, or add the attribute " - "`dtype='timedelta64[ns]'` to this variable on disk.\n" - "To opt-in to future behavior, set " - "`decode_timedelta=False`.", - FutureWarning, - ) if self.time_unit is None: time_unit = "ns" else: diff --git a/xarray/conventions.py b/xarray/conventions.py index 8f804923f30..d3ee05e5da1 100644 --- a/xarray/conventions.py +++ b/xarray/conventions.py @@ -173,12 +173,13 @@ def decode_cf_variable( original_dtype = var.dtype - decode_timedelta_was_none = decode_timedelta is None if decode_timedelta is None: if isinstance(decode_times, CFDatetimeCoder): decode_timedelta = CFTimedeltaCoder(time_unit=decode_times.time_unit) + elif decode_times: + decode_timedelta = CFTimedeltaCoder() else: - decode_timedelta = bool(decode_times) + decode_timedelta = False if concat_characters: if stack_char_dim: @@ -208,9 +209,6 @@ def decode_cf_variable( decode_timedelta = CFTimedeltaCoder( decode_via_units=decode_timedelta, decode_via_dtype=decode_timedelta ) - decode_timedelta._emit_decode_timedelta_future_warning = ( - decode_timedelta_was_none - ) var = decode_timedelta.decode(var, name=name) if decode_times: # remove checks after end of deprecation cycle diff --git a/xarray/tests/test_coding_times.py b/xarray/tests/test_coding_times.py index 1bb533dc7d8..15e38233c4b 100644 --- a/xarray/tests/test_coding_times.py +++ b/xarray/tests/test_coding_times.py @@ -1825,67 +1825,49 @@ def test_encode_cf_timedelta_small_dtype_missing_value(use_dask) -> None: _DECODE_TIMEDELTA_VIA_UNITS_TESTS = { - "default": (True, None, np.dtype("timedelta64[ns]"), True), - "decode_timedelta=True": (True, True, np.dtype("timedelta64[ns]"), False), - "decode_timedelta=False": (True, False, np.dtype("int64"), False), - "inherit-time_unit-from-decode_times": ( - CFDatetimeCoder(time_unit="s"), - None, - np.dtype("timedelta64[s]"), - True, - ), + "default": (True, None, np.dtype("int64")), + "decode_timedelta=True": (True, True, np.dtype("timedelta64[ns]")), + "decode_timedelta=False": (True, False, np.dtype("int64")), "set-time_unit-via-CFTimedeltaCoder-decode_times=True": ( True, - CFTimedeltaCoder(time_unit="s"), + CFTimedeltaCoder(decode_via_units=True, time_unit="s"), np.dtype("timedelta64[s]"), - False, ), "set-time_unit-via-CFTimedeltaCoder-decode_times=False": ( False, - CFTimedeltaCoder(time_unit="s"), + CFTimedeltaCoder(decode_via_units=True, time_unit="s"), np.dtype("timedelta64[s]"), - False, ), "override-time_unit-from-decode_times": ( CFDatetimeCoder(time_unit="ns"), - CFTimedeltaCoder(time_unit="s"), + CFTimedeltaCoder(decode_via_units=True, time_unit="s"), np.dtype("timedelta64[s]"), - False, ), } @pytest.mark.parametrize( - ("decode_times", "decode_timedelta", "expected_dtype", "warns"), + ("decode_times", "decode_timedelta", "expected_dtype"), list(_DECODE_TIMEDELTA_VIA_UNITS_TESTS.values()), ids=list(_DECODE_TIMEDELTA_VIA_UNITS_TESTS.keys()), ) def test_decode_timedelta_via_units( - decode_times, decode_timedelta, expected_dtype, warns + decode_times, decode_timedelta, expected_dtype ) -> None: timedeltas = pd.timedelta_range(0, freq="D", periods=3) attrs = {"units": "days"} var = Variable(["time"], timedeltas, encoding=attrs) encoded = Variable(["time"], np.array([0, 1, 2]), attrs=attrs) - if warns: - with pytest.warns( - FutureWarning, - match="xarray will not decode the variable 'foo' into a timedelta64 dtype", - ): - decoded = conventions.decode_cf_variable( - "foo", - encoded, - decode_times=decode_times, - decode_timedelta=decode_timedelta, - ) + decoded = conventions.decode_cf_variable( + "foo", encoded, decode_times=decode_times, decode_timedelta=decode_timedelta + ) + if decode_timedelta is True or ( + isinstance(decode_timedelta, CFTimedeltaCoder) + and decode_timedelta.decode_via_units + ): + assert_equal(var, decoded) else: - decoded = conventions.decode_cf_variable( - "foo", encoded, decode_times=decode_times, decode_timedelta=decode_timedelta - ) - if decode_timedelta is False: assert_equal(encoded, decoded) - else: - assert_equal(var, decoded) assert decoded.dtype == expected_dtype @@ -1954,7 +1936,7 @@ def test_decode_timedelta_via_dtype( @pytest.mark.parametrize("dtype", [np.uint64, np.int64, np.float64]) def test_decode_timedelta_dtypes(dtype) -> None: encoded = Variable(["time"], np.arange(10), {"units": "seconds"}) - coder = CFTimedeltaCoder(time_unit="s") + coder = CFTimedeltaCoder(decode_via_units=True, time_unit="s") decoded = coder.decode(encoded) assert decoded.dtype.kind == "m" assert_equal(coder.encode(decoded), encoded) @@ -1963,8 +1945,9 @@ def test_decode_timedelta_dtypes(dtype) -> None: def test_lazy_decode_timedelta_unexpected_dtype() -> None: attrs = {"units": "seconds"} encoded = Variable(["time"], [0, 0.5, 1], attrs=attrs) + decode_timedelta = CFTimedeltaCoder(decode_via_units=True, time_unit="s") decoded = conventions.decode_cf_variable( - "foo", encoded, decode_timedelta=CFTimedeltaCoder(time_unit="s") + "foo", encoded, decode_timedelta=decode_timedelta ) expected_dtype_upon_lazy_decoding = np.dtype("timedelta64[s]") @@ -1978,8 +1961,9 @@ def test_lazy_decode_timedelta_unexpected_dtype() -> None: def test_lazy_decode_timedelta_error() -> None: attrs = {"units": "seconds"} encoded = Variable(["time"], [0, np.iinfo(np.int64).max, 1], attrs=attrs) + decode_timedelta = CFTimedeltaCoder(decode_via_units=True, time_unit="ms") decoded = conventions.decode_cf_variable( - "foo", encoded, decode_timedelta=CFTimedeltaCoder(time_unit="ms") + "foo", encoded, decode_timedelta=decode_timedelta ) with pytest.raises(OutOfBoundsTimedelta, match="overflow"): decoded.load() diff --git a/xarray/tests/test_conventions.py b/xarray/tests/test_conventions.py index 461bd727e64..38b835fd3d5 100644 --- a/xarray/tests/test_conventions.py +++ b/xarray/tests/test_conventions.py @@ -553,7 +553,9 @@ def test_decode_cf_time_kwargs(self, time_unit) -> None: dsc = conventions.decode_cf( ds, decode_times=CFDatetimeCoder(time_unit=time_unit), - decode_timedelta=CFTimedeltaCoder(time_unit=time_unit), + decode_timedelta=CFTimedeltaCoder( + decode_via_units=True, time_unit=time_unit + ), ) assert dsc.timedelta.dtype == np.dtype(f"m8[{time_unit}]") assert dsc.time.dtype == np.dtype(f"M8[{time_unit}]") @@ -679,15 +681,3 @@ def test_encode_cf_variable_with_vlen_dtype() -> None: encoded_v = conventions.encode_cf_variable(v) assert encoded_v.data.dtype.kind == "O" assert coding.strings.check_vlen_dtype(encoded_v.data.dtype) is str - - -def test_decode_cf_variables_decode_timedelta_warning() -> None: - v = Variable(["time"], [1, 2], attrs={"units": "seconds"}) - variables = {"a": v} - - with warnings.catch_warnings(): - warnings.filterwarnings("error", "decode_timedelta", FutureWarning) - conventions.decode_cf_variables(variables, {}, decode_timedelta=True) - - with pytest.warns(FutureWarning, match="decode_timedelta"): - conventions.decode_cf_variables(variables, {})