Merge remote-tracking branch 'origin/feature/fix_knutson_scaling' into feature/tc-synth-new-model-euler

Author: spjuhel

CHANGELOG.md

@@ -6,6 +6,22 @@ Release date: 2024-07-19
 
 ### Dependency Changes
 
+### Added
+
+### Changed
+
+### Fixed
+
+### Deprecated
+
+### Removed
+
+## 5.0.0
+
+Release date: 2024-07-19
+
+### Dependency Changes
+
 Added:
 
 - `bayesian-optimization`

climada/_version.py

@@ -1 +1 @@
-__version__ = '5.0.0'
+__version__ = '5.0.1-dev'

climada/hazard/tc_clim_change.py

@@ -334,11 +334,11 @@ def from_ibtracs_netcdf(cls, provider=None, rescale_windspeeds=True, storm_id=No
 
         When using data from IBTrACS, make sure to be familiar with the scope and limitations of
         IBTrACS, e.g. by reading the official documentation
-        (https://www.ncdc.noaa.gov/ibtracs/pdf/IBTrACS_version4_Technical_Details.pdf). Reading the
-        CLIMADA documentation can't replace a thorough understanding of the underlying data. This
-        function only provides a (hopefully useful) interface for the data input, but cannot
-        provide any guidance or make recommendations about if and how to use IBTrACS data for your
-        particular project.
+        (https://www.ncei.noaa.gov/sites/default/files/2021-07/IBTrACS_version4_Technical_Details.pdf).
+        Reading the CLIMADA documentation can't replace a thorough understanding of the underlying
+        data. This function only provides a (hopefully useful) interface for the data input, but
+        cannot provide any guidance or make recommendations about if and how to use IBTrACS data
+        for your particular project.
 
         Resulting tracks are required to have both pressure and wind speed information at all time
         steps. Therefore, all track positions where one of wind speed or pressure are missing are
@@ -377,8 +377,8 @@ def from_ibtracs_netcdf(cls, provider=None, rescale_windspeeds=True, storm_id=No
         rescale_windspeeds : bool, optional
             If True, all wind speeds are linearly rescaled to 1-minute sustained winds.
             Note however that the IBTrACS documentation (Section 5.2,
-            https://www.ncdc.noaa.gov/ibtracs/pdf/IBTrACS_version4_Technical_Details.pdf) includes
-            a warning about this kind of conversion: "While a multiplicative factor can describe
+            https://www.ncei.noaa.gov/sites/default/files/2021-07/IBTrACS_version4_Technical_Details.pdf)
+            includes a warning about this kind of conversion: "While a multiplicative factor can
             the numerical differences, there are procedural and observational differences between
             agencies that can change through time, which confounds the simple multiplicative
             factor." Default: True

climada/hazard/tc_tracks.py

@@ -21,45 +21,127 @@
 
 import unittest
 
+import unittest
+from math import log
+import pandas as pd
+import numpy as np
 import climada.hazard.tc_clim_change as tc_cc
 
 class TestKnutson(unittest.TestCase):
-    """Test loading funcions from the TropCyclone class"""
 
-    def test_get_pass(self):
+    def test_get_knutson_scaling_calculations(self):
+
+        basin = 'NA'
+        variable = 'cat05'
+        percentile = '5/10'
+        base_start, base_end = 1950, 2018
+        yearly_steps = 5
+
+        predicted_changes = tc_cc.get_knutson_scaling_factor(
+                    percentile=percentile,
+                    variable=variable,
+                    basin=basin,
+                    baseline=(base_start, base_end),
+                    yearly_steps=yearly_steps
+        )
+
+        ## Test computations of future changes
+        # Load data
+        gmst_info = tc_cc.get_gmst_info()
+
+        var_id, basin_id, perc_id = (tc_cc.MAP_VARS_NAMES[variable],
+                                    tc_cc.MAP_BASINS_NAMES[basin],
+                                    tc_cc.MAP_PERC_NAMES[percentile])
+
+        knutson_data = tc_cc.get_knutson_data()
+        knutson_value = knutson_data[var_id, basin_id, perc_id]
+    
+        start_ind = base_start - gmst_info['gmst_start_year']
+        end_ind = base_end - gmst_info['gmst_start_year']
+
+        # Apply model
+        beta = 0.5 * log(0.01 * knutson_value + 1)
+        tc_properties = np.exp(beta * gmst_info['gmst_data'])
+
+        # Assess baseline value
+        baseline = np.mean(tc_properties[:, start_ind:end_ind + 1], 1)
+
+        # Assess future value and test predicted change from baseline is
+        # the same as given by function
+        smoothing = 5
+
+        for target_year in [2030, 2050, 2070, 2090]:
+            target_year_ind = target_year - gmst_info['gmst_start_year']
+            ind1 = target_year_ind - smoothing
+            ind2 = target_year_ind + smoothing + 1
+
+            prediction = np.mean(tc_properties[:, ind1:ind2], 1)
+            predicted_change = ((prediction - baseline) / baseline) * 100
+
+            np.testing.assert_array_almost_equal(predicted_changes.loc[target_year, '2.6'], predicted_change[0])
+            np.testing.assert_array_almost_equal(predicted_changes.loc[target_year, '4.5'], predicted_change[1])
+            np.testing.assert_array_almost_equal(predicted_changes.loc[target_year, '6.0'], predicted_change[2])
+            np.testing.assert_array_almost_equal(predicted_changes.loc[target_year, '8.5'], predicted_change[3])
+
+    def test_get_knutson_scaling_structure(self):
         """Test get_knutson_criterion function."""
-        criterion = tc_cc.get_knutson_criterion()
-        self.assertTrue(len(criterion), 20)
-        for crit_val in criterion:
-            self.assertTrue('year' in crit_val)
-            self.assertTrue('change' in crit_val)
-            self.assertTrue('variable' in crit_val)
-        self.assertEqual(criterion[0]['variable'], "frequency")
-        self.assertEqual(criterion[0]['change'], 1)
-        self.assertEqual(criterion[4]['variable'], "intensity")
-        self.assertEqual(criterion[4]['change'], 1.045)
-        self.assertEqual(criterion[-10]['basin'], "SP")
-        self.assertEqual(criterion[-10]['variable'], "frequency")
-        self.assertEqual(criterion[-10]['change'], 1 - 0.583)
-
-    def test_scale_pass(self):
-        """Test calc_scale_knutson function."""
-        self.assertAlmostEqual(tc_cc.calc_scale_knutson(ref_year=2050, rcp_scenario=45),
-                               0.759630751756698)
-        self.assertAlmostEqual(tc_cc.calc_scale_knutson(ref_year=2070, rcp_scenario=45),
-                               0.958978483788876)
-        self.assertAlmostEqual(tc_cc.calc_scale_knutson(ref_year=2060, rcp_scenario=60),
-                               0.825572149523299)
-        self.assertAlmostEqual(tc_cc.calc_scale_knutson(ref_year=2080, rcp_scenario=60),
-                               1.309882943406079)
-        self.assertAlmostEqual(tc_cc.calc_scale_knutson(ref_year=2090, rcp_scenario=85),
-                               2.635069196605717)
-        self.assertAlmostEqual(tc_cc.calc_scale_knutson(ref_year=2100, rcp_scenario=85),
-                               2.940055236533517)
-        self.assertAlmostEqual(tc_cc.calc_scale_knutson(ref_year=2066, rcp_scenario=26),
-                               0.341930203294547)
-        self.assertAlmostEqual(tc_cc.calc_scale_knutson(ref_year=2078, rcp_scenario=26),
-                               0.312383928930456)
+        
+        yearly_steps = 8
+        predicted_changes = tc_cc.get_knutson_scaling_factor(yearly_steps=yearly_steps)
+
+        np.testing.assert_equal(predicted_changes.columns, np.array(['2.6', '4.5', '6.0', '8.5']))
+
+        simulated_years = np.arange(tc_cc.YEAR_WINDOWS_PROPS['start'],
+                                    tc_cc.YEAR_WINDOWS_PROPS['end']+1,
+                                    yearly_steps)
+        np.testing.assert_equal(predicted_changes.index, simulated_years)
+
+    def test_get_knutson_scaling_valid_inputs(self):
+        df = tc_cc.get_knutson_scaling_factor()
+        self.assertIsInstance(df, pd.DataFrame)
+        np.testing.assert_equal(df.shape, (21, 4))
+
+    def test_get_knutson_scaling_invalid_baseline_start_year(self):
+        with self.assertRaises(ValueError):
+            tc_cc.get_knutson_scaling_factor(baseline=(1870, 2022))
+
+    def test_get_knutson_scaling_invalid_baseline_end_year(self):
+        with self.assertRaises(ValueError):
+            tc_cc.get_knutson_scaling_factor(baseline=(1982, 2110))
+
+    def test_get_knutson_scaling_no_scaling_factors_for_unknonw_basin(self):
+        df = tc_cc.get_knutson_scaling_factor(basin='ZZZZZ')
+        self.assertIsInstance(df, pd.DataFrame)
+        np.testing.assert_equal(df.values, np.ones_like(df.values))
+
+    def test_get_gmst(self):
+        """Test get_gmst_info function."""
+        gmst_info = tc_cc.get_gmst_info()
+
+        self.assertAlmostEqual(gmst_info['gmst_start_year'], 1880)
+        self.assertAlmostEqual(gmst_info['gmst_end_year'], 2100)
+        self.assertAlmostEqual(len(gmst_info['rcps']), 4)
+
+        self.assertAlmostEqual(gmst_info['gmst_data'].shape,
+                               (len(gmst_info['rcps']),
+                                gmst_info['gmst_end_year']-gmst_info['gmst_start_year']+1))
+        self.assertAlmostEqual(gmst_info['gmst_data'][0,0], -0.16)
+        self.assertAlmostEqual(gmst_info['gmst_data'][0,-1], 1.27641, 4)
+        self.assertAlmostEqual(gmst_info['gmst_data'][-1,0], -0.16)
+        self.assertAlmostEqual(gmst_info['gmst_data'][-1,-1], 4.477764, 4)
+
+    def test_get_knutson_data_pass(self):
+        """Test get_knutson_data function."""
+
+        data_knutson = tc_cc.get_knutson_data()
+
+        self.assertAlmostEqual(data_knutson.shape, (4,6,5))
+        self.assertAlmostEqual(data_knutson[0,0,0], -34.49)
+        self.assertAlmostEqual(data_knutson[-1,-1,-1], 15.419)
+        self.assertAlmostEqual(data_knutson[0,-1,-1], 4.689)
+        self.assertAlmostEqual(data_knutson[-1,0,0], 5.848)
+        self.assertAlmostEqual(data_knutson[-1,0,-1], 22.803)
+        self.assertAlmostEqual(data_knutson[2,3,2], 4.324)
 
 if __name__ == "__main__":
     TESTS = unittest.TestLoader().loadTestsFromTestCase(TestKnutson)