Any way to specify on disk that an array with timedelta-like units should not be decoded as timedelta?

[delgadom]

I'm working with xclim climate indicators, and trying to stick to their specification as much as possible to maintain compatibility with their operations, which can be chained. I know this sounds like an xclim question, but it's actually an xarray encoding question I promise 😄

Some of their transforms (e.g. see xclim.indicators.atmos.tx_days_above) return DataArrays with units 'days'. This is a totally reasonable description of the data contents - count of days above a threshold - but it also falls into the list of units keywords parsed by the decode_timedeltas argument of xr.open_zarr. So, if I save a dataset containing this variable to zarr (or to netcdf or any other backend as far as I can tell), the data will be read in as a timedelta unless I specify decode_timedeltas=False.

My question is - is there any way to specify in the encodings attribute or in saving the array that the data should not be interpreted as timedelta? I've tried explicitly setting the encoding attribute ds.tx_days_above.encoding = {'units': 'day_counts'} or some other unit but this doesn't help. The only way I can get the array to be read correctly using the default behavior is by changing the units in the attributes.

Here's an MRE:

In [1]: import xarray as xr, pandas as pd, numpy as np

In [2]: ds = xr.Dataset({
   ...:     'tx_days_above': xr.DataArray(
   ...:         np.random.randint(0, 365, size=(20, 20, 20)).astype('float32'),
   ...:         dims=['x', 'y', 'year'],
   ...:         coords=[range(23, 43), range(-110, -90), range(2020, 2040)],
   ...:         attrs={'units': 'days'},
   ...:     ),
   ...: })

In [3]: ds
Out[3]:
<xarray.Dataset>
Dimensions:        (x: 20, y: 20, year: 20)
Coordinates:
  * x              (x) int64 23 24 25 26 27 28 29 30 ... 35 36 37 38 39 40 41 42
  * y              (y) int64 -110 -109 -108 -107 -106 ... -95 -94 -93 -92 -91
  * year           (year) int64 2020 2021 2022 2023 2024 ... 2036 2037 2038 2039
Data variables:
    tx_days_above  (x, y, year) float32 241.0 264.0 202.0 ... 274.0 272.0 133.0

In [4]: ds.to_zarr('days_above.zarr')
Out[4]: <xarray.backends.zarr.ZarrStore at 0x10f8ba2d0>

In [5]: ds_zarr = xr.open_zarr('days_above.zarr').load()

In [6]: ds_zarr
Out[6]:
<xarray.Dataset>
Dimensions:        (x: 20, y: 20, year: 20)
Coordinates:
  * x              (x) int64 23 24 25 26 27 28 29 30 ... 35 36 37 38 39 40 41 42
  * y              (y) int64 -110 -109 -108 -107 -106 ... -95 -94 -93 -92 -91
  * year           (year) int64 2020 2021 2022 2023 2024 ... 2036 2037 2038 2039
Data variables:
    tx_days_above  (x, y, year) timedelta64[ns] 241 days 264 days ... 133 days

[spencerkclark]

Indeed there are no solutions outside the workarounds you mentioned currently. This is something that has come up in the past -- see discussion in #1621 -- and spawned the addition of the decode_timedelta argument. A possible solution you might be interested in, however, is the one @shoyer recommended in #1621 (comment). Would you be up for implementing that approach?

[delgadom]

thanks for the archival research! 🥇 I didn't see any issues or discussion items on this but I guess I didn't go back far enough 😆 will definitely check it out!