Hazard: winter windstorms / extratropical cyclones in Europe#
Or: The StormEurope hazard subclass of CLIMADA#
Auth: Jan Hartman & Thomas Röösli
Date: 2018-04-26 & 2020-03-03
This notebook will give a quick tour of the capabilities of the StormEurope hazard class. This includes functionalities to apply probabilistic alterations to historical storms.
%matplotlib inline
import matplotlib.pyplot as plt
plt.rcParams["figure.figsize"] = [15, 10]
Reading Data#
StormEurope was written under the presumption that you’d start out with WISC storm footprint data in netCDF format. This notebook works with a demo dataset. If you would like to work with the real data: (1) Please follow the link and download the file C3S_WISC_FOOTPRINT_NETCDF_0100.tgz from the Copernicus Windstorm Information Service, (2) unzip it (3) uncomment the last two lines in the following codeblock and (4) adjust the variable “WISC_files”.
We first construct an instance and then point the reader at a directory containing compatible .nc
files. Since there are other files in there, we must be explicit and use a globbing pattern; supplying incompatible files will make the reader fail.
The reader actually calls climada.util.files_handler.get_file_names
, so it’s also possible to hand it an explicit list of filenames, or a dirname, or even a list of glob patterns or directories.
from climada.hazard import StormEurope
from climada.util.constants import WS_DEMO_NC
storm_instance = StormEurope.from_footprints(WS_DEMO_NC)
# WISC_files = '/path/to/folder/C3S_WISC_FOOTPRINT_NETCDF_0100/fp_era[!er5]*_0.nc'
# storm_instance = StormEurope.from_footprints(WISC_files)
$CLIMADA_SRC/climada/hazard/centroids/centr.py:822: UserWarning: Geometry is in a geographic CRS. Results from 'buffer' are likely incorrect. Use 'GeoSeries.to_crs()' to re-project geometries to a projected CRS before this operation.
xy_pixels = self.geometry.buffer(res / 2).envelope
Introspection#
Let’s quickly see what attributes this class brings with it:
?storm_instance
Type: StormEurope
String form: <climada.hazard.storm_europe.StormEurope object at 0x7f2a986b4c70>
File: ~/code/climada_python/climada/hazard/storm_europe.py
Docstring:
A hazard set containing european winter storm events. Historic storm
events can be downloaded at http://wisc.climate.copernicus.eu/ and read
with `from_footprints`. Weather forecasts can be automatically downloaded from
https://opendata.dwd.de/ and read with from_icon_grib(). Weather forecast
from the COSMO-Consortium http://www.cosmo-model.org/ can be read with
from_cosmoe_file().
Attributes
----------
ssi_wisc : np.array, float
Storm Severity Index (SSI) as recorded in
the footprint files; apparently not reproducible from the footprint
values only.
ssi : np.array, float
SSI as set by set_ssi; uses the Dawkins
definition by default.
Init docstring: Calls the Hazard init dunder. Sets unit to 'm/s'.
You could also try listing all permissible methods with dir(storm_instance)
, but since that would include the methods from the Hazard
base class, you wouldn’t know what’s special. The best way is to read the source: uncomment the following statement to read more.
# StormEurope??
Into the Storm Severity Index (SSI)#
The SSI, according to Dawkins et al. 2016 or Lamb and Frydendahl, 1991, can be set using set_ssi
. For demonstration purposes, I show the default arguments. (Check also the defaults using storm_instance.calc_ssi?
, the method for which set_ssi is a wrapper.)
We won’t be using the plot_ssi
functionality just yet, because we only have two events; the graph really isn’t informative. After this, we’ll generate some more storms to make that plot more aesthetically pleasing.
storm_instance.set_ssi(
method="wind_gust",
intensity=storm_instance.intensity,
# the above is just a more explicit way of passing the default
on_land=True,
threshold=25,
sel_cen=None,
# None is default. sel_cen could be used to subset centroids
)
Probabilistic Storms#
This class allows generating probabilistic storms from historical ones according to a method outlined in Schwierz et al. 2010. This means that per historical event, we generate 29 new ones with altered intensities. Since it’s just a bunch of vector operations, this is pretty fast.
However, we should not return the entire probabilistic dataset in-memory: in trials, this used up 60 GB of RAM, thus requiring a great amount of swap space. Instead, we must select a country by setting the reg_id
parameter to an ISO_N3 country code used in the Natural Earth dataset. It is also possible to supply a list of ISO codes. If your machine is up for the job of handling the whole dataset, set the reg_id
parameter to None
.
Since assigning each centroid a country ID is a rather inefficient affair, you may need to wait a minute or two for the entire WISC dataset to be processed. For the small demo dataset, it runs pretty quickly.
%%time
storm_prob = storm_instance.generate_prob_storms(reg_id=528)
storm_prob.plot_intensity(0);
2020-03-05 10:29:31,845 - climada.hazard.centroids.centr - INFO - Setting geometry points.
2020-03-05 10:29:32,248 - climada.hazard.centroids.centr - DEBUG - Setting region_id 9944 points.
2020-03-05 10:29:32,466 - climada.util.coordinates - DEBUG - Setting region_id 9944 points.
2020-03-05 10:29:33,506 - climada.hazard.storm_europe - INFO - Commencing probabilistic calculations
2020-03-05 10:29:33,620 - climada.hazard.storm_europe - INFO - Generating new StormEurope instance
2020-03-05 10:29:33,663 - climada.util.checker - DEBUG - Hazard.ssi not set.
2020-03-05 10:29:33,664 - climada.util.checker - DEBUG - Hazard.ssi_wisc not set.
2020-03-05 10:29:33,665 - climada.util.checker - DEBUG - Hazard.event_name not set. Default values set.
C:\shortpaths\GitHub\climada_python\climada\util\plot.py:311: UserWarning: Tight layout not applied. The left and right margins cannot be made large enough to accommodate all axes decorations.
fig.tight_layout()
Wall time: 2.24 s
<cartopy.mpl.geoaxes.GeoAxesSubplot at 0x1dafba69940>
We can get much more fancy in our calls to generate_prob_storms
; the keyword arguments after ssi_args
are passed on to _hist2prob
, allowing us to tweak the probabilistic permutations.
ssi_args = {
"on_land": True,
"threshold": 25,
}
storm_prob_xtreme = storm_instance.generate_prob_storms(
reg_id=[56, 528], # BEL and NLD
spatial_shift=2,
ssi_args=ssi_args,
power=1.5,
scale=0.3,
)
We can now check out the SSI plots of both these calculations. The comparison between the historic and probabilistic ssi values, only makes sense for the full dataset.
storm_prob_xtreme.plot_ssi(full_area=True)
storm_prob.plot_ssi(full_area=True);
(<Figure size 1080x720 with 1 Axes>,
<AxesSubplot:xlabel='Exceedance Frequency [1/a]', ylabel='Storm Severity Index'>)