issue_comments

I assume you figured this out somehow. If not, feel free to reopen.

is this still an issue, @mada0007? If it is, I tried a bit more (it is still incomplete since I don't really understand what happens, especially with the last groupby-agg) python In [7]: def get_grps(s, axis, thresh=-1, Nmin=2): ...: s = xr.DataArray(s, dims=("time", "lon", "lat", "window")) ...: m = ( ...: xr.concat( ...: [ ...: s.shift(window=-i).pipe(lambda arr: arr <= thresh) ...: for i in range(Nmin) ...: ], ...: dim="shifted", ...: ) ...: .reduce(lambda arr, axis: np.logical_and.reduce(arr, axis=axis), dim="shifted") ...: ) ...: if Nmin > 1: ...: m = ( ...: m.where(lambda arr: arr == False, np.nan) ...: .ffill(limit=Nmin - 1, dim="window") ...: .fillna(False) ...: ) ...: ...: gps = (m != m.shift(window=1)).cumsum(dim="window").where(m) ...: ...: groups = s.groupby(gps.rename("gps_groups")) ...: return xr.where( ...: gps.isnull().all(dim="window"), ...: 0, ...: 1, # use `groups` here once it has the same size as `gps` ...: ) ...: ...: spi = xr.DataArray(dims=("time", "lon", "lat"), data=np.random.randn(118, 18, 36)) ...: spi.rolling(time=59).reduce(get_grps) which mirrors your function except it does the calculation for all points in your space-time-grid at once.

Of course, once again I can't really check for errors since I don't understand what I'm computing :grin:

the reduce function needs to be a function that takes a 4D array (time x lat x lon x window in your case). It should use the elements of the window dimension (1 x 1 x 1 x 59) to calculate a time x lat x lon array.

The only thing wrong with your function is that it ignores the axis parameter: it can only work on 1D arrays (with a shape of 59 in this case).

A simple albeit horribly slow way to use your function without rewriting it would be numpy.apply_along_axis (slow because the function is called time * lat * lon == 20995200 times): python spi.rolling(time=59).reduce(lambda x, axis: np.apply_along_axis(get_grps, axis, x))

I'm not too familiar with writing numpy reduce functions so you might want to wait for the opinion of someone more experienced.

I started to rewrite your function using xarray.DataArray instead of pandas.Series, but I don't have time to finish that right now. Hopefully this points you in the right direction. python In [32]: axis = 3 ...: thresh = -1 ...: Nmin = 2 ...: # function code starts here ...: dim_names = ["time", "lat", "lon"] ...: # fails if axis == -1 but axis from reduce is always a positive number ...: dim_names.insert(axis, "window") ...: # use s = xr.DataArray(s, dims=dim_names) in the function ...: # so we can interactively write and try the function code ...: s = spi.rolling(time=59).construct("window") ...: shifted = xr.concat([s.shift({"window": i}) for i in range(Nmin)], dim="shifted") ...: thresholded = shifted <= thresh ...: # xarray does not support `ufunc.reduce`, so we need to fall back to numpy arrays ...: m = r.copy(data=np.logical_and.reduce(thresholded.values, axis=0)) ...: if Nmin > 1: ...: m = m.where(~m).ffill(limit=Nmin - 1, dim="window").fillna(False) ...: gps = m != m.shift({"window": 1}).cumsum().where(m)

the reason for the exception is that you are trying to convert the input data (a 4D numpy.ndarray of shape 324 x 180 x 360 x 59) to a pandas.Series which only takes 1D data. You will have to rewrite your function to accept a 4D array and reduce it to a 3D array along the given axis (see numpy functions like std or mean).

We should definitely improve the documentation, right now the docstring of DataArrayRolling.reduce and DatasetRolling.reduce is not really clear on the requirements for the reduce function and the examples don't show the use of custom functions.