diff --git a/src/xclim/indices/_threshold.py b/src/xclim/indices/_threshold.py index 9f11be29c..4da2ac250 100644 --- a/src/xclim/indices/_threshold.py +++ b/src/xclim/indices/_threshold.py @@ -9,7 +9,7 @@ import xarray from xclim.core import DayOfYearStr, Quantified -from xclim.core.calendar import doy_from_string, get_calendar +from xclim.core.calendar import doy_from_string, get_calendar, select_time from xclim.core.missing import at_least_n_valid from xclim.core.units import ( convert_units_to, @@ -67,6 +67,8 @@ "growing_season_start", "heat_wave_index", "heating_degree_days", + "holiday_snow_and_snowfall_days", + "holiday_snow_days", "hot_spell_frequency", "hot_spell_max_length", "hot_spell_max_magnitude", @@ -376,7 +378,8 @@ def snd_season_end( window : int Minimum number of days with snow depth below threshold. freq : str - Resampling frequency. The default value is chosen for the northern hemisphere. + Resampling frequency. Default: "YS-JUL". + The default value is chosen for the northern hemisphere. Returns ------- @@ -2766,7 +2769,7 @@ def maximum_consecutive_frost_days( Let :math:`\mathbf{t}=t_0, t_1, \ldots, t_n` be a minimum daily temperature series and :math:`thresh` the threshold below which a day is considered a frost day. Let :math:`\mathbf{s}` be the sorted vector of indices :math:`i` where :math:`[t_i < thresh] \neq [t_{i+1} < thresh]`, that is, the days where the temperature crosses the threshold. - Then the maximum number of consecutive frost days is given by + Then the maximum number of consecutive frost days is given by: .. math:: @@ -2821,7 +2824,7 @@ def maximum_consecutive_dry_days( Let :math:`\mathbf{p}=p_0, p_1, \ldots, p_n` be a daily precipitation series and :math:`thresh` the threshold under which a day is considered dry. Then let :math:`\mathbf{s}` be the sorted vector of indices :math:`i` where :math:`[p_i < thresh] \neq [p_{i+1} < thresh]`, that is, the days where the precipitation crosses the threshold. - Then the maximum number of consecutive dry days is given by + Then the maximum number of consecutive dry days is given by: .. math:: @@ -3633,3 +3636,142 @@ def wet_spell_max_length( resample_before_rl=resample_before_rl, **indexer, ) + + +@declare_units( + snd="[length]", + prsn="[precipitation]", + snd_thresh="[length]", + prsn_thresh="[length]", +) +def holiday_snow_days( + snd: xarray.DataArray, + snd_thresh: Quantified = "20 mm", + date_start: str = "12-25", + date_end: str | None = None, + freq: str = "YS", +) -> xarray.DataArray: # numpydoc ignore=SS05 + r""" + Christmas Days. + + Whether there is a significant amount of snow on the ground on December 25th (or around that time). + + Parameters + ---------- + snd : xarray.DataArray + Surface snow depth. + snd_thresh : Quantified + Threshold snow amount. Default: 20 mm. + date_start : str + Beginning of analysis period. Default: "12-25" (December 25th). + date_end : str, optional + End of analysis period. If not provided, `date_start` is used. + Default: None. + freq : str + Resampling frequency. Default: "YS". + The default value is chosen for the northern hemisphere. + + Returns + ------- + xarray.DataArray, [bool] + Boolean array of years with Christmas Days. + + References + ---------- + https://www.canada.ca/en/environment-climate-change/services/weather-general-tools-resources/historical-christmas-snowfall-data.html + """ + snow_depth = convert_units_to(snd, "mm", context="hydro") + snow_depth_thresh = convert_units_to(snd_thresh, "mm", context="hydro") + snow_depth_constrained = select_time( + snow_depth, + drop=True, + date_bounds=(date_start, date_start if date_end is None else date_end), + ) + + xmas_days = ( + (snow_depth_constrained >= snow_depth_thresh) + .resample(time=freq) + .map(lambda x: x.all(dim="time")) + ) + + xmas_days = xmas_days.assign_attrs({"units": "1"}) + return xmas_days + + +@declare_units( + snd="[length]", + prsn="[precipitation]", + snd_thresh="[length]", + prsn_thresh="[length]", +) +def holiday_snow_and_snowfall_days( + snd: xarray.DataArray, + prsn: xarray.DataArray | None = None, + snd_thresh: Quantified = "20 mm", + prsn_thresh: Quantified = "10 mm", + date_start: str = "12-25", + date_end: str | None = None, + freq: str = "YS-JUL", +) -> xarray.DataArray: + r""" + Perfect Christmas Days. + + Whether there is a significant amount of snow on the ground and snowfall occurring on December 25th. + + Parameters + ---------- + snd : xarray.DataArray + Surface snow depth. + prsn : xarray.DataArray + Snowfall flux. + snd_thresh : Quantified + Threshold snow amount. Default: 20 mm. + prsn_thresh : Quantified + Threshold snowfall flux. Default: 10 mm. + date_start : str + Beginning of analysis period. Default: "12-25" (December 25th). + date_end : str, optional + End of analysis period. If not provided, `date_start` is used. + Default: None. + freq : str + Resampling frequency. Default: "YS-JUL". + The default value is chosen for the northern hemisphere. + + Returns + ------- + xarray.DataArray, [bool] + Boolean array of years with Perfect Christmas Days. + + References + ---------- + https://www.canada.ca/en/environment-climate-change/services/weather-general-tools-resources/historical-christmas-snowfall-data.html + """ + snowfall_rate = rate2amount( + convert_units_to(prsn, "mm/d", context="hydro"), out_units="mm" + ) + snowfall_thresh = convert_units_to(prsn_thresh, "mm", context="hydro") + snowfall_constrained = select_time( + snowfall_rate, + drop=True, + date_bounds=(date_start, date_start if date_end is None else date_end), + ) + + snow_depth = convert_units_to(snd, "mm", context="hydro") + snow_depth_thresh = convert_units_to(snd_thresh, "mm", context="hydro") + snow_depth_constrained = select_time( + snow_depth, + drop=True, + date_bounds=(date_start, date_start if date_end is None else date_end), + ) + + perfect_xmas_days = ( + ( + (snow_depth_constrained >= snow_depth_thresh) + & (snowfall_constrained >= snowfall_thresh) + ) + .resample(time=freq) + .map(lambda x: x.all(dim="time")) + ) + + perfect_xmas_days = perfect_xmas_days.assign_attrs({"units": "1"}) + return perfect_xmas_days