Repository: xCDAT/xcdat

Issue:
[Feature]: Improve support for user generated time coordinates and temporal operations
### Is your feature request related to a problem?

Related to https://github.com/xCDAT/xcdat/discussions/391#discussioncomment-4828490

Generating time coordinates using pandas/xarray and using them with xCDAT temporal APIs seem to be a desired functionality.
We should look more closely at how we can better support this workflow (e.g., better documentation, alternative solutions, etc.). I'll open up an issue and link it here.

> At some point, we should investigate whether the first KeyError is needed. I think this was probably used to figure out if the time axis could be used for this computation, but the calendar attribute probably isn't actually required. So the logic could possibly be updated here.

When decoding time, [xarray.decode_cf_datetime](https://github.com/pydata/xarray/blob/b4515582ffc8b7f63632bfccd109d19889d00384/xarray/coding/times.py#L281-L283) defaults `calendar = "standard"` if the `calendar` encoding attribute is not set.

We also adopted this behavior in `xcdat.decode_time`: https://github.com/xCDAT/xcdat/blob/16a807cb8c885700efc85d3223b8392f3627bbc3/xcdat/dataset.py#L318-L324

We can consider doing the same thing in the temporal APIs and raise a warning for a missing calendar attribute (instead of `KeyError`).

_Originally posted by @tomvothecoder in https://github.com/xCDAT/xcdat/discussions/391#discussioncomment-4832383_
    

### Describe the solution you'd like

_**Work in progress**_

1. Raise a warning instead of a `KeyError` if `calendar` attribute is missing
  - https://github.com/xCDAT/xcdat/blob/16a807cb8c885700efc85d3223b8392f3627bbc3/xcdat/temporal.py#L731-L740
  - https://github.com/xCDAT/xcdat/blob/16a807cb8c885700efc85d3223b8392f3627bbc3/xcdat/temporal.py#L889-L896
2. Add explanation on why leap days are dropped for high frequency (`day`) climatology and departure calculations
  - https://github.com/xCDAT/xcdat/blob/16a807cb8c885700efc85d3223b8392f3627bbc3/xcdat/temporal.py#L513-L515
  - https://github.com/xCDAT/xcdat/blob/16a807cb8c885700efc85d3223b8392f3627bbc3/xcdat/temporal.py#L513-L515
3. Improve the `KeyError` raised when getting bounds -- point user to use `xcdat.add_missing_bounds` or `xcdat.add_bounds`
  - https://github.com/xCDAT/xcdat/blob/16a807cb8c885700efc85d3223b8392f3627bbc3/xcdat/bounds.py#L217-L222

### Describe alternatives you've considered

_No response_

### Additional context

_No response_

Discussion:
codecov-commenter: # [Codecov](https://codecov.io/gh/xCDAT/xcdat/pull/415?src=pr&el=h1&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT) Report
Base: **100.00**% // Head: **100.00**% // No change to project coverage :thumbsup:
> Coverage data is based on head [(`9911ad1`)](https://codecov.io/gh/xCDAT/xcdat/pull/415?src=pr&el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT) compared to base [(`c1989c9`)](https://codecov.io/gh/xCDAT/xcdat/commit/c1989c950b0362f55636d7d276af2083e97dcdb7?el=desc&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT).
> Patch coverage: 100.00% of modified lines in pull request are covered.

:mega: This organization is not using Codecov’s [GitHub App Integration](https://github.com/apps/codecov). We recommend you install it so Codecov can continue to function properly for your repositories. [Learn more](https://about.codecov.io/blog/codecov-is-updating-its-github-integration/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT)

<details><summary>Additional details and impacted files</summary>


```diff
@@            Coverage Diff            @@
##              main      #415   +/-   ##
=========================================
  Coverage   100.00%   100.00%           
=========================================
  Files           14        14           
  Lines         1260      1266    +6     
=========================================
+ Hits          1260      1266    +6     
```


| [Impacted Files](https://codecov.io/gh/xCDAT/xcdat/pull/415?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT) | Coverage Δ | |
|---|---|---|
| [xcdat/bounds.py](https://codecov.io/gh/xCDAT/xcdat/pull/415?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT#diff-eGNkYXQvYm91bmRzLnB5) | `100.00% <ø> (ø)` | |
| [xcdat/temporal.py](https://codecov.io/gh/xCDAT/xcdat/pull/415?src=pr&el=tree&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT#diff-eGNkYXQvdGVtcG9yYWwucHk=) | `100.00% <100.00%> (ø)` | |

Help us with your feedback. Take ten seconds to tell us [how you rate us](https://about.codecov.io/nps?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT). Have a feature suggestion? [Share it here.](https://app.codecov.io/gh/feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT)

</details>

[:umbrella: View full report at Codecov](https://codecov.io/gh/xCDAT/xcdat/pull/415?src=pr&el=continue&utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT).   
:loudspeaker: Do you have feedback about the report comment? [Let us know in this issue](https://about.codecov.io/codecov-pr-comment-feedback/?utm_medium=referral&utm_source=github&utm_content=comment&utm_campaign=pr+comments&utm_term=xCDAT).

tomvothecoder: > This looks good to me (and solves the initial issue / discussion item).
> 
> The only thing I got caught up on was the changes to `_set_data_var_attrs()`. Do you think the top line documentation for this function should be augmented, e.g.,
> 
> Set data variable metadata as object attributes **and checks whether the time axis is decoded**.
> 
> It seems like this function is partially used to enforce this (prior to calling key temporal methods).

Thanks for that catch. I just updated the docstring based on your feedback.

I was considering splitting out the check for decoded time as a separate method to keep the implementation cleaner. However, it requires looking at the data variable values which means it's not a big deal to keep it inside `_set_data_var_attrs()`.


Patch:
diff --git a/xcdat/bounds.py b/xcdat/bounds.py
index 7b5aff1..f2b772c 100644
--- a/xcdat/bounds.py
+++ b/xcdat/bounds.py
@@ -216,9 +216,11 @@ class BoundsAccessor:
 
         if len(bounds_keys) == 0:
             raise KeyError(
-                f"No bounds were found for the '{axis}' axis. Make sure bounds vars "
-                "exist in the Dataset with names that match the 'bounds' keys, or try "
-                "adding bounds."
+                f"No bounds data variables were found for the '{axis}' axis. Make sure "
+                "the dataset has bound data vars and their names match the 'bounds' "
+                "attributes found on their related time coordinate variables. "
+                "Alternatively, you can add bounds with `xcdat.add_missing_bounds` "
+                "or `xcdat.add_bounds`."
             )
 
         bounds: Union[xr.Dataset, xr.DataArray] = self._dataset[
diff --git a/xcdat/temporal.py b/xcdat/temporal.py
index 7d06bcd..f4f0cc7 100644
--- a/xcdat/temporal.py
+++ b/xcdat/temporal.py
@@ -3,6 +3,7 @@ from itertools import chain
 from typing import Dict, List, Literal, Optional, Tuple, TypedDict, get_args
 
 import cf_xarray  # noqa: F401
+import cftime
 import numpy as np
 import pandas as pd
 import xarray as xr
@@ -354,8 +355,11 @@ class TemporalAccessor:
             * "season": groups by season for the seasonal cycle climatology.
             * "month": groups by month for the annual cycle climatology.
             * "day": groups by (month, day) for the daily cycle climatology.
-              Leap days (if present) are dropped if the CF calendar type is
-              ``"gregorian"``, ``"proleptic_gregorian"``, or ``"standard"``
+              If the CF calendar type is ``"gregorian"``,
+              ``"proleptic_gregorian"``, or ``"standard"``, leap days (if
+              present) are dropped to avoid inconsistencies when calculating
+              climatologies. Refer to [1]_ for more details on this
+              implementation decision.
 
         weighted : bool, optional
             Calculate averages using weights, by default True.
@@ -416,6 +420,10 @@ class TemporalAccessor:
         xr.Dataset
             Dataset with the climatology of a data variable.
 
+        References
+        ----------
+        .. [1] https://github.com/xCDAT/xcdat/discussions/332
+
         Examples
         --------
 
@@ -510,9 +518,12 @@ class TemporalAccessor:
 
             * "season": groups by season for the seasonal cycle departures.
             * "month": groups by month for the annual cycle departures.
-            * "day": groups by (month, day) for the daily cycle departures. Leap
-              days (if present) are dropped if the CF calendar type is
-              ``"gregorian"``, ``"proleptic_gregorian"``, or ``"standard"``
+            * "day": groups by (month, day) for the daily cycle departures.
+              If the CF calendar type is ``"gregorian"``,
+              ``"proleptic_gregorian"``, or ``"standard"``, leap days (if
+              present) are dropped to avoid inconsistencies when calculating
+              climatologies. Refer to [2]_ for more details on this
+              implementation decision.
 
         weighted : bool, optional
             Calculate averages using weights, by default True.
@@ -575,11 +586,12 @@ class TemporalAccessor:
 
         Notes
         -----
-        Refer to [1]_ to learn more about how xarray's grouped arithmetic works.
+        Refer to [3]_ to learn more about how xarray's grouped arithmetic works.
 
         References
         ----------
-        .. [1] https://xarray.pydata.org/en/stable/user-guide/groupby.html#grouped-arithmetic
+        .. [2] https://github.com/xCDAT/xcdat/discussions/332
+        .. [3] https://xarray.pydata.org/en/stable/user-guide/groupby.html#grouped-arithmetic
 
         Examples
         --------
@@ -708,7 +720,9 @@ class TemporalAccessor:
         return ds
 
     def _set_data_var_attrs(self, data_var: str):
-        """Set data variable metadata as object attributes.
+        """
+        Set data variable metadata as object attributes and checks whether the
+        time axis is decoded.
 
         This includes the name of the data variable, the time axis dimension
         name, the calendar type and its corresponding cftime object (date type).
@@ -728,17 +742,34 @@ class TemporalAccessor:
         self.data_var = data_var
         self.dim = get_dim_coords(dv, "T").name
 
+        first_time_coord = dv[self.dim].values[0]
+        is_decoded = isinstance(first_time_coord, cftime.datetime) or isinstance(
+            first_time_coord, np.datetime64
+        )
+        if not is_decoded:
+            raise ValueError(
+                f"The time coordinates type is {type(first_time_coord)}, not "
+                "'cftime.datetime' or 'np.datetime64'. Time coordinates must be "
+                "decoded as datetime before using TemporalAccessor methods. "
+            )
+
+        # Get the `cftime` date type based on the CF calendar attribute.
+        # The date type is used to get the correct cftime.datetime sub-class
+        # type for creating new grouped time coordinates for averaging.
         try:
             self.calendar = dv[self.dim].encoding["calendar"]
-            self.date_type = get_date_type(self.calendar)
         except KeyError:
-            raise KeyError(
-                f"The 'calendar' encoding attribute is not set on the '{data_var}' "
-                f"time coordinate variable ({self.dim}). This might indicate that the "
-                "time coordinates were not decoded, which is required for temporal "
-                "averaging operations. "
+            self.calendar = "standard"
+            logger.warning(
+                f"'{self.dim}' does not have a calendar encoding attribute set, "
+                "which is used to determine the `cftime.datetime` object type for the "
+                "output time coordinates. Defaulting to CF 'standard' calendar. "
+                "Otherwise, set the calendar type (e.g., "
+                "ds['time'].encoding['calendar'] = 'noleap') and try again."
             )
 
+        self.date_type = get_date_type(self.calendar)
+
     def _set_arg_attrs(
         self,
         mode: Mode,
@@ -1084,11 +1115,11 @@ class TemporalAccessor:
 
         Notes
         -----
-        Refer to [5]_ for the supported CF convention calendar types.
+        Refer to [4]_ for the supported CF convention calendar types.
 
         References
         ----------
-        .. [5] https://cfconventions.org/cf-conventions/cf-conventions.html#calendar
+        .. [4] https://cfconventions.org/cf-conventions/cf-conventions.html#calendar
         """
         with xr.set_options(keep_attrs=True):
             time_lengths: xr.DataArray = time_bounds[:, 1] - time_bounds[:, 0]
@@ -1228,11 +1259,11 @@ class TemporalAccessor:
 
         Notes
         -----
-        Refer to [2]_ for information on xarray datetime accessor components.
+        Refer to [5]_ for information on xarray datetime accessor components.
 
         References
         ----------
-        .. [2] https://xarray.pydata.org/en/stable/user-guide/time-series.html#datetime-components
+        .. [5] https://xarray.pydata.org/en/stable/user-guide/time-series.html#datetime-components
         """
         df = pd.DataFrame()
 
@@ -1457,14 +1488,13 @@ class TemporalAccessor:
 
         Notes
         -----
-        Refer to [3]_ and [4]_ for more information on Timestamp-valid range.
+        Refer to [6]_ and [7]_ for more information on Timestamp-valid range.
         We use cftime.datetime objects to avoid these time range issues.
 
         References
         ----------
-        .. [3] https://docs.xarray.dev/en/stable/user-guide/weather-climate.html#non-standard-calendars-and-dates-outside-the-timestamp-valid-range
-
-        .. [4] https://pandas.pydata.org/pandas-docs/stable/user_guide/timeseries.html#timestamp-limitations
+        .. [6] https://docs.xarray.dev/en/stable/user-guide/weather-climate.html#non-standard-calendars-and-dates-outside-the-timestamp-valid-range
+        .. [7] https://pandas.pydata.org/pandas-docs/stable/user_guide/timeseries.html#timestamp-limitations
         """
         df_new = df.copy()
 


Test patch:
diff --git a/tests/test_temporal.py b/tests/test_temporal.py
index fc25ef3..cf48116 100644
--- a/tests/test_temporal.py
+++ b/tests/test_temporal.py
@@ -6,8 +6,11 @@ from xarray.coding.cftime_offsets import get_date_type
 from xarray.tests import requires_dask
 
 from tests.fixtures import generate_dataset
+from xcdat.logger import setup_custom_logger
 from xcdat.temporal import TemporalAccessor
 
+logger = setup_custom_logger("xcdat.temporal", propagate=True)
+
 
 class TestTemporalAccessor:
     def test__init__(self):
@@ -26,16 +29,31 @@ class TestTemporalAccessor:
 
 
 class TestAverage:
-    def test_raises_error_if_calendar_encoding_attr_not_found_on_data_var_time_coords(
-        self,
+    def test_raises_error_if_time_coords_are_not_decoded(self):
+        ds: xr.Dataset = generate_dataset(
+            decode_times=False, cf_compliant=False, has_bounds=True
+        )
+
+        with pytest.raises(ValueError):
+            ds.temporal.average("ts")
+
+    def test_raises_warning_if_calendar_encoding_attr_not_found_on_data_var_time_coords(
+        self, caplog
     ):
         ds: xr.Dataset = generate_dataset(
             decode_times=True, cf_compliant=False, has_bounds=True
         )
         ds.ts.time.encoding = {}
 
-        with pytest.raises(KeyError):
-            ds.temporal.average("ts")
+        ds.temporal.average("ts")
+
+        assert (
+            "'time' does not have a calendar encoding attribute set, "
+            "which is used to determine the `cftime.datetime` object type for the "
+            "output time coordinates. Defaulting to CF 'standard' calendar. "
+            "Otherwise, set the calendar type (e.g., "
+            "ds['time'].encoding['calendar'] = 'noleap') and try again."
+        ) in caplog.text
 
     def test_averages_for_yearly_time_series(self):
         ds = xr.Dataset(
@@ -430,16 +448,31 @@ class TestGroupAverage:
             attrs={"test_attr": "test"},
         )
 
-    def test_raises_error_if_calendar_encoding_attr_not_found_on_data_var_time_coords(
-        self,
+    def test_raises_error_if_time_coords_are_not_decoded(self):
+        ds: xr.Dataset = generate_dataset(
+            decode_times=False, cf_compliant=False, has_bounds=True
+        )
+
+        with pytest.raises(ValueError):
+            ds.temporal.group_average("ts", freq="year")
+
+    def test_raises_warning_if_calendar_encoding_attr_not_found_on_data_var_time_coords(
+        self, caplog
     ):
         ds: xr.Dataset = generate_dataset(
             decode_times=True, cf_compliant=False, has_bounds=True
         )
         ds.ts.time.encoding = {}
 
-        with pytest.raises(KeyError):
-            ds.temporal.group_average("ts", freq="year")
+        ds.temporal.group_average("ts", freq="year")
+
+        assert (
+            "'time' does not have a calendar encoding attribute set, "
+            "which is used to determine the `cftime.datetime` object type for the "
+            "output time coordinates. Defaulting to CF 'standard' calendar. "
+            "Otherwise, set the calendar type (e.g., "
+            "ds['time'].encoding['calendar'] = 'noleap') and try again."
+        ) in caplog.text
 
     def test_weighted_annual_averages(self):
         ds = self.ds.copy()
@@ -994,16 +1027,31 @@ class TestClimatology:
             decode_times=True, cf_compliant=False, has_bounds=True
         )
 
-    def test_raises_error_if_calendar_encoding_attr_not_found_on_data_var_time_coords(
-        self,
+    def test_raises_error_if_time_coords_are_not_decoded(self):
+        ds: xr.Dataset = generate_dataset(
+            decode_times=False, cf_compliant=False, has_bounds=True
+        )
+
+        with pytest.raises(ValueError):
+            ds.temporal.climatology("ts", freq="year")
+
+    def test_raises_warning_if_calendar_encoding_attr_not_found_on_data_var_time_coords(
+        self, caplog
     ):
         ds: xr.Dataset = generate_dataset(
             decode_times=True, cf_compliant=False, has_bounds=True
         )
         ds.ts.time.encoding = {}
 
-        with pytest.raises(KeyError):
-            ds.temporal.climatology("ts", freq="season")
+        ds.temporal.climatology("ts", freq="season")
+
+        assert (
+            "'time' does not have a calendar encoding attribute set, "
+            "which is used to determine the `cftime.datetime` object type for the "
+            "output time coordinates. Defaulting to CF 'standard' calendar. "
+            "Otherwise, set the calendar type (e.g., "
+            "ds['time'].encoding['calendar'] = 'noleap') and try again."
+        ) in caplog.text
 
     def test_weighted_seasonal_climatology_with_DJF(self):
         ds = self.ds.copy()
@@ -1580,16 +1628,31 @@ class TestDepartures:
 
         self.seasons = ["JJA", "MAM", "SON", "DJF"]
 
-    def test_raises_error_if_calendar_encoding_attr_not_found_on_data_var_time_coords(
-        self,
+    def test_raises_error_if_time_coords_are_not_decoded(self):
+        ds: xr.Dataset = generate_dataset(
+            decode_times=False, cf_compliant=False, has_bounds=True
+        )
+
+        with pytest.raises(ValueError):
+            ds.temporal.departures("ts", freq="season")
+
+    def test_raises_warning_if_calendar_encoding_attr_not_found_on_data_var_time_coords(
+        self, caplog
     ):
         ds: xr.Dataset = generate_dataset(
             decode_times=True, cf_compliant=False, has_bounds=True
         )
         ds.ts.time.encoding = {}
 
-        with pytest.raises(KeyError):
-            ds.temporal.departures("ts", freq="year")
+        ds.temporal.departures("ts", freq="season")
+
+        assert (
+            "'time' does not have a calendar encoding attribute set, "
+            "which is used to determine the `cftime.datetime` object type for the "
+            "output time coordinates. Defaulting to CF 'standard' calendar. "
+            "Otherwise, set the calendar type (e.g., "
+            "ds['time'].encoding['calendar'] = 'noleap') and try again."
+        ) in caplog.text
 
     def test_weighted_seasonal_departures_with_DJF(self):
         ds = self.ds.copy()