diff --git a/AUTHORS.md b/AUTHORS.md index 49ecd1e246..45ab0b88bc 100644 --- a/AUTHORS.md +++ b/AUTHORS.md @@ -39,3 +39,4 @@ * Dahyann Araya * Giovanni Cozzolongo * Thomas Struys +* Dylan Pulver diff --git a/CHANGELOG.md b/CHANGELOG.md index 6b449f2225..5707d3fe66 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -31,6 +31,7 @@ Code freeze date: YYYY-MM-DD - `Hazard.from_raster_xarray` now returns a sparse matrix instead of a sparse array [#1261](https://github.com/CLIMADA-project/climada_python/pull/1261). - `ImpactCalc.impact` now raises a clear `ValueError` when the supplied `Hazard` contains no events, instead of failing later inside `np.array_split` with an obscure message [#814](https://github.com/CLIMADA-project/climada_python/issues/814). - Fix TCTracks.from_FAST duplicate loading from year loop [#1269](github.com/CLIMADA-project/climada_python/pull/1269) +- `compute_angular_windspeeds` no longer overrides `model_kwargs["cyclostrophic"]` with its own default, and no longer emits a `DeprecationWarning` when the deprecated `cyclostrophic` argument is not passed. The caller's `model_kwargs` dict is no longer modified in place. [#1209](https://github.com/CLIMADA-project/climada_python/issues/1209) ### Deprecated - `Impact.calc_freq_curve()` should not be given the parameter `return_per`. Use the parameter `return_periods` in `Impact.calc_freq_curve().interpolate()` instead. diff --git a/climada/hazard/test/test_trop_cyclone_windfields.py b/climada/hazard/test/test_trop_cyclone_windfields.py index 418e52867e..62e7947a54 100644 --- a/climada/hazard/test/test_trop_cyclone_windfields.py +++ b/climada/hazard/test/test_trop_cyclone_windfields.py @@ -1,4 +1,5 @@ import unittest +import warnings import numpy as np import xarray as xr @@ -9,6 +10,7 @@ H_TO_S, KM_TO_M, MBAR_TO_PA, + MODEL_VANG, _B_holland_1980, _bs_holland_2008, _stat_er_2011, @@ -17,6 +19,7 @@ _v_max_s_holland_2008, _vtrans, _x_holland_2010, + compute_angular_windspeeds, get_close_centroids, tctrack_to_si, ) @@ -26,6 +29,50 @@ class TestWindfieldHelpers(unittest.TestCase): """Test helper functions of TC wind field model""" + def setUp(self): + # Emanuel and Rotunno 2011 model test inputs + self.er_2011_inputs = self._er_2011_inputs() + + # Reference wind speeds for the er_2011 inputs above, with the Coriolis term + # (the default). + self.ER11_EXPECTED = [ + [28.258025, 36.782418, 36.869995, 22.521237], + [39.670883, 0.0, 3.300626, 10.827206], + ] + + # Same inputs under the cyclostrophic approximation, which drops the Coriolis term. + # From equation (36) of Emanuel and Rotunno 2011, + # M = M_max * 2 * (r/r_max)^2 / (1 + (r/r_max)^2) and v = M / r, with + # M_max = r_max * v_max = 1.6e6 m^2/s: at r = 30 km, + # v = 1.6e6 * 2 * 0.5625 / 1.5625 / 30e3 = 38.4 m/s. + self.ER11_NODE1_CYCLOSTROPHIC = [ + 38.4, + 0.0, + 3.194888178913738, + 10.480349344978167, + ] + + def _er_2011_inputs(self): + """Track/centroid inputs for the Emanuel and Rotunno 2011 tests. + + Centroids sit both within and outside the radius of maximum wind. + """ + d_centr = KM_TO_M * np.array( + [[35, 70, 75, 220], [30, 150, 1000, 300]], dtype=float + ) + si_track = xr.Dataset( + { + "rad": ("time", KM_TO_M * np.array([75.0, 40.0])), + "vmax": ("time", [35.0, 40.0]), + "lat": ("time", [20.0, 27.0]), + "cp": ("time", [4.98665369e-05, 6.61918149e-05]), + } + ) + mask = np.array( + [[True, True, True, True], [True, False, True, True]], dtype=bool + ) + return si_track, d_centr, mask + def test_get_close_centroids_pass(self): """Test get_close_centroids function.""" si_track = xr.Dataset( @@ -278,30 +325,66 @@ def test_stat_holland_1980(self): def test_er_2011_pass(self): """Test Emanuel and Rotunno 2011 wind field model.""" - # test at centroids within and outside of radius of max wind - d_centr = KM_TO_M * np.array( - [[35, 70, 75, 220], [30, 150, 1000, 300]], dtype=float - ) - si_track = xr.Dataset( - { - "rad": ("time", KM_TO_M * np.array([75.0, 40.0])), - "vmax": ("time", [35.0, 40.0]), - "lat": ("time", [20.0, 27.0]), - "cp": ("time", [4.98665369e-05, 6.61918149e-05]), - } + si_track, d_centr, mask = self.er_2011_inputs + v_ang_norm = _stat_er_2011(si_track, d_centr, mask) + np.testing.assert_allclose(v_ang_norm, self.ER11_EXPECTED) + + def test_er_2011_pass_cyclostrophic(self): + """Test Emanuel and Rotunno 2011 wind field model with cyclostrophic approx.""" + si_track, d_centr, mask = self.er_2011_inputs + v_ang_norm = _stat_er_2011(si_track, d_centr, mask, cyclostrophic=True) + np.testing.assert_allclose(v_ang_norm[1], self.ER11_NODE1_CYCLOSTROPHIC) + + def test_compute_angular_windspeeds_cyclostrophic_model_kwarg(self): + """``cyclostrophic`` passed via ``model_kwargs`` must reach the wind model.""" + si_track, d_centr, mask = self.er_2011_inputs + v_ang_norm = compute_angular_windspeeds( + si_track, + d_centr, + mask, + MODEL_VANG["ER11"], + model_kwargs={"cyclostrophic": True}, ) - mask = np.array( - [[True, True, True, True], [True, False, True, True]], dtype=bool + np.testing.assert_allclose( + v_ang_norm[1], self.ER11_NODE1_CYCLOSTROPHIC, atol=1e-6 ) - v_ang_norm = _stat_er_2011(si_track, d_centr, mask) - np.testing.assert_array_almost_equal( - v_ang_norm, - [ - [28.258025, 36.782418, 36.869995, 22.521237], - [39.670883, 0, 3.300626, 10.827206], - ], + + def test_compute_angular_windspeeds_no_spurious_deprecation(self): + """No DeprecationWarning unless the deprecated argument is actually passed.""" + si_track, d_centr, mask = self.er_2011_inputs + with warnings.catch_warnings(record=True) as caught: + warnings.simplefilter("always") + compute_angular_windspeeds( + si_track, d_centr, mask, MODEL_VANG["ER11"], model_kwargs={} + ) + self.assertEqual( + [w for w in caught if issubclass(w.category, DeprecationWarning)], [] ) + def test_compute_angular_windspeeds_cyclostrophic_deprecation(self): + """Passing the deprecated ``cyclostrophic`` argument warns.""" + si_track, d_centr, mask = self.er_2011_inputs + with self.assertWarns(DeprecationWarning): + compute_angular_windspeeds( + si_track, d_centr, mask, MODEL_VANG["ER11"], cyclostrophic=True + ) + + def test_compute_angular_windspeeds_does_not_mutate_model_kwargs(self): + """The caller's ``model_kwargs`` dict must not be modified in place.""" + si_track, d_centr, mask = self.er_2011_inputs + model_kwargs = {} + with warnings.catch_warnings(): + warnings.simplefilter("ignore", DeprecationWarning) + compute_angular_windspeeds( + si_track, + d_centr, + mask, + MODEL_VANG["ER11"], + cyclostrophic=True, + model_kwargs=model_kwargs, + ) + self.assertEqual(model_kwargs, {}) + def test_vtrans_pass(self): """Test _vtrans function. Compare to MATLAB reference.""" tc_track = TCTracks.from_processed_ibtracs_csv(TEST_TRACK) diff --git a/climada/hazard/trop_cyclone/trop_cyclone_windfields.py b/climada/hazard/trop_cyclone/trop_cyclone_windfields.py index 86a2e144a2..aa59d96c49 100644 --- a/climada/hazard/trop_cyclone/trop_cyclone_windfields.py +++ b/climada/hazard/trop_cyclone/trop_cyclone_windfields.py @@ -100,7 +100,7 @@ def compute_angular_windspeeds( d_centr: np.ndarray, mask_centr_close: np.ndarray, model: int, - cyclostrophic: Optional[bool] = False, + cyclostrophic: Optional[bool] = None, model_kwargs: Optional[dict] = None, ): """Compute (absolute) angular wind speeds according to a parametric wind profile @@ -120,14 +120,16 @@ def compute_angular_windspeeds( If given, forward these kwargs to the selected model. Default: None cyclostrophic: bool, optional, deprecated This argument is deprecated and will be removed in a future release. - Include `cyclostrophic` as `model_kwargs` instead. + Include `cyclostrophic` as `model_kwargs` instead. If given, it takes precedence + over the value in `model_kwargs`. Default: None (not set). Returns ------- ndarray of shape (npositions, ncentroids) containing the magnitude of the angular windspeed per track position per centroid location """ - model_kwargs = {} if model_kwargs is None else model_kwargs + # copy, so that the deprecation shim below never mutates the caller's dict + model_kwargs = {} if model_kwargs is None else dict(model_kwargs) if cyclostrophic is not None: warnings.warn(