issues

Like pandas.date_range, cftime_range would previously return dates outside the range of the specified start and end dates if provided a negative frequency: ```

start = cftime.DatetimeGregorian(2023, 10, 31) end = cftime.DatetimeGregorian(2021, 11, 1) xr.cftime_range(start, end, freq="-1YE") CFTimeIndex([2023-12-31 00:00:00, 2022-12-31 00:00:00, 2021-12-31 00:00:00], dtype='object', length=3, calendar='standard', freq='-1YE-DEC') ```

This PR ports a bug fix from pandas (https://github.com/pandas-dev/pandas/issues/56147) to prevent this from happening. The above example now produces: ```

start = cftime.DatetimeGregorian(2023, 10, 31) end = cftime.DatetimeGregorian(2021, 11, 1) xr.cftime_range(start, end, freq="-1YE") CFTimeIndex([2022-12-31 00:00:00, 2021-12-31 00:00:00], dtype='object', length=2, calendar='standard', freq=None) ```

Since this is a bug fix, we do not make any attempt to preserve the old behavior if an earlier version of pandas is installed. In the testing context this means we skip some tests for pandas versions less than 3.0.

• [x] User visible changes (including notable bug fixes) are documented in whats-new.rst

Since cftime version 1.3.0, the base cftime.datetime object can be calendar-aware, obviating the need for calendar-specific subclasses like cftime.DatetimeNoLeap. This PR aims to finally enable the use of these objects in xarray. We can also use this moment to remove cruft around accommodating inexact cftime datetime arithmetic, since that has been fixed since cftime version 1.2.0.

The idea will be to support both for a period of time and eventually drop support for the calendar-specific subclasses. I do not think too much should need to change within xarray—the main challenge will be to see if we can maintain adequate test coverage without multiplying the number of cftime tests by two. This draft PR is at least a start towards that.

• [ ] Closes #4336
• [ ] Closes #4853
• [ ] Closes #5551
• [ ] Closes #8298
• [ ] Closes #8941
• [ ] Tests added
• [ ] User visible changes (including notable bug fixes) are documented in whats-new.rst
• [ ] New functions/methods are listed in api.rst

By default xarray will exactly roundtrip times saved to disk by encoding them using int64 values. However, if a user specifies time encoding units that prevent this, float64 values will be used, and this has the potential to cause roundtripping differences due to roundoff error. Recently, cftime added the ability to encode times using longdouble values (https://github.com/Unidata/cftime/pull/284). On some platforms this offers greater precision than float64 values (though typically not full quad precision). Nevertheless some users might be interested in encoding their times using such values.

The main thing that longdouble values have going for them is that they enable greater precision when using arbitrary units to encode the dates (with int64 we are constrained to using units that allow for time intervals to be expressed with integers). That said, the more I think about this, the more I feel it may not be the best idea:

• Since the meaning of longdouble can vary from platform to platform, I wonder what happens if you encode times using longdouble values on one machine and decode them on another?
• longdouble values cannot be stored with all backends; for example zarr supports it, but netCDF does not.
• We already provide a robust way to exactly roundtrip any dates--i.e. encode them with int64 values--so adding a less robust (if slightly more flexible in terms of units) option might just cause confusion.

It's perhaps still worth opening this issue for discussion in case others have thoughts that might allay those concerns.

cc: @jswhit @dcherian

What happened?

If I construct the following DataArray, and try to chunk its "x" coordinate, I get back a NumPy-backed DataArray: ``` In [2]: a = xr.DataArray([1, 2, 3], dims=["x"], coords=[[4, 5, 6]])

In [3]: a.x.chunk() Out[3]: <xarray.DataArray 'x' (x: 3)> array([4, 5, 6]) Coordinates: * x (x) int64 4 5 6 If I construct a copy of the `"x"` coordinate, things work as I would expect: In [4]: x = xr.DataArray(a.x, dims=a.x.dims, coords=a.x.coords, name="x")

In [5]: x.chunk() Out[5]: <xarray.DataArray 'x' (x: 3)> dask.array<xarray-\<this-array>, shape=(3,), dtype=int64, chunksize=(3,), chunktype=numpy.ndarray> Coordinates: * x (x) int64 4 5 6 ```

What did you expect to happen?

I would expect the following to happen: ``` In [2]: a = xr.DataArray([1, 2, 3], dims=["x"], coords=[[4, 5, 6]])

In [3]: a.x.chunk() Out[3]: <xarray.DataArray 'x' (x: 3)> dask.array<xarray-\<this-array>, shape=(3,), dtype=int64, chunksize=(3,), chunktype=numpy.ndarray> Coordinates: * x (x) int64 4 5 6 ```

Minimal Complete Verifiable Example

No response

Relevant log output

No response

Anything else we need to know?

No response

Environment

INSTALLED VERSIONS

commit: None python: 3.7.10 | packaged by conda-forge | (default, Feb 19 2021, 15:59:12) [Clang 11.0.1 ] python-bits: 64 OS: Darwin OS-release: 21.2.0 machine: x86_64 processor: i386 byteorder: little LC_ALL: None LANG: en_US.UTF-8 LOCALE: ('en_US', 'UTF-8') libhdf5: 1.10.5 libnetcdf: 4.6.3

xarray: 0.20.1 pandas: 1.3.5 numpy: 1.19.4 scipy: 1.5.4 netCDF4: 1.5.5 pydap: None h5netcdf: 0.8.1 h5py: 2.10.0 Nio: None zarr: 2.7.0 cftime: 1.2.1 nc_time_axis: 1.2.0 PseudoNetCDF: None rasterio: None cfgrib: None iris: None bottleneck: None dask: 2.22.0 distributed: None matplotlib: 3.2.2 cartopy: 0.19.0.post1 seaborn: None numbagg: None fsspec: 2021.06.0 cupy: None pint: 0.15 sparse: None setuptools: 49.6.0.post20210108 pip: 20.2.4 conda: 4.10.1 pytest: 6.0.1 IPython: 7.27.0 sphinx: 3.2.1

The discussion of the use of the indexes property in #5102 got me thinking about this StackOverflow answer. For a while I have thought that my answer there isn't very satisfying, not only because it relies on this somewhat obscure indexes property, but also because it only works on dimension coordinates -- i.e. something that would be backed by an index.

Describe the solution you'd like

It would be better if we could do this conversion with astype, e.g. da.astype("datetime64[ns]"). This would allow conversion to datetime64 values for all cftime.datetime DataArrays -- dask-backed or NumPy-backed, 1D or ND -- through a fairly standard and well-known method. To my surprise, while you do not get the nice calendar-switching warning that CFTimeIndex.to_datetimeindex provides, this actually already kind of seems to work (?!):

``` In [1]: import xarray as xr

In [2]: times = xr.cftime_range("2000", periods=6, calendar="noleap")

In [3]: da = xr.DataArray(times.values.reshape((2, 3)), dims=["a", "b"])

In [4]: da.astype("datetime64[ns]") Out[4]: <xarray.DataArray (a: 2, b: 3)> array([['2000-01-01T00:00:00.000000000', '2000-01-02T00:00:00.000000000', '2000-01-03T00:00:00.000000000'], ['2000-01-04T00:00:00.000000000', '2000-01-05T00:00:00.000000000', '2000-01-06T00:00:00.000000000']], dtype='datetime64[ns]') Dimensions without coordinates: a, b ```

NumPy obviously does not officially support this -- nor would I expect it to -- so I would be wary of simply documenting this behavior as is. Would it be reasonable for us to modify xarray.core.duck_array_ops.astype to explicitly implement this conversion ourselves for cftime.datetime arrays? This way we could ensure this was always supported, and we could include appropriate errors for out-of-bounds times (the NumPy method currently overflows in that case) and warnings for switching from non-standard calendars.