Skip to content

Remove intensity_thres attribute from Hazard class #1065

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.

Sign up for GitHub

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

Merged
ValentinGebhart merged 11 commits into from
Jul 31, 2025
Merged

Remove intensity_thres attribute from Hazard class #1065

Show file tree
Hide file tree
Changes from 5 commits
Commits
Show all changes
11 commits
Select commit Hold shift + click to select a range
08cfe7b
removed attribute in hazard.base
ValentinGebhart May 15, 2025
4f21ae0
remove attribute in hazard.test_base
ValentinGebhart Jul 1, 2025
658d5ec
adapted storm_europe.py to similar structure as trop_cyclones
ValentinGebhart Jul 1, 2025
6f0969d
removed global hazard intensity_thres attribute
ValentinGebhart Jul 1, 2025
c61b780
Merge branch 'develop' into feature/remove_intensity_thres
ValentinGebhart Jul 1, 2025
de7fc0d
include attribute in class docstrings
ValentinGebhart Jul 1, 2025
b0189ba
update docstring
ValentinGebhart Jul 17, 2025
cf1a5c4
update changelog
ValentinGebhart Jul 17, 2025
3468a7e
removed class attributes from StormEurope and TropCyclone and removed…
ValentinGebhart Jul 22, 2025
2eb02af
Apply suggestions from code review
ValentinGebhart Jul 28, 2025
9470b83
Merge branch 'develop' into feature/remove_intensity_thres
ValentinGebhart Jul 31, 2025
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
18 changes: 4 additions & 14 deletions climada/hazard/base.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -101,10 +101,6 @@ class Hazard(HazardIO, HazardPlot):
(i.e., is equivalent to fraction is 1 everywhere).
"""

intensity_thres = 10
"""Intensity threshold per hazard used to filter lower intensities. To be
set for every hazard type"""

vars_oblig = {
"units",
"centroids",
Expand Down Expand Up @@ -485,7 +481,7 @@ def local_exceedance_intensity(
self,
return_periods=(25, 50, 100, 250),
method="interpolate",
min_intensity=None,
min_intensity=0,

@chahank chahank Jul 1, 2025

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Does this mean that all sparse data is projected onto a dense matrix? This might lead to problems no?

@ValentinGebhart ValentinGebhart Jul 1, 2025

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The calculation of local exceedance intensity is centroid wise, so yes, the sparse matrix is converted to a dense matrix column by column. The min_intensity parameter being zero or positive does not change this.

It might be that a positive min_intensity parameter here would result in larger geodataframes (ie. centroids with only very small intensities would still lead to a computation and nonzero values, instead of faster handling if all intensities are considered zero). However, as we do not have a minimum intensity for a general hazard object, I think this should be set to zero per default, and can be adapted by the user if wanted.

Does this answer your question?

@chahank chahank Jul 1, 2025

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, thanks!

log_frequency=True,
log_intensity=True,
bin_decimals=None,
Expand All @@ -512,8 +508,7 @@ def local_exceedance_intensity(
periods and extends the interpolation between these points to the given return period
(similar for small return periods). Defauls to "interpolate".
min_intensity : float, optional
Minimum threshold to filter the hazard intensity. If set to None, self.intensity_thres
will be used. Defaults to None.
Minimum threshold to filter the hazard intensity. Defaults to 0.
log_frequency : bool, optional
If set to True, (cummulative) frequency values are converted to log scale before
inter- and extrapolation. Defaults to True.
Expand Down Expand Up @@ -554,8 +549,6 @@ def local_exceedance_intensity(
intensities range from 1e6 to 1e9, you could use bin_decimals=-5, if your intensities
range from 0.0001 to .01, you could use bin_decimals=5.
"""
if not min_intensity and min_intensity != 0:
min_intensity = self.intensity_thres
# check frequency unit
return_period_unit = u_dt.convert_frequency_unit_to_time_unit(
self.frequency_unit
Expand Down Expand Up @@ -639,7 +632,7 @@ def local_return_period(
self,
threshold_intensities=(10.0, 20.0),
method="interpolate",
min_intensity=None,
min_intensity=0,
log_frequency=True,
log_intensity=True,
bin_decimals=None,
Expand Down Expand Up @@ -667,8 +660,7 @@ def local_return_period(
points to the given threshold intensity (similar for small threshold intensites).
Defaults to "interpolate".
min_intensity : float, optional
Minimum threshold to filter the hazard intensity. If set to None, self.intensity_thres
will be used. Defaults to None.
Minimum threshold to filter the hazard intensity. Defaults to 0.
log_frequency : bool, optional
If set to True, (cummulative) frequency values are converted to log scale before
inter- and extrapolation. Defaults to True.
Expand Down Expand Up @@ -710,8 +702,6 @@ def local_return_period(
intensities range from 1e6 to 1e9, you could use bin_decimals=-5, if your intensities
range from 0.0001 to .01, you could use bin_decimals=5.
"""
if not min_intensity and min_intensity != 0:
min_intensity = self.intensity_thres
# check frequency unit
return_period_unit = u_dt.convert_frequency_unit_to_time_unit(
self.frequency_unit
Expand Down
37 changes: 17 additions & 20 deletions climada/hazard/storm_europe.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -57,6 +57,12 @@
N_PROB_EVENTS = 5 * 6
"""Number of events per historic event in probabilistic dataset"""

DEF_INTENSITY_THRES = 14.7
"""
Default value for the threshold below which wind speeds (in m/s) are stored as 0.
Same as used by WISC SSI calculations.
"""


class StormEurope(Hazard):
"""A hazard set containing european winter storm events. Historic storm
Expand All @@ -81,8 +87,8 @@
SSI as set by set_ssi; uses the Dawkins definition by default.
"""

intensity_thres = 14.7
"""Intensity threshold for storage in m/s; same as used by WISC SSI calculations."""
intensity_thres = DEF_INTENSITY_THRES
"""Intensity threshold for storage in m/s."""

@chahank chahank Jul 21, 2025

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sorry for the late extra comment, but is the class attribute then still needed? It looks like the intensity_thresh is now a parameter of all relevant methods no?

@ValentinGebhart ValentinGebhart Jul 21, 2025

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, intensity_thresh is an optional parameter of all relevant methods. This class attribute is used as the default values for these methods. If you change the attribute for a specific instance, the default value of all corresponding methods changes, and you do not have to change it manually each time.

Does this make sense? If there is a smarter implementation, happy to change this.

@chahank chahank Jul 21, 2025

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It makes sense, although I am not sure the current implementation is doing what you want.

The code now has an class attribute, but actually never uses it, since the value of this attribute is always passed as parameter to all relevant methods. Thus, the class attribute seems redundant to me. What do you think?

@ValentinGebhart ValentinGebhart Jul 22, 2025

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, almost agree. There is one place in StormEurope and one in TropCyclone where the correpsonding class attributes are used (called on class instances though):
In StormEurope: in the function calc_ssi
In TropCyclone: in the function set_from_tracks which has a deprecation warning.

It looks to me as if we can remove these uses and then also remove the class attributes. I have to think about the calc_ssi though. If you confirm that it is better to remove them, I can try to implement that :)

@chahank chahank Jul 22, 2025

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

How about we add the parameter to the methods calc_ssi, remove the method set_from_tracks, adjust the method from_tracks (remove the haz.intensity_thres = intensity_thres), and remove the class parameter?

@ValentinGebhart ValentinGebhart Jul 22, 2025 •
edited
Loading

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I removed the class attribute from StormEurope and TropCyclone, and the deprecated TropCyclone.set_from_tracks. In the StormEurope.calc_ssi method, there was already an optional parameter threshold for the calculation of SSI. Before, the method checked if this parameter was larger than the class attribute intensity_thres and, if not, gave an error "threshold cannot be below threshold upon read_footprint". I don't understand why it should be a problem though if the calculation threshold is below the read-in threshold...

Let me know what you think, for me this looks good.


vars_opt = Hazard.vars_opt.union({"ssi_wisc", "ssi", "ssi_full_area"})
"""Name of the variables that aren't need to compute the impact."""
Expand Down Expand Up @@ -139,7 +145,7 @@
centroids=None,
files_omit="fp_era20c_1990012515_701_0.nc",
combine_threshold=None,
intensity_thres=None,
intensity_thres=DEF_INTENSITY_THRES,
):
"""Create new StormEurope object from WISC footprints.

Expand Down Expand Up @@ -171,18 +177,15 @@
events are combined into one.
Default is None, Advised for WISC is 2
intensity_thres : float, optional
Intensity threshold for storage in m/s. Default: class attribute
StormEurope.intensity_thres (same as used by WISC SSI calculations)
Intensity threshold for storage in m/s. Default: 14.7
(same as used by WISC SSI calculations)

Returns
-------
haz : StormEurope
StormEurope object with data from WISC footprints.
"""
# pylint: disable=protected-access
intensity_thres = (
cls.intensity_thres if intensity_thres is None else intensity_thres
)
file_names = get_file_names(path)

if ref_raster is not None and centroids is not None:
Expand Down Expand Up @@ -300,7 +303,7 @@
event_date=None,
model_name="COSMO-2E",
description=None,
intensity_thres=None,
intensity_thres=DEF_INTENSITY_THRES,
):
"""Create a new StormEurope object with gust footprint from weather forecast.

Expand Down Expand Up @@ -331,17 +334,14 @@
description of the events, defaults
to a combination of model_name and run_datetime
intensity_thres : float, optional
Intensity threshold for storage in m/s. Default: class attribute
StormEurope.intensity_thres (same as used by WISC SSI calculations)
Intensity threshold for storage in m/s. Default: 14.7
(same as used by WISC SSI calculations)

Returns
-------
haz : StormEurope
StormEurope object with data from COSMO ensemble file.
"""
intensity_thres = (
cls.intensity_thres if intensity_thres is None else intensity_thres
)

# read intensity from file
with xr.open_dataset(fp_file) as ncdf:
Expand Down Expand Up @@ -435,7 +435,7 @@
description=None,
grib_dir=None,
delete_raw_data=True,
intensity_thres=None,
intensity_thres=DEF_INTENSITY_THRES,
):
"""Create new StormEurope object from DWD icon weather forecast footprints.

Expand Down Expand Up @@ -472,20 +472,17 @@
.grib.bz2 file format should be stored on the computer or
removed
intensity_thres : float, optional
Intensity threshold for storage in m/s. Default: class attribute
StormEurope.intensity_thres (same as used by WISC SSI calculations)
Intensity threshold for storage in m/s. Default: 14.7
(same as used by WISC SSI calculations)

Returns
-------
haz : StormEurope
StormEurope object with data from DWD icon weather forecast footprints.
"""
# pylint: disable=protected-access
intensity_thres = (
cls.intensity_thres if intensity_thres is None else intensity_thres
)

if not (run_datetime.hour == 0 or run_datetime.hour == 12):

Check warning on line 485 in climada/hazard/storm_europe.py

View check run for this annotation

Jenkins - WCR / Pylint

consider-using-in 

LOW:
Consider merging these comparisons with 'in' by using 'run_datetime.hour in (0, 12)'. Use a set instead if elements are hashable.
Raw output 
no description found
LOGGER.warning(
"The event definition is inaccuratly implemented "
"for starting times, which are not 00H or 12H."
Expand Down
21 changes: 13 additions & 8 deletions climada/hazard/test/test_base.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1156,7 +1156,6 @@ def test_local_exceedance_intensity(self):
"""Test local exceedance frequencies with lin lin interpolation"""
haz = dummy_hazard()
haz.intensity = sparse.csr_matrix([[1.0, 3.0, 1.0], [2.0, 3.0, 0.0]])
haz.intensity_thres = 0.5
haz.frequency = np.full(2, 1.0)
return_period = np.array([0.5, 2.0 / 3.0, 1.0])
# first centroid has intensities 1,2 with cum frequencies 2,1
Expand All @@ -1180,17 +1179,19 @@ def test_local_exceedance_intensity_methods(self):
haz.intensity = sparse.csr_matrix(
[[0, 0, 1e1], [0.2, 1e1, 1e2], [1e3, 1e3, 1e3]]
)
haz.intensity_thres = 0.5
min_intensity = 0.5
haz.frequency = np.array([1.0, 0.1, 0.01])
return_period = (1000, 30, 0.1)
return_periods = (1000, 30, 0.1)
# first centroid has intensities 1e3 with frequencies .01, cum freq .01
# second centroid has intensities 1e1, 1e3 with cum frequencies .1, .01, cum freq .11, .01
# third centroid has intensities 1e1, 1e2, 1e3 with cum frequencies 1., .1, .01, cum freq 1.11, .11, .01
# testing at frequencies .001, .033, 10.

# test stepfunction
inten_stats, _, _ = haz.local_exceedance_intensity(
return_periods=(1000, 30, 0.1), method="stepfunction"
return_periods=return_periods,
method="stepfunction",
min_intensity=min_intensity,
)
np.testing.assert_allclose(
inten_stats.values[:, 1:].astype(float),
Expand All @@ -1199,7 +1200,9 @@ def test_local_exceedance_intensity_methods(self):

# test log log extrapolation
inten_stats, _, _ = haz.local_exceedance_intensity(
return_periods=(1000, 30, 0.1), method="extrapolate"
return_periods=return_periods,
method="extrapolate",
min_intensity=min_intensity,
)
np.testing.assert_allclose(
inten_stats.values[:, 1:].astype(float),
Expand All @@ -1209,7 +1212,9 @@ def test_local_exceedance_intensity_methods(self):

# test log log interpolation and extrapolation with constant
inten_stats, _, _ = haz.local_exceedance_intensity(
return_periods=(1000, 30, 0.1), method="extrapolate_constant"
return_periods=return_periods,
method="extrapolate_constant",
min_intensity=min_intensity,
)
np.testing.assert_allclose(
inten_stats.values[:, 1:].astype(float),
Expand All @@ -1219,7 +1224,7 @@ def test_local_exceedance_intensity_methods(self):

# test log log interpolation and no extrapolation
inten_stats, _, _ = haz.local_exceedance_intensity(
return_periods=(1000, 30, 0.1)
return_periods=return_periods, min_intensity=min_intensity
)
np.testing.assert_allclose(
inten_stats.values[:, 1:].astype(float),
Expand All @@ -1235,6 +1240,7 @@ def test_local_exceedance_intensity_methods(self):
log_frequency=False,
log_intensity=False,
method="extrapolate_constant",
min_intensity=min_intensity,
)
np.testing.assert_allclose(
inten_stats.values[:, 1:].astype(float),
Expand Down Expand Up @@ -1270,7 +1276,6 @@ def test_local_return_period_methods(self):
haz.intensity = sparse.csr_matrix(
[[0, 0, 1e1], [0.0, 1e1, 1e2], [1e3, 1e3, 1e3]]
)
haz.intensity_thres = 0.5
haz.frequency = np.array([1.0, 0.1, 0.01])
# first centroid has intensities 1e3 with frequencies .01, cum freq .01
# second centroid has intensities 1e1, 1e3 with cum frequencies .1, .01, cum freq .11, .01
Expand Down
3 changes: 1 addition & 2 deletions doc/user-guide/climada_util_local_exceedance_values.ipynb
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -40,7 +40,7 @@
},
{
"cell_type": "code",
"execution_count": 1,
"execution_count": null,
"metadata": {},
"outputs": [
{
Expand Down Expand Up @@ -86,7 +86,6 @@
" frequency_unit=\"1/year\",\n",
" units=\"m/s\",\n",
")\n",
"hazard.intensity_thres = 0\n",
"\n",
"# plot first event of Hazard object\n",
"hazard.plot_intensity(event=1, smooth=False, figsize=(4, 4));"
Expand Down
Loading