diff --git a/pyrealm/pmodel/scaler.py b/pyrealm/pmodel/scaler.py index 4fa8d0a7..c7e9e7e2 100644 --- a/pyrealm/pmodel/scaler.py +++ b/pyrealm/pmodel/scaler.py @@ -407,6 +407,7 @@ def get_daily_means( def fill_daily_to_subdaily( self, values: NDArray, + previous_value: NDArray | None = None, update_point: str = "max", kind: str = "previous", fill_from: np.timedelta64 | None = None, @@ -419,45 +420,53 @@ def fill_daily_to_subdaily( axis of the `values` must be the same length as the number of days used to create the instance. + The update point defaults to the maximum time of day during the acclimation + window. It can also be set to the mean time of day, but note that this implies + that the plant predicts the daily values between the mean and max observation + time. The ``fill_from`` argument can be used to set the update point to an + arbitrary time of day. + Two interpolation kinds are currently implemented: * ``previous`` interpolates the daily value as a constant, until updating to the next daily value. This option will fill values until the end of the time - series. + series. The * ``linear`` interpolates linearly between the update points of the daily values. The interpolated values are held constant for the first day and then interpolated linearly: this is to avoid plants adapting optimally to future conditions. - The update point defaults to the maximum time of day during the acclimation - window. It can also be set to the mean time of day, but note that this implies - that the plant predicts the daily values between the mean and max observation - time. The ``fill_from`` argument can be used to set the update point to an - arbitrary time of day. + Subdaily observations before the update point on the first day of the time + series are filled with ``np.nan``. The ``previous_value`` argument can be used + to provide an alternative value, allowing time series to be processed in blocks, + but this option is only currently implemented for the ``previous`` + interpolation method. Args: values: An array with the first dimension matching the number of days in the - instances :class:`~pyrealm.pmodel.scaler.SubdailyScaler` object. + instances :class:`~pyrealm.pmodel.scaler.SubdailyScaler` object. + previous_value: An array with dimensions equal to a slice across the first + axis of the values array. update_point: The point in the acclimation window at which the plant updates - to the new daily value: one of 'mean' or 'max'. + to the new daily value: one of 'mean' or 'max'. kind: The kind of interpolation to use to fill between daily values: one of - 'previous' or 'linear', + 'previous' or 'linear', fill_from: As an alternative to ``update_point``, an - :class:`numpy.timedelta64` value giving the time of day from which to fill - values forward. + :class:`numpy.timedelta64` value giving the time of day from which to + fill values forward. """ if values.shape[0] != self.n_days: - raise ValueError("Values is not of length n_days on its first axis.") + raise ValueError("Values is not of length n_days on its first axis") if fill_from is not None: if not isinstance(fill_from, np.timedelta64): - raise ValueError("The fill_from argument must be a timedelta64 value.") + raise ValueError("The fill_from argument must be a timedelta64 value") # Convert to seconds and check it is in range _fill_from = fill_from.astype("timedelta64[s]") if not (_fill_from >= 0 and _fill_from < 24 * 60 * 60): - raise ValueError("The fill_from argument is not >= 0 and < 24 hours.") + raise ValueError("The fill_from argument is not >= 0 and < 24 hours") update_time = self.observation_dates + _fill_from @@ -468,8 +477,23 @@ def fill_daily_to_subdaily( else: raise ValueError("Unknown update point") - # Note that interp1d cannot handle datetime64 inputs, so need to interpolate - # using datetimes cast to integer types + # Check previous value settings - only allow with previous interpolation and + # check the previous value shape matches + if previous_value is not None: + if kind == "linear": + raise NotImplementedError( + "Using previous value with kind='linear' is not implemented" + ) + + # Use np.broadcast_shapes here to handle checking array shapes. This is + # mostly to catch the fact that () and (1,) are equivalent. + try: + np.broadcast_shapes(previous_value.shape, values.shape) + except ValueError: + raise ValueError( + "The input to previous_value is not congruent with " + "the shape of the observed data" + ) # Use fill_value to handle extrapolation before or after update point: # - previous will fill the last value forward to the end of the time series, @@ -479,8 +503,16 @@ def fill_daily_to_subdaily( # value until _after_ the update point. if kind == "previous": - fill_value = (None, values[-1]) + # The fill values here are used to extend the last daily value out to the + # end of the subdaily observations but also to fill any provided previous + # values for subdaily observations _before_ the first daily value. If + # the default previous value of None is supplied, this inserts np.nan as + # expected. + fill_value = (previous_value, values[-1]) elif kind == "linear": + # Shift the values forward by a day, inserting a copy of the first day at + # the start. This then avoids plants seeing the future and provides values + # up until the last observation. values = np.insert(values, 0, values[0], axis=0) update_time = np.append( update_time, update_time[-1] + np.timedelta64(1, "D") @@ -489,6 +521,8 @@ def fill_daily_to_subdaily( else: raise ValueError("Unsupported interpolation option") + # Note that interp1d cannot handle datetime64 inputs, so need to interpolate + # using datetimes cast to integer types interp_fun = interp1d( update_time.astype("int"), values, @@ -496,6 +530,7 @@ def fill_daily_to_subdaily( kind=kind, bounds_error=False, fill_value=fill_value, + assume_sorted=True, ) # TODO - The kind "previous" might be replaceable with bottleneck.push diff --git a/pyrealm/pmodel/subdaily.py b/pyrealm/pmodel/subdaily.py index 0367ce38..24d871d9 100644 --- a/pyrealm/pmodel/subdaily.py +++ b/pyrealm/pmodel/subdaily.py @@ -41,7 +41,10 @@ def memory_effect( - values: NDArray, alpha: float = 0.067, allow_holdover: bool = False + values: NDArray, + previous_values: NDArray | None = None, + alpha: float = 0.067, + allow_holdover: bool = False, ) -> NDArray: r"""Apply a memory effect to a variable. @@ -89,6 +92,8 @@ def memory_effect( Args: values: The values to apply the memory effect to. + previous_values: Last available realised value used if model is fitted in + chunks and value at t=0 is not optimal. alpha: The relative weight applied to the most recent observation. allow_holdover: Allow missing values to be filled by holding over earlier values. @@ -105,7 +110,10 @@ def memory_effect( # Initialise the output storage and set the first values to be a slice along the # first axis of the input values memory_values = np.empty_like(values, dtype=np.float32) - memory_values[0] = values[0] + if previous_values is None: + memory_values[0] = values[0] + else: + memory_values[0] = previous_values * (1 - alpha) + values[0] * alpha # Handle the data if there are no missing data, if not nan_present: @@ -177,6 +185,10 @@ class SubdailyPModel: more rapid acclimation: :math:`\alpha=1` results in immediate acclimation and :math:`\alpha=0` results in no acclimation at all, with values pinned to the initial estimates. + * By default, the initial realised value :math:`R_1` for each of the three slowly + acclimating variables is assumed to be the first optimal value :math:`O_1`, but + the `previous_realised` argument can be used to provide values of :math:`R_0` from + which to calculate :math:`R_{1} = R_{0}(1 - \alpha) + O_{1} \alpha`. * The realised values are then filled back onto the original subdaily timescale, with :math:`V_{cmax}` and :math:`J_{max}` then being calculated from the slowly responding :math:`V_{cmax25}` and :math:`J_{max25}` and the actual subdaily @@ -221,10 +233,12 @@ class SubdailyPModel: allow_partial_data: Should estimates of daily optimal conditions be calculated with missing values in the acclimation window. reference_kphio: An optional alternative reference value for the quantum yield - efficiency of photosynthesis (:math:`\phi_0`, -) to be passed to the kphio - calculation method. + efficiency of photosynthesis (:math:`\phi_0`, -) to be passed to the kphio + calculation method. fill_kind: The approach used to fill daily realised values to the subdaily timescale, currently one of 'previous' or 'linear'. + previous_realised: A tuple of previous realised values of three NumPy arrays + (xi_real, vcmax25_real, jmax25_real). """ def __init__( @@ -241,6 +255,7 @@ def __init__( allow_holdover: bool = False, allow_partial_data: bool = False, fill_kind: str = "previous", + previous_realised: tuple[NDArray, NDArray, NDArray] | None = None, ) -> None: # Warn about the API warn( @@ -379,25 +394,73 @@ def __init__( 1 / calc_ftemp_arrh(tk_acclim, self.env.pmodel_const.subdaily_jmax25_ha) ) + """Instantaneous optimal :math:`x_{i}`, :math:`V_{cmax}` and :math:`J_{max}`""" + # Check the shape of previous realised values are congruent with a slice across + # the time axis + if previous_realised is not None: + if fill_kind != "previous": + raise NotImplementedError( + "Using previous_realised is only implemented for " + "fill_kind = 'previous'" + ) + + # All variables should share the shape of a slice along the first axis of + # the environmental forcings + expected_shape = self.env.tc[0].shape + if not ( + (previous_realised[0].shape == expected_shape) + and (previous_realised[1].shape == expected_shape) + and (previous_realised[2].shape == expected_shape) + ): + raise ValueError( + "`previous_realised` entries have wrong shape in Subdaily PModel" + ) + else: + previous_xi_real, previous_vcmax25_real, previous_jmax25_real = ( + previous_realised + ) + else: + previous_xi_real, previous_vcmax25_real, previous_jmax25_real = [ + None, + None, + None, + ] + # 5) Calculate the realised daily values from the instantaneous optimal values self.xi_real: NDArray = memory_effect( - self.pmodel_acclim.optchi.xi, alpha=alpha, allow_holdover=allow_holdover + self.pmodel_acclim.optchi.xi, + previous_values=previous_xi_real, + alpha=alpha, + allow_holdover=allow_holdover, ) r"""Realised daily slow responses in :math:`\xi`""" self.vcmax25_real: NDArray = memory_effect( - self.vcmax25_opt, alpha=alpha, allow_holdover=allow_holdover + self.vcmax25_opt, + previous_values=previous_vcmax25_real, + alpha=alpha, + allow_holdover=allow_holdover, ) r"""Realised daily slow responses in :math:`V_{cmax25}`""" self.jmax25_real: NDArray = memory_effect( - self.jmax25_opt, alpha=alpha, allow_holdover=allow_holdover + self.jmax25_opt, + previous_values=previous_jmax25_real, + alpha=alpha, + allow_holdover=allow_holdover, ) + r"""Realised daily slow responses in :math:`J_{max25}`""" # 6) Fill the realised xi, jmax25 and vcmax25 from daily values back to the # subdaily timescale. - self.subdaily_vcmax25 = fs_scaler.fill_daily_to_subdaily(self.vcmax25_real) - self.subdaily_jmax25 = fs_scaler.fill_daily_to_subdaily(self.jmax25_real) - self.subdaily_xi = fs_scaler.fill_daily_to_subdaily(self.xi_real) + self.subdaily_xi = fs_scaler.fill_daily_to_subdaily( + self.xi_real, previous_value=previous_xi_real + ) + self.subdaily_vcmax25 = fs_scaler.fill_daily_to_subdaily( + self.vcmax25_real, previous_value=previous_vcmax25_real + ) + self.subdaily_jmax25 = fs_scaler.fill_daily_to_subdaily( + self.jmax25_real, previous_value=previous_jmax25_real + ) # 7) Adjust subdaily jmax25 and vcmax25 back to jmax and vcmax given the # actual subdaily temperatures. diff --git a/tests/unit/pmodel/test_memory_effect.py b/tests/unit/pmodel/test_memory_effect.py index b4d49cc2..d5ab4402 100644 --- a/tests/unit/pmodel/test_memory_effect.py +++ b/tests/unit/pmodel/test_memory_effect.py @@ -54,6 +54,43 @@ def test_memory_effect(inputs, alpha): assert np.allclose(result, expected) +@pytest.mark.parametrize( + argnames="inputs_whole", + argvalues=[ + pytest.param(np.arange(0, 10), id="1D"), + pytest.param( + np.column_stack([np.arange(0, 10)] * 4) + np.arange(4), + id="2D", + ), + pytest.param( + np.dstack([np.column_stack([np.arange(0, 10)] * 4)] * 4) + + np.arange(16).reshape(4, 4), + id="3D", + ), + ], +) +@pytest.mark.parametrize(argnames="alpha", argvalues=(0.0, 0.5, 1.0)) +def test_memory_effect_chunked(inputs_whole, alpha): + """Test that the memory effect works when chunking the time series up. + + This compares the output of `test_memory_effect` with the output of + `memory_effect` which gets the two time chunks fed sequentially., + """ + from pyrealm.pmodel import memory_effect + + result_whole = memory_effect(inputs_whole, alpha=alpha) + + [inputs_chunk1, inputs_chunk2] = np.split(inputs_whole, [5], axis=0) + + result_chunk1 = memory_effect(inputs_chunk1, alpha=alpha) + + result_chunk2 = memory_effect( + inputs_chunk2, previous_values=result_chunk1[-1], alpha=alpha + ) + + assert np.allclose(result_whole[-1], result_chunk2[-1]) + + @pytest.mark.parametrize( argnames="inputs,allow_holdover,context_manager,expected", argvalues=[ diff --git a/tests/unit/pmodel/test_subdaily.py b/tests/unit/pmodel/test_subdaily.py index 0796086d..471f60dd 100644 --- a/tests/unit/pmodel/test_subdaily.py +++ b/tests/unit/pmodel/test_subdaily.py @@ -114,6 +114,90 @@ def test_FSPModel_corr(be_vie_data_components, data_args): assert np.all(r_vals > 0.995) +def test_SubdailyPModel_previous_realised(be_vie_data_components): + """Test the functionality that allows the subdaily model to restart in blocks.""" + + from pyrealm.pmodel import SubdailyScaler + from pyrealm.pmodel.subdaily import SubdailyPModel + + # Run all in one model + env, ppfd, fapar, datetime, _ = be_vie_data_components.get( + mode="crop", start=0, end=17520 + ) + + # Get the fast slow scaler and set window + fsscaler = SubdailyScaler(datetime) + fsscaler.set_window( + window_center=np.timedelta64(12, "h"), + half_width=np.timedelta64(30, "m"), + ) + + # Run as a subdaily model using the kphio used in the reference implementation. + all_in_one_subdaily_pmodel = SubdailyPModel( + env=env, + ppfd=ppfd, + fapar=fapar, + reference_kphio=1 / 8, + fs_scaler=fsscaler, + allow_holdover=True, + ) + + # Run first half of year + env1, ppfd1, fapar1, datetime1, _ = be_vie_data_components.get( + mode="crop", start=0, end=182 * 48 + ) + + # Get the fast slow scaler and set window + fsscaler1 = SubdailyScaler(datetime1) + fsscaler1.set_window( + window_center=np.timedelta64(12, "h"), + half_width=np.timedelta64(30, "m"), + ) + + # Run as a subdaily model using the kphio used in the reference implementation. + part_1_subdaily_pmodel = SubdailyPModel( + env=env1, + ppfd=ppfd1, + fapar=fapar1, + reference_kphio=1 / 8, + fs_scaler=fsscaler1, + allow_holdover=True, + ) + + # Run second year + env2, ppfd2, fapar2, datetime2, _ = be_vie_data_components.get( + mode="crop", start=182 * 48, end=365 * 48 + ) + + # Get the fast slow scaler and set window + fsscaler2 = SubdailyScaler(datetime2) + fsscaler2.set_window( + window_center=np.timedelta64(12, "h"), + half_width=np.timedelta64(30, "m"), + ) + + # Run as a subdaily model using the kphio used in the reference implementation. + part_2_subdaily_pmodel = SubdailyPModel( + env=env2, + ppfd=ppfd2, + fapar=fapar2, + reference_kphio=1 / 8, + fs_scaler=fsscaler2, + allow_holdover=True, + previous_realised=( + part_1_subdaily_pmodel.optimal_chi.xi[-1], + part_1_subdaily_pmodel.vcmax25_real[-1], + part_1_subdaily_pmodel.jmax25_real[-1], + ), + ) + + assert np.allclose( + all_in_one_subdaily_pmodel.gpp, + np.concat([part_1_subdaily_pmodel.gpp, part_2_subdaily_pmodel.gpp]), + equal_nan=True, + ) + + @pytest.mark.parametrize("ndims", [2, 3, 4]) def test_FSPModel_dimensionality(be_vie_data, ndims): """Tests that the SubdailyPModel handles dimensions correctly. diff --git a/tests/unit/pmodel/test_fast_slow_scaler.py b/tests/unit/pmodel/test_subdailyscaler.py similarity index 79% rename from tests/unit/pmodel/test_fast_slow_scaler.py rename to tests/unit/pmodel/test_subdailyscaler.py index 1dfc1758..977095dc 100644 --- a/tests/unit/pmodel/test_fast_slow_scaler.py +++ b/tests/unit/pmodel/test_subdailyscaler.py @@ -11,7 +11,7 @@ @pytest.fixture -def fixture_FSS(): +def fixture_SubdailyScaler(): """A fixture providing a SubdailyScaler object.""" from pyrealm.pmodel import SubdailyScaler @@ -114,7 +114,7 @@ def fixture_FSS(): ), ], ) -def test_FSS_init(ctext_mngr, msg, datetimes): +def test_SubdailyScaler_init(ctext_mngr, msg, datetimes): """Test the SubdailyScaler init handling of date ranges.""" from pyrealm.pmodel import SubdailyScaler @@ -162,11 +162,13 @@ def test_FSS_init(ctext_mngr, msg, datetimes): ), ], ) -def test_FSS_set_window(fixture_FSS, ctext_mngr, msg, kwargs, samp_mean, samp_max): +def test_SubdailyScaler_set_window( + fixture_SubdailyScaler, ctext_mngr, msg, kwargs, samp_mean, samp_max +): """Test the SubdailyScaler set_window method.""" with ctext_mngr as cman: - fixture_FSS.set_window(**kwargs) + fixture_SubdailyScaler.set_window(**kwargs) if msg is not None: assert str(cman.value) == msg @@ -174,8 +176,8 @@ def test_FSS_set_window(fixture_FSS, ctext_mngr, msg, kwargs, samp_mean, samp_ma # Check that _set_times has run correctly. Can't use allclose directly on # datetimes and since these are integers under the hood, don't need float # testing - assert np.all(fixture_FSS.sample_datetimes_mean == samp_mean) - assert np.all(fixture_FSS.sample_datetimes_max == samp_max) + assert np.all(fixture_SubdailyScaler.sample_datetimes_mean == samp_mean) + assert np.all(fixture_SubdailyScaler.sample_datetimes_max == samp_max) @pytest.mark.parametrize( @@ -227,10 +229,12 @@ def test_FSS_set_window(fixture_FSS, ctext_mngr, msg, kwargs, samp_mean, samp_ma ), ], ) -def test_FSS_set_include(fixture_FSS, ctext_mngr, msg, include, samp_mean, samp_max): +def test_SubdailyScaler_set_include( + fixture_SubdailyScaler, ctext_mngr, msg, include, samp_mean, samp_max +): """Test the SubdailyScaler set_include method.""" with ctext_mngr as cman: - fixture_FSS.set_include(include) + fixture_SubdailyScaler.set_include(include) if msg is not None: assert str(cman.value) == msg @@ -239,8 +243,8 @@ def test_FSS_set_include(fixture_FSS, ctext_mngr, msg, include, samp_mean, samp_ # Check that _set_times has run correctly. Can't use allclose directly on # datetimes and since these are integers under the hood, don't need float # testing - assert np.all(fixture_FSS.sample_datetimes_mean == samp_mean) - assert np.all(fixture_FSS.sample_datetimes_max == samp_max) + assert np.all(fixture_SubdailyScaler.sample_datetimes_mean == samp_mean) + assert np.all(fixture_SubdailyScaler.sample_datetimes_max == samp_max) @pytest.mark.parametrize( @@ -290,10 +294,12 @@ def test_FSS_set_include(fixture_FSS, ctext_mngr, msg, include, samp_mean, samp_ ), ], ) -def test_FSS_set_nearest(fixture_FSS, ctext_mngr, msg, time, samp_mean, samp_max): +def test_SubdailyScaler_set_nearest( + fixture_SubdailyScaler, ctext_mngr, msg, time, samp_mean, samp_max +): """Test the SubdailyScaler set_nearest method.""" with ctext_mngr as cman: - fixture_FSS.set_nearest(time) + fixture_SubdailyScaler.set_nearest(time) if msg is not None: assert str(cman.value) == msg @@ -302,8 +308,8 @@ def test_FSS_set_nearest(fixture_FSS, ctext_mngr, msg, time, samp_mean, samp_max # Check that _set_times has run correctly. Can't use allclose directly on # datetimes and since these are integers under the hood, don't need float # testing - assert np.all(fixture_FSS.sample_datetimes_mean == samp_mean) - assert np.all(fixture_FSS.sample_datetimes_max == samp_max) + assert np.all(fixture_SubdailyScaler.sample_datetimes_mean == samp_mean) + assert np.all(fixture_SubdailyScaler.sample_datetimes_max == samp_max) @pytest.mark.parametrize( @@ -317,15 +323,15 @@ def test_FSS_set_nearest(fixture_FSS, ctext_mngr, msg, time, samp_mean, samp_max ), ], ) -def test_FSS_get_wv_errors(fixture_FSS, ctext_mngr, msg, values): +def test_SubdailyScaler_get_wv_errors(fixture_SubdailyScaler, ctext_mngr, msg, values): """Test errors arising in the SubdailyScaler get_window_value method.""" - fixture_FSS.set_window( + fixture_SubdailyScaler.set_window( window_center=np.timedelta64(12, "h"), half_width=np.timedelta64(2, "h"), ) with ctext_mngr as cman: - _ = fixture_FSS.get_window_values(values) + _ = fixture_SubdailyScaler.get_window_values(values) assert str(cman.value) == msg @@ -488,8 +494,8 @@ def test_FSS_get_wv_errors(fixture_FSS, ctext_mngr, msg, values): ), ], ) -class Test_FSS_get_vals_window_and_include: - """Test FSS get methods for set_window and set_include. +class Test_SubdailyScaler_get_vals_window_and_include: + """Test SubdailyScaler get methods for set_window and set_include. The daily values extracted using the set_window and set_include methods can be the same, by setting the window and the include to cover the same observations, so these @@ -509,30 +515,30 @@ class Test_FSS_get_vals_window_and_include: are np.nan - this should revert to setting np.nan in the first day. """ - def test_FSS_get_vals_window( - self, fixture_FSS, values, expected_means, allow_partial_data + def test_SubdailyScaler_get_vals_window( + self, fixture_SubdailyScaler, values, expected_means, allow_partial_data ): """Test a window.""" - fixture_FSS.set_window( + fixture_SubdailyScaler.set_window( window_center=np.timedelta64(12, "h"), half_width=np.timedelta64(2, "h"), ) - calculated_means = fixture_FSS.get_daily_means( + calculated_means = fixture_SubdailyScaler.get_daily_means( values, allow_partial_data=allow_partial_data ) assert np.allclose(calculated_means, expected_means, equal_nan=True) - def test_FSS_get_vals_include( - self, fixture_FSS, values, expected_means, allow_partial_data + def test_SubdailyScaler_get_vals_include( + self, fixture_SubdailyScaler, values, expected_means, allow_partial_data ): """Test include.""" # This duplicates the selection of the window test but using direct include inc = np.zeros(48, dtype=np.bool_) inc[20:29] = True - fixture_FSS.set_include(inc) - calculated_means = fixture_FSS.get_daily_means( + fixture_SubdailyScaler.set_include(inc) + calculated_means = fixture_SubdailyScaler.get_daily_means( values, allow_partial_data=allow_partial_data ) @@ -607,7 +613,9 @@ def test_FSS_get_vals_include( ), ], ) -def test_FSS_get_vals_nearest(fixture_FSS, values, expected_means): +def test_SubdailyScaler_get_vals_nearest( + fixture_SubdailyScaler, values, expected_means +): """Test get_daily_values. This tests the specific behaviour when set_nearest is used and a single observation @@ -616,8 +624,8 @@ def test_FSS_get_vals_nearest(fixture_FSS, values, expected_means): """ # Select the 11:30 observation, which is missing in PARTIAL_ONES and PARTIAL_VARYING - fixture_FSS.set_nearest(np.timedelta64(11 * 60 + 29, "m")) - calculated_means = fixture_FSS.get_daily_means(values) + fixture_SubdailyScaler.set_nearest(np.timedelta64(11 * 60 + 29, "m")) + calculated_means = fixture_SubdailyScaler.get_daily_means(values) assert np.allclose(calculated_means, expected_means, equal_nan=True) @@ -678,27 +686,37 @@ def test_FSS_get_vals_nearest(fixture_FSS, values, expected_means): ], ) @pytest.mark.parametrize( - argnames=["ctext_mngr", "msg", "input_values", "exp_values", "fill_from"], + argnames=["input_values", "exp_values", "fill_from", "previous_value"], argvalues=[ pytest.param( - does_not_raise(), - None, np.array([1, 2, 3]), np.repeat([np.nan, 1, 2, 3], (26, 48, 48, 22)), None, + None, id="1D test", ), pytest.param( - does_not_raise(), - None, np.array([1, 2, 3]), np.repeat([1, 2, 3], (48, 48, 48)), np.timedelta64(0, "h"), + None, id="1D test - fill from", ), pytest.param( - does_not_raise(), + np.array([1, 2, 3]), + np.repeat([0, 1, 2, 3], (26, 48, 48, 22)), + None, + np.array([0]), + id="1D test - previous value 1D", + ), + pytest.param( + np.array([1, 2, 3]), + np.repeat([0, 1, 2, 3], (26, 48, 48, 22)), None, + np.array(0), + id="1D test - previous value 0D", + ), + pytest.param( np.array([[[1, 4], [7, 10]], [[2, 5], [8, 11]], [[3, 6], [9, 12]]]), np.repeat( a=[ @@ -711,11 +729,10 @@ def test_FSS_get_vals_nearest(fixture_FSS, values, expected_means): axis=0, ), None, - id="2D test", + None, + id="3D test", ), pytest.param( - does_not_raise(), - None, np.array([[[1, 4], [7, 10]], [[2, 5], [8, 11]], [[3, 6], [9, 12]]]), np.repeat( a=[ @@ -728,11 +745,26 @@ def test_FSS_get_vals_nearest(fixture_FSS, values, expected_means): axis=0, ), np.timedelta64(2, "h"), - id="2D test - fill from", + None, + id="3D test - fill from", ), pytest.param( - does_not_raise(), + np.array([[[1, 4], [7, 10]], [[2, 5], [8, 11]], [[3, 6], [9, 12]]]), + np.repeat( + a=[ + [[0, 3], [6, 9]], + [[1, 4], [7, 10]], + [[2, 5], [8, 11]], + [[3, 6], [9, 12]], + ], + repeats=[26, 48, 48, 22], + axis=0, + ), None, + np.array([[0, 3], [6, 9]]), + id="3D test - previous value 2D", + ), + pytest.param( np.array([[1, 4], [2, 5], [3, 6]]), np.repeat( a=[[np.nan, np.nan], [1, 4], [2, 5], [3, 6]], @@ -740,37 +772,41 @@ def test_FSS_get_vals_nearest(fixture_FSS, values, expected_means): axis=0, ), None, - id="3D test", + None, + id="2D test", ), ], ) -def test_FSS_resample_subdaily( - fixture_FSS, +def test_SubdailyScaler_fill_daily_to_subdaily_previous( + fixture_SubdailyScaler, method_name, kwargs, update_point, - ctext_mngr, - msg, input_values, exp_values, fill_from, + previous_value, ): - """Test the calculation of subdaily samples using SubdailyScaler.""" + """Test fill_daily_to_subdaily using SubdailyScale with method previous. + + The first parameterisation sets the exact same acclimation windows in a bunch of + different ways. The second paramaterisation provides inputs with different + dimensionality. + """ # Set the included observations - the different parameterisations here and for # the update point should all select the same update point. - func = getattr(fixture_FSS, method_name) + func = getattr(fixture_SubdailyScaler, method_name) func(**kwargs) - with ctext_mngr as cman: - res = fixture_FSS.fill_daily_to_subdaily( - input_values, update_point=update_point, fill_from=fill_from - ) + res = fixture_SubdailyScaler.fill_daily_to_subdaily( + input_values, + update_point=update_point, + fill_from=fill_from, + previous_value=previous_value, + ) - if cman is not None: - assert str(cman.value) == msg - else: - assert np.allclose(res, exp_values, equal_nan=True) + assert np.allclose(res, exp_values, equal_nan=True) @pytest.mark.parametrize( @@ -829,21 +865,91 @@ def test_FSS_resample_subdaily( ), ], ) -def test_FSS_resample_subdaily_linear( - fixture_FSS, +def test_SubdailyScaler_fill_daily_to_subdaily_linear( + fixture_SubdailyScaler, update_point, input_values, exp_values, ): - """Test SubdailyScaler resampling to subdaily timescale by linear interpolation.""" + """Test fill_daily_to_subdaily using SubdailyScaler with method linear.""" # Set the included observations - fixture_FSS.set_window( + fixture_SubdailyScaler.set_window( window_center=np.timedelta64(13, "h"), half_width=np.timedelta64(1, "h") ) - res = fixture_FSS.fill_daily_to_subdaily( + res = fixture_SubdailyScaler.fill_daily_to_subdaily( input_values, update_point=update_point, kind="linear" ) assert np.allclose(res, exp_values, equal_nan=True) + + +@pytest.mark.parametrize( + argnames="inputs, outcome, msg", + argvalues=[ + pytest.param( + {"values": np.arange(12)}, + pytest.raises(ValueError), + "Values is not of length n_days on its first axis", + id="values wrong shape", + ), + pytest.param( + {"values": np.arange(3), "fill_from": 3}, + pytest.raises(ValueError), + "The fill_from argument must be a timedelta64 value", + id="fill_from not timedelta64", + ), + pytest.param( + {"values": np.arange(3), "fill_from": np.timedelta64(12, "D")}, + pytest.raises(ValueError), + "The fill_from argument is not >= 0 and < 24 hours", + id="fill_from too large", + ), + pytest.param( + {"values": np.arange(3), "fill_from": np.timedelta64(-1, "s")}, + pytest.raises(ValueError), + "The fill_from argument is not >= 0 and < 24 hours", + id="fill_from negative", + ), + pytest.param( + {"values": np.arange(3), "update_point": "noon"}, + pytest.raises(ValueError), + "Unknown update point", + id="unknown update point", + ), + pytest.param( + {"values": np.arange(3), "previous_value": np.array(1), "kind": "linear"}, + pytest.raises(NotImplementedError), + "Using previous value with kind='linear' is not implemented", + id="previous_value with linear", + ), + pytest.param( + {"values": np.arange(3), "previous_value": np.ones(4)}, + pytest.raises(ValueError), + "The input to previous_value is not congruent with " + "the shape of the observed data", + id="previous_value shape issue", + ), + pytest.param( + {"values": np.arange(3), "kind": "quadratic"}, + pytest.raises(ValueError), + "Unsupported interpolation option", + id="unsupported interpolation", + ), + ], +) +def test_SubdailyScaler_fill_daily_to_subdaily_failure_modes( + fixture_SubdailyScaler, inputs, outcome, msg +): + """Test fill_daily_to_subdaily using SubdailyScaler with method linear.""" + + # Set the included observations + fixture_SubdailyScaler.set_window( + window_center=np.timedelta64(13, "h"), half_width=np.timedelta64(1, "h") + ) + + with outcome as excep: + _ = fixture_SubdailyScaler.fill_daily_to_subdaily(**inputs) + + assert str(excep.value) == msg