Template for unit tests

environment-dev.yml

@@ -26,7 +26,7 @@ dependencies:
   - rechunker
   - shapely
   - xarray
-  - xclim >=0.37
+  - xclim >=0.43  # Harder requirement because of xclim.testing.helpers.test_timeseries introduced in 0.43
   - xesmf >=0.7
   - zarr
   # To avoid bugs

tests/test_aggregate.py

@@ -0,0 +1,152 @@
+import cftime
+import numpy as np
+import pytest
+import xarray as xr
+from xclim.core.calendar import convert_calendar
+from xclim.testing.helpers import test_timeseries as timeseries
+
+import xscen as xs
+
+
+class TestClimatologicalMean:
+    def test_daily(self):
+        ds = timeseries(
+            np.tile(np.arange(1, 13), 3),
+            variable="tas",
+            start="2001-01-01",
+            freq="D",
+            as_dataset=True,
+        )
+        with pytest.raises(NotImplementedError):
+            xs.climatological_mean(ds)
+
+    def test_all_default(self):
+        for xrfreq in ["MS", "AS-JAN"]:
+            o = 12 if xrfreq == "MS" else 1
+
+            ds = timeseries(
+                np.tile(np.arange(1, o + 1), 30),
+                variable="tas",
+                start="2001-01-01",
+                freq=xrfreq,
+                as_dataset=True,
+            )
+            out = xs.climatological_mean(ds)
+
+            np.testing.assert_array_equal(out.tas, np.arange(1, o + 1))
+            np.testing.assert_array_equal(
+                len(out.time), o * len(np.unique(out.horizon.values))
+            )
+            np.testing.assert_array_equal(out.time[0], ds.time[0])
+            assert (
+                out.tas.attrs["description"]
+                == f"30-year mean of {ds.tas.attrs['description']}"
+            )
+            assert (
+                "30-year rolling average (non-centered) with a minimum of 30 years of data"
+                in out.tas.attrs["history"]
+            )
+            assert all(out.horizon.values == "2001-2030")
+            assert out.attrs["cat:processing_level"] == "climatology"
+
+    def test_options(self):
+        for xrfreq in ["MS", "AS-JAN"]:
+            o = 12 if xrfreq == "MS" else 1
+
+            ds = timeseries(
+                np.tile(np.arange(1, o + 1), 30),
+                variable="tas",
+                start="2001-01-01",
+                freq=xrfreq,
+                as_dataset=True,
+            )
+            out = xs.climatological_mean(
+                ds, window=15, interval=5, to_level="for_testing"
+            )
+
+            np.testing.assert_array_equal(
+                out.tas,
+                np.tile(np.arange(1, o + 1), len(np.unique(out.horizon.values))),
+            )
+            np.testing.assert_array_equal(
+                len(out.time), o * len(np.unique(out.horizon.values))
+            )
+            np.testing.assert_array_equal(out.time[0], ds.time[0])
+            assert (
+                out.tas.attrs["description"]
+                == f"15-year mean of {ds.tas.attrs['description']}"
+            )
+            assert (
+                "15-year rolling average (non-centered) with a minimum of 15 years of data"
+                in out.tas.attrs["history"]
+            )
+            assert all(
+                h in np.unique(out.horizon.values)
+                for h in ["2001-2015", "2006-2020", "2011-2025", "2016-2030"]
+            )
+            assert out.attrs["cat:processing_level"] == "for_testing"
+
+    def test_minperiods(self):
+        ds = timeseries(
+            np.tile(np.arange(1, 5), 30),
+            variable="tas",
+            start="2000-12-01",
+            freq="QS-DEC",
+            as_dataset=True,
+        )
+        ds = ds.where(ds.time.dt.year >= 2001)
+        out = xs.climatological_mean(ds, window=30)
+        assert all(np.isreal(out.tas))
+        assert len(out.time) == 8
+        out = xs.climatological_mean(ds, window=30, min_periods=30)
+        assert np.sum(np.isreal(out.tas))
+
+        with pytest.raises(ValueError):
+            xs.climatological_mean(ds, window=5, min_periods=6)
+
+    def test_periods(self):
+        ds1 = timeseries(
+            np.tile(np.arange(1, 2), 10),
+            variable="tas",
+            start="2001-01-01",
+            freq="AS-JAN",
+            as_dataset=True,
+        )
+        ds2 = timeseries(
+            np.tile(np.arange(1, 2), 10),
+            variable="tas",
+            start="2021-01-01",
+            freq="AS-JAN",
+            as_dataset=True,
+        )
+        ds = xr.concat([ds1, ds2], dim="time")
+        with pytest.raises(ValueError):
+            xs.climatological_mean(ds)
+
+        out = xs.climatological_mean(ds, periods=[["2001", "2010"], ["2021", "2030"]])
+        np.testing.assert_array_equal(len(out.time), 2)
+        assert all(
+            h in np.unique(out.horizon.values) for h in ["2001-2010", "2021-2030"]
+        )
+
+    def test_calendars(self):
+        ds = timeseries(
+            np.tile(np.arange(1, 2), 30),
+            variable="tas",
+            start="2001-01-01",
+            freq="AS-JAN",
+            as_dataset=True,
+        )
+
+        assert isinstance(
+            xs.climatological_mean(convert_calendar(ds, "standard")),
+            cftime.Datetime360Day,
+        )
+        assert isinstance(
+            xs.climatological_mean(convert_calendar(ds, "noleap")),
+            cftime.DatetimeNoLeap,
+        )
+        assert isinstance(
+            xs.climatological_mean(convert_calendar(ds, "360_day", align_on="random")),
+            cftime.Datetime360Day,
+        )

xscen/aggregate.py

@@ -73,6 +73,11 @@ def climatological_mean(
         Returns a Dataset of the climatological mean
 
     """
+    if xr.infer_freq(ds.time) == "D":
+        raise NotImplementedError(
+            "xs.climatological_mean does not currently support daily data."
+        )
+
     # there is one less occurrence when a period crosses years
     freq_across_year = [
         f"{f}-{mon}"
@@ -100,9 +105,17 @@ def climatological_mean(
 
     window = window or int(periods[0][1]) - int(periods[0][0]) + 1
 
-    if ds.attrs.get("cat:xrfreq") in freq_across_year and min_periods is None:
+    if (
+        any(
+            x in freq_across_year
+            for x in [ds.attrs.get("cat:xrfreq"), xr.infer_freq(ds.time)]
+        )
+        and min_periods is None
+    ):
         min_periods = window - 1
     min_periods = min_periods or window
+    if min_periods > window:
+        raise ValueError("'min_periods' should be smaller or equal to 'window'")
 
     for period in periods:
         # Rolling average
@@ -173,9 +186,7 @@ def climatological_mean(
             else new_history
         )
         ds_rolling[vv].attrs["history"] = history
-
-    if to_level is not None:
-        ds_rolling.attrs["cat:processing_level"] = to_level
+    ds_rolling.attrs["cat:processing_level"] = to_level
 
     return ds_rolling
 
@@ -233,6 +244,11 @@ def compute_deltas(
         )
 
     if "time" in ds:
+        if xr.infer_freq(ds.time) == "D":
+            raise NotImplementedError(
+                "xs.climatological_mean does not currently support daily data."
+            )
+
         # Remove references to 'year' in REF
         ind = pd.MultiIndex.from_arrays(
             [ref.time.dt.month.values, ref.time.dt.day.values], names=["month", "day"]