Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Thresholded events - Quantified rate2amount #1778

Merged
merged 26 commits into from
Oct 10, 2024
Merged
Show file tree
Hide file tree
Changes from 6 commits
Commits
Show all changes
26 commits
Select commit Hold shift + click to select a range
48efef8
black ice events
coxipi Jun 13, 2024
35b0f75
black ice -> freezing rain
coxipi Jun 13, 2024
dab33b9
Merge branch 'main' of github.com:Ouranosinc/xclim into black_ice
coxipi Jun 13, 2024
44d8681
update CHANGES
coxipi Jun 13, 2024
d036823
event dimension with padding
coxipi Jun 14, 2024
d61182c
merge main
aulemahal Sep 6, 2024
d556d95
first pass to generalize
aulemahal Sep 6, 2024
8ab9442
Merge branch 'main' into black_ice
aulemahal Oct 2, 2024
d27ba5d
cleaner?
aulemahal Oct 2, 2024
c5ff6d3
merge resample-map
aulemahal Oct 3, 2024
3a4ebbb
Working find events
aulemahal Oct 3, 2024
74b8cc9
Add quantified support to rate2amount
aulemahal Oct 3, 2024
9d80208
Fix bad merge...
aulemahal Oct 3, 2024
29a8271
Merge branch 'resample-map' into black_ice
aulemahal Oct 3, 2024
cb6d2cc
Merge branch 'resample-map' into black_ice
aulemahal Oct 7, 2024
f411144
fix fmt, fix tests
aulemahal Oct 7, 2024
cc02b19
Add min_gap to spell length stats
aulemahal Oct 7, 2024
385d59b
Merge branch 'resample-map' into black_ice
aulemahal Oct 7, 2024
1241551
Merge branch 'resample-map' into black_ice
aulemahal Oct 8, 2024
6080aa3
fix test
aulemahal Oct 8, 2024
4457d11
Merge branch 'resample-map' into black_ice
aulemahal Oct 8, 2024
12a4738
Merge branch 'resample-map' into black_ice
aulemahal Oct 9, 2024
d32cf58
Apply suggestions from code review
aulemahal Oct 9, 2024
24e2f7a
Merge branch 'black_ice' of github.com:Ouranosinc/xclim into black_ice
aulemahal Oct 9, 2024
fb6cf7d
Suggestions from review
aulemahal Oct 9, 2024
0f03cdd
Merge branch 'main' into black_ice
aulemahal Oct 10, 2024
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 3 additions & 1 deletion CHANGELOG.rst
Original file line number Diff line number Diff line change
Expand Up @@ -14,13 +14,15 @@ Announcements
New indicators
^^^^^^^^^^^^^^
* New ``heat_spell_frequency``, ``heat_spell_max_length`` and ``heat_spell_total_length`` : spell length statistics on a bivariate condition that uses the average over a window by default. (:pull:`1885`).
* New indicator ``freezing_rain_events`` gives statistics about freezing rain sequences. (:pull:`1778`).

New features and enhancements
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
* New generic ``xclim.indices.generic.spell_mask`` that returns a mask of which days are part of a spell. Supports multivariate conditions and weights. Used in new generic index ``xclim.indices.generic.bivariate_spell_length_statistics`` that extends ``spell_length_statistics`` to two variables. (:pull:`1885`).
* Indicator parameters can now be assigned a new name, different from the argument name in the compute function. (:pull:`1885`).
* New global option ``resample_map_blocks`` to wrap all ``resample().map()`` code inside a ``xr.map_blocks`` to lower the number of dask tasks. Uses utility ``xclim.indices.helpers.resample_map`` and requires ``flox`` to ensure the chunking allows such block-mapping. Defaults to False. (:pull:`1848`).
* ``xclim.indices.run_length.runs_with_holes`` allows to input a condition that must be met for a run to start and a second condition that must be met for the run to stop. (:pull:`1778`).
* New generic compute function ``xclim.indices.generic.thresholded_events`` that finds events based on a threshold condition and returns basic stats for each. See also ``xclim.indices.run_length.find_events``. (:pull:`1778`).
* ``xclim.core.units.rate2amount`` and ``xclim.core.units.amount2rate`` can now also accept quantities (pint objects or strings), in which case the ``dim`` argument must be the ``time`` coordinate through which we can find the sampling rate. (:pull:`1778`).

Bug fixes
^^^^^^^^^
Expand Down
304 changes: 304 additions & 0 deletions tests/test_generic.py
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,7 @@

from xclim.core.calendar import doy_to_days_since, select_time
from xclim.indices import generic
from xclim.testing.helpers import assert_lazy

K2C = 273.15

Expand Down Expand Up @@ -768,3 +769,306 @@ def test_spell_length_statistics_multi(tasmin_series, tasmax_series):
)
xr.testing.assert_equal(outs, outm)
np.testing.assert_allclose(outc, 1)


class TestThresholdedEvents:

@pytest.mark.parametrize("use_dask", [True, False])
def test_simple(self, pr_series, use_dask):
arr = np.array(
[
0,
0,
0,
1,
2,
3,
0,
3,
3,
10,
0,
0,
0,
0,
0,
1,
2,
2,
2,
0,
0,
0,
0,
0,
0,
1,
3,
3,
2,
0,
0,
0,
2,
0,
0,
0,
0,
]
) # noqa
pr = pr_series(arr, start="2000-01-01", units="mm")
if use_dask:
pr = pr.chunk(-1)

with assert_lazy:
out = generic.thresholded_events(
pr,
thresh="1 mm",
op=">=",
window=3,
)

assert out.event.size == np.ceil(arr.size / 6)
out = out.load().dropna("event", how="all")

np.testing.assert_array_equal(out.event_length, [7, 4, 4])
np.testing.assert_array_equal(out.event_effective_length, [6, 4, 4])
np.testing.assert_array_equal(out.event_sum, [22, 7, 9])
np.testing.assert_array_equal(
out.event_start,
np.array(
["2000-01-04", "2000-01-16", "2000-01-26"], dtype="datetime64[ns]"
),
)

@pytest.mark.parametrize("use_dask", [True, False])
def test_diff_windows(self, pr_series, use_dask):
arr = np.array(
[
0,
0,
0,
1,
2,
3,
0,
3,
3,
10,
0,
0,
0,
0,
0,
1,
2,
2,
2,
0,
0,
0,
0,
0,
0,
1,
3,
3,
2,
0,
0,
0,
2,
0,
0,
0,
0,
]
) # noqa
pr = pr_series(arr, start="2000-01-01", units="mm")
if use_dask:
pr = pr.chunk(-1)

# different window stop
out = (
generic.thresholded_events(
pr, thresh="2 mm", op=">=", window=3, window_stop=4
)
.load()
.dropna("event", how="all")
)

np.testing.assert_array_equal(out.event_length, [3, 3, 7])
np.testing.assert_array_equal(out.event_effective_length, [3, 3, 4])
np.testing.assert_array_equal(out.event_sum, [16, 6, 10])
np.testing.assert_array_equal(
out.event_start,
np.array(
["2000-01-08", "2000-01-17", "2000-01-27"], dtype="datetime64[ns]"
),
)

@pytest.mark.parametrize("use_dask", [True, False])
def test_cftime(self, pr_series, use_dask):
arr = np.array(
[
0,
0,
0,
1,
2,
3,
0,
3,
3,
10,
0,
0,
0,
0,
0,
1,
2,
2,
2,
0,
0,
0,
0,
0,
0,
1,
3,
3,
2,
0,
0,
0,
2,
0,
0,
0,
0,
]
) # noqa
pr = pr_series(arr, start="2000-01-01", units="mm").convert_calendar("noleap")
if use_dask:
pr = pr.chunk(-1)

# cftime
with assert_lazy:
out = generic.thresholded_events(
pr,
thresh="1 mm",
op=">=",
window=3,
)
out = out.load().dropna("event", how="all")

np.testing.assert_array_equal(out.event_length, [7, 4, 4])
np.testing.assert_array_equal(out.event_effective_length, [6, 4, 4])
np.testing.assert_array_equal(out.event_sum, [22, 7, 9])
exp = xr.DataArray(
[1, 2, 3],
dims=("time",),
coords={
"time": np.array(
["2000-01-04", "2000-01-16", "2000-01-26"], dtype="datetime64[ns]"
)
},
)
np.testing.assert_array_equal(
out.event_start, exp.convert_calendar("noleap").time
)

@pytest.mark.parametrize("use_dask", [True, False])
def test_freq(self, pr_series, use_dask):
jan = [
0,
0,
0,
1,
2,
3,
0,
3,
3,
10,
0,
0,
0,
0,
0,
0,
2,
2,
2,
2,
2,
2,
0,
0,
0,
0,
0,
3,
2,
3,
2,
] # noqa
fev = [
2,
2,
1,
0,
0,
0,
3,
3,
4,
5,
2,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
0,
] # noqa
pr = pr_series(np.array(jan + fev), start="2000-01-01", units="mm")
if use_dask:
pr = pr.chunk(-1)

with assert_lazy:
out = generic.thresholded_events(
pr, thresh="1 mm", op=">=", window=3, freq="MS"
)
assert out.event_length.shape == (2, 6)
out = out.load().dropna("event", how="all")

np.testing.assert_array_equal(out.event_length, [[7, 6, 4], [3, 5, np.nan]])
np.testing.assert_array_equal(
out.event_effective_length, [[6, 6, 4], [3, 5, np.nan]]
)
np.testing.assert_array_equal(out.event_sum, [[22, 12, 10], [5, 17, np.nan]])
np.testing.assert_array_equal(
out.event_start,
np.array(
[
["2000-01-04", "2000-01-17", "2000-01-28"],
["2000-02-01", "2000-02-07", "NaT"],
],
dtype="datetime64[ns]",
),
)
28 changes: 1 addition & 27 deletions tests/test_indices.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@
from xclim.core import ValidationError
from xclim.core.calendar import percentile_doy
from xclim.core.options import set_options
from xclim.core.units import convert_units_to, rate2amount, units
from xclim.core.units import convert_units_to, units

K2C = 273.15

Expand Down Expand Up @@ -346,32 +346,6 @@ def test_bedd(self, method, end_date, deg_days, max_deg_days):
if method == "icclim":
np.testing.assert_array_equal(bedd, bedd_high_lat)

def test_freezing_rain_events(self, open_dataset):
times = pd.date_range("1950-01-01", "1950-01-31", freq="D")
da = xr.DataArray(
np.zeros(len(times)),
dims={"time"},
coords={"time": times},
attrs={"units": "kg m-2 s-1"},
)

# Two sequences separated by 3 days, which will be split with window_stop = 3
da[0:10] = 1
da[13:20] = 1
out = xci.freezing_rain_events(
da, thresh="0.5 kg m-2 s-1", window_start=3, window_stop=3
)
assert (out.run_lengths.values[0:2] == [10, 7]).all()

# Two sequences separated by 2 days, which will form a single large event
da[12] = 1
out = xci.freezing_rain_events(
da, thresh="0.5 kg m-2 s-1", window_start=3, window_stop=3
)
pram = rate2amount(da)
assert out.run_lengths.values[0] == 20
assert (out.cumulative_precipitation - pram.sum()).sum() == 0

def test_cool_night_index(self, open_dataset):
ds = open_dataset("cmip5/tas_Amon_CanESM2_rcp85_r1i1p1_200701-200712.nc")
ds = ds.rename(dict(tas="tasmin"))
Expand Down
Loading