Hazard class#

What is a hazard?#

A hazard describes weather events such as storms, floods, droughts, or heat waves both in terms of probability of occurrence as well as physical intensity.


How are hazards embedded in the CLIMADA architecture?#

Hazards are defined by the base class Hazard which gathers the required attributes that enable the impact computation (such as centroids, frequency per event, and intensity per event and centroid) and common methods such as readers and visualization functions. Each hazard class collects historical data or model simulations and transforms them, if necessary, in order to construct a coherent event database. Stochastic events can be generated taking into account the frequency and main intensity characteristics (such as local water depth for floods or gust speed for storms) of historical events, producing an ensemble of probabilistic events for each historical event. CLIMADA provides therefore an event-based probabilistic approach which does not depend on a hypothesis of a priori general probability distribution choices. Note that one can also reduce the probabilistic approach to a deterministic approach (e.g., story-line or forecasting) by defining the frequency to be 1. The source of the historical data (e.g. inventories or satellite images) or model simulations (e.g. synthetic tropical cyclone tracks) and the methodologies used to compute the hazard attributes and its stochastic events depend on each hazard type and are defined in its corresponding Hazard-derived class (e.g. TropCylcone for tropical cyclones, explained in the tutorial TropCyclone). This procedure provides a solid and homogeneous methodology to compute impacts worldwide. In the case where the risk analysis comprises a specific region where good quality data or models describing the hazard intensity and frequency are available, these can be directly ingested by the platform through the reader functions, skipping the hazard modelling part (in total or partially), and allowing us to easily and seamlessly combine CLIMADA with external sources. Hence the impact model can be used for a wide variety of applications, e.g. deterministically to assess the impact of a single (past or future) event or to quantify risk based on a (large) set of probabilistic events. Note that since the Hazard class is not an abstract class, any hazard that is not defined in CLIMADA can still be used by providing the Hazard attributes.


What do hazards look like in CLIMADA?#

A Hazard contains events of some hazard type defined at centroids. There are certain variables in a Hazard instance that are needed to compute the impact, while others are descriptive and can therefore be set with default values. The full list of looks like this:

Mandatory variables             

Data Type                                    

Description                                        

units

(str)

units of the intensity

centroids

Centroids()

centroids of the events

event_id

(np.array)

id (>0) of each event

frequency

(np.array)

frequency of each event in years

intensity

(sparse.csr_matrix)

intensity of the events at centroids

fraction

(sparse.csr_matrix)

fraction of affected exposures for each event at each centroid



Descriptive variables

Data Type

Description

date

(np.array)

integer date corresponding to the proleptic Gregorian ordinal, where January 1 of year 1 has
ordinal 1 (ordinal format of datetime library)

orig

(np.array)

flags indicating historical events (True) or probabilistic (False)

event_name

(list(str))

name of each event (default: event_id)


Note that intensity and fraction are scipy.sparse matrices of size num_events x num_centroids. The fraction attribute is optional. The Centroids class contains the geographical coordinates where the hazard is defined. A Centroids instance provides the coordinates either as points or raster data together with their Coordinate Reference System (CRS). The default CRS used in climada is the usual EPSG:4326. Centroids provides moreover methods to compute centroids areas, on land mask, country iso mask or distance to coast.

How is this tutorial structured?#

Part 1: Read hazards from raster data

Part 2: Read hazards from other data

Part 3: Define hazards manually

Part 4: Analyse hazards

Part 5: Visualize hazards

Part 6: Write (=save) hazards

## Part 1: Read hazards from raster data

Raster data can be read in any format accepted by rasterio using Hazard’s from_raster() method. The raster information might refer to the intensity or fractionof the hazard. Different configuration options such as transforming the coordinates, changing the CRS and reading only a selected area or band are available through the from_raster() arguments as follows:

%matplotlib inline
import numpy as np
from climada.hazard import Hazard
from climada.util.constants import HAZ_DEMO_FL
# to hide the warnings
import warnings
warnings.filterwarnings('ignore')

# read intensity from raster file HAZ_DEMO_FL and set frequency for the contained event
haz_ven = Hazard.from_raster([HAZ_DEMO_FL], attrs={'frequency':np.ones(1)/2}, haz_type='FL')
haz_ven.check()

# The masked values of the raster are set to 0
# Sometimes the raster file does not contain all the information, as in this case the mask value -9999
# We mask it manuall and plot it using plot_intensity()
haz_ven.intensity[haz_ven.intensity==-9999] = 0
haz_ven.plot_intensity(1, smooth=False) # if smooth=True (default value) is used, the computation time might increase

# per default the following attributes have been set
print('event_id: ', haz_ven.event_id)
print('event_name: ', haz_ven.event_name)
print('date: ', haz_ven.date)
print('frequency: ', haz_ven.frequency)
print('orig: ', haz_ven.orig)
print('min, max fraction: ', haz_ven.fraction.min(), haz_ven.fraction.max())
2022-03-16 22:09:56,965 - climada.util.coordinates - INFO - Reading C:\Users\yyljy\climada\demo\data\SC22000_VE__M1.grd.gz
2022-03-16 22:10:01,099 - climada.util.coordinates - INFO - Reading C:\Users\yyljy\climada\demo\data\SC22000_VE__M1.grd.gz
event_id:  [1]
event_name:  ['1']
date:  [1.]
frequency:  [0.5]
orig:  [ True]
min, max fraction:  0.0 1.0
../_images/d23565a8141d024fe33117e9c91b6a06c1ad6cd8ad2e8b7d298387681c4c807e.png

EXERCISE:#

  1. Read raster data in EPSG 2201 Coordinate Reference System (CRS)

  2. Read raster data in its given CRS and transform it to the affine transformation Affine(0.009000000000000341, 0.0, -69.33714959699981, 0.0, -0.009000000000000341, 10.42822096697894), height=500, width=501)

  3. Read raster data in window Window(10, 10, 20, 30)

# Put your code here
# Solution:

# 1. The CRS can be reprojected using dst_crs option
haz = Hazard.from_raster([HAZ_DEMO_FL], dst_crs='epsg:2201', haz_type='FL')
haz.check()
print('\n Solution 1:')
print('centroids CRS:', haz.centroids.crs)
print('raster info:', haz.centroids.meta)

# 2. Transformations of the coordinates can be set using the transform option and Affine
from rasterio import Affine
haz = Hazard.from_raster([HAZ_DEMO_FL], haz_type='FL',
                         transform=Affine(0.009000000000000341, 0.0, -69.33714959699981, \
                                          0.0, -0.009000000000000341, 10.42822096697894),
                         height=500, width=501)
haz.check()
print('\n Solution 2:')
print('raster info:', haz.centroids.meta)
print('intensity size:', haz.intensity.shape)

# 3. A partial part of the raster can be loaded using the window or geometry
from rasterio.windows import Window
haz = Hazard.from_raster([HAZ_DEMO_FL], haz_type='FL', window=Window(10, 10, 20, 30))
haz.check()
print('\n Solution 3:')
print('raster info:', haz.centroids.meta)
print('intensity size:', haz.intensity.shape)
2022-03-16 22:11:00,548 - climada.util.coordinates - INFO - Reading C:\Users\yyljy\climada\demo\data\SC22000_VE__M1.grd.gz
2022-03-16 22:11:04,941 - climada.util.coordinates - INFO - Reading C:\Users\yyljy\climada\demo\data\SC22000_VE__M1.grd.gz

 Solution 1:
centroids CRS: epsg:2201
raster info: {'driver': 'GSBG', 'dtype': 'float32', 'nodata': 1.701410009187828e+38, 'width': 978, 'height': 1091, 'count': 1, 'crs': 'epsg:2201', 'transform': Affine(1011.5372910988809, 0.0, 1120744.5486664253,
       0.0, -1011.5372910988809, 1189133.7652687666)}
2022-03-16 22:11:09,308 - climada.util.coordinates - INFO - Reading C:\Users\yyljy\climada\demo\data\SC22000_VE__M1.grd.gz
2022-03-16 22:11:13,503 - climada.util.coordinates - INFO - Reading C:\Users\yyljy\climada\demo\data\SC22000_VE__M1.grd.gz

 Solution 2:
raster info: {'driver': 'GSBG', 'dtype': 'float32', 'nodata': 1.701410009187828e+38, 'width': 501, 'height': 500, 'count': 1, 'crs': CRS.from_epsg(4326), 'transform': Affine(0.009000000000000341, 0.0, -69.33714959699981,
       0.0, -0.009000000000000341, 10.42822096697894)}
intensity size: (1, 250500)
2022-03-16 22:11:17,739 - climada.util.coordinates - INFO - Reading C:\Users\yyljy\climada\demo\data\SC22000_VE__M1.grd.gz
2022-03-16 22:11:17,780 - climada.util.coordinates - INFO - Reading C:\Users\yyljy\climada\demo\data\SC22000_VE__M1.grd.gz

 Solution 3:
raster info: {'driver': 'GSBG', 'dtype': 'float32', 'nodata': 1.701410009187828e+38, 'width': 20, 'height': 30, 'count': 1, 'crs': CRS.from_epsg(4326), 'transform': Affine(0.009000000000000341, 0.0, -69.2471495969998,
       0.0, -0.009000000000000341, 10.338220966978936)}
intensity size: (1, 600)

## Part 2: Read hazards from other data

  • excel: Hazards can be read from Excel files following the template in climada_python/climada/data/system/hazard_template.xlsx using the from_excel() method.

  • MATLAB: Hazards generated with CLIMADA’s MATLAB version (.mat format) can be read using from_mat().

  • vector data: Use Hazard’s from_vector-constructor to read shape data (all formats supported by fiona).

  • hdf5: Hazards generated with the CLIMADA in Python (.h5 format) can be read using from_hdf5().

from climada.hazard import Hazard, Centroids
from climada.util import HAZ_DEMO_H5 # CLIMADA's Python file
# Hazard needs to know the acronym of the hazard type to be constructed!!! Use 'NA' if not known.
haz_tc_fl = Hazard.from_hdf5(HAZ_DEMO_H5)  # Historic tropical cyclones in Florida from 1990 to 2004
haz_tc_fl.check() # Use always the check() method to see if the hazard has been loaded correctly
2022-03-16 22:11:26,805 - climada.hazard.base - INFO - Reading C:\Users\yyljy\climada\demo\data\tc_fl_1990_2004.h5

## Part 3: Define hazards manually A Hazard can be defined by filling its values one by one, as follows:

# setting points
import numpy as np
from scipy import sparse

lat = np.array([26.933899, 26.957203, 26.783846, 26.645524, 26.897796, 26.925359, \
       26.914768, 26.853491, 26.845099, 26.82651 , 26.842772, 26.825905, \
       26.80465 , 26.788649, 26.704277, 26.71005 , 26.755412, 26.678449, \
       26.725649, 26.720599, 26.71255 , 26.6649  , 26.664699, 26.663149, \
       26.66875 , 26.638517, 26.59309 , 26.617449, 26.620079, 26.596795, \
       26.577049, 26.524585, 26.524158, 26.523737, 26.520284, 26.547349, \
       26.463399, 26.45905 , 26.45558 , 26.453699, 26.449999, 26.397299, \
       26.4084  , 26.40875 , 26.379113, 26.3809  , 26.349068, 26.346349, \
       26.348015, 26.347957])

lon = np.array([-80.128799, -80.098284, -80.748947, -80.550704, -80.596929, \
       -80.220966, -80.07466 , -80.190281, -80.083904, -80.213493, \
       -80.0591  , -80.630096, -80.075301, -80.069885, -80.656841, \
       -80.190085, -80.08955 , -80.041179, -80.1324  , -80.091746, \
       -80.068579, -80.090698, -80.1254  , -80.151401, -80.058749, \
       -80.283371, -80.206901, -80.090649, -80.055001, -80.128711, \
       -80.076435, -80.080105, -80.06398 , -80.178973, -80.110519, \
       -80.057701, -80.064251, -80.07875 , -80.139247, -80.104316, \
       -80.188545, -80.21902 , -80.092391, -80.1575  , -80.102028, \
       -80.16885 , -80.116401, -80.08385 , -80.241305, -80.158855])

n_cen = lon.size # number of centroids
n_ev = 10 # number of events

intensity = sparse.csr_matrix(np.random.random((n_ev, n_cen)))
fraction = intensity.copy()
fraction.data.fill(1)

haz = Hazard(haz_type='TC',
             intensity=intensity,
             fraction=fraction,
             centroids=Centroids.from_lat_lon(lat, lon),  # default crs used
             units='m',
             event_id=np.arange(n_ev, dtype=int),
             event_name=['ev_12', 'ev_21', 'Maria', 'ev_35',
                         'Irma', 'ev_16', 'ev_15', 'Edgar', 'ev_1', 'ev_9'],
             date=np.array([721166, 734447, 734447, 734447, 721167,
                            721166, 721167, 721200, 721166, 721166]),
             orig=np.zeros(n_ev, bool),
             frequency=np.ones(n_ev)/n_ev,)

haz.check()
haz.centroids.plot();
../_images/e4c471cadde480157884fc432041d7df65448b6b694328a2ee3b4fd29eaa773f.png
# setting raster
import numpy as np
from scipy import sparse

# raster info:
# border upper left corner (of the pixel, not of the center of the pixel)
xf_lat = 22
xo_lon = -72
# resolution in lat and lon
d_lat = -0.5 # negative because starting in upper corner
d_lon = 0.5 # same step as d_lat
# number of points
n_lat = 50
n_lon = 40

n_ev = 10 # number of events
centroids = Centroids.from_pix_bounds(xf_lat, xo_lon, d_lat, d_lon, n_lat, n_lon) # default crs used
intensity = sparse.csr_matrix(np.random.random((n_ev, centroids.size)))
fraction = intensity.copy()
fraction.data.fill(1)

haz = Hazard('TC',
             centroids=centroids,
             intensity=intensity,
             fraction=fraction,
             units='m',
             event_id=np.arange(n_ev, dtype=int),
             event_name=['ev_12', 'ev_21', 'Maria', 'ev_35',
                         'Irma', 'ev_16', 'ev_15', 'Edgar', 'ev_1', 'ev_9'],
             date=np.array([721166, 734447, 734447, 734447, 721167,
                            721166, 721167, 721200, 721166, 721166]),
             orig=np.zeros(n_ev, bool),
             frequency=np.ones(n_ev)/n_ev,)

haz.check()
print('Check centroids borders:', haz.centroids.total_bounds)
haz.centroids.plot()

# using from_pnt_bounds, the bounds refer to the bounds of the center of the pixel
left, bottom, right, top = xo_lon, -3.0, -52.0, xf_lat
haz.centroids = Centroids.from_pnt_bounds((left, bottom, right, top), 0.5) # default crs used
print('Check centroids borders:', haz.centroids.total_bounds)
Check centroids borders: (-72.0, -3.0, -52.0, 22.0)
Check centroids borders: (-72.25, -3.25, -51.75, 22.25)
../_images/a5a3f3744f3573f4d742be6e07e8d5c4ed4b106cf1acb7b09d8d8f0279ac74fe.png

Part 4: Analyse Hazards#

The following methods can be used to analyse the data in Hazard:

  • calc_year_set() method returns a dictionary with all the historical (not synthetic) event ids that happened at each year.

  • get_event_date() returns strings of dates in ISO format.

  • To obtain the relation between event ids and event names, two methods can be used get_event_name() and get_event_id().

Other methods to handle one or several Hazards are:

  • the property size returns the number of events contained.

  • append() is used to expand events with data from another Hazard (and same centroids).

  • select() returns a new hazard with the selected region, date and/or synthetic or historical filter.

  • remove_duplicates() removes events with same name and date.

  • local_exceedance_inten() returns a matrix with the exceedence frequency at every frequency and provided return periods. This is the one used in plot_rp_intensity().

  • reproject_raster(), reproject_vector(), raster_to_vector(), vector_to_raster() are methods to change centroids’ CRS and between raster and vector data.

Centroids methods:

  • centroids properties such as area per pixel, distance to coast, country ISO code, on land mask or elevation are available through different set_XX()methods.

  • set_lat_lon_to_meta() computes the raster meta dictionary from present lat and lon. set_meta_to_lat_lon() computes lat and lon of the center of the pixels described in attribute meta. The raster meta information contains at least: width, height, crs and transform data (use help(Centroids) for more info). Using raster centroids can increase computing performance for several computations.

  • when using lats and lons (vector data) the geopandas.GeoSeries geometry attribute contains the CRS information and can be filled with point shapes to perform different computation. The geometry points can be then released using empty_geometry_points().

EXERCISE:#

Using the previous hazard haz_tc_fl answer these questions:

  1. How many synthetic events are contained?

  2. Generate a hazard with historical hurricanes ocurring between 1995 and 2001.

  3. How many historical hurricanes occured in 1999? Which was the year with most hurricanes between 1995 and 2001?

  4. What is the number of centroids with distance to coast smaller than 1km?

# Put your code here:
#help(hist_tc.centroids) # If you want to run it, do it after you execute the next block
# SOLUTION:

# 1.How many synthetic events are contained?
print('Number of total events:', haz_tc_fl.size)
print('Number of synthetic events:', np.logical_not(haz_tc_fl.orig).astype(int).sum())

# 2. Generate a hazard with historical hurricanes ocurring between 1995 and 2001.
hist_tc = haz_tc_fl.select(date=('1995-01-01', '2001-12-31'), orig=True)
print('Number of historical events between 1995 and 2001:', hist_tc.size)

# 3. How many historical hurricanes occured in 1999? Which was the year with most hurricanes between 1995 and 2001?
ev_per_year = hist_tc.calc_year_set() # events ids per year
print('Number of events in 1999:', ev_per_year[1999].size)
max_year = 1995
max_ev = ev_per_year[1995].size
for year, ev in ev_per_year.items():
    if ev.size > max_ev:
        max_year = year
print('Year with most hurricanes between 1995 and 2001:', max_year)

# 4. What is the number of centroids with distance to coast smaller than 1km?
hist_tc.centroids.set_dist_coast()
num_cen_coast = np.argwhere(hist_tc.centroids.dist_coast < 1000).size
print('Number of centroids close to coast: ', num_cen_coast)
Number of total events: 216
Number of synthetic events: 0
Number of historical events between 1995 and 2001: 109
Number of events in 1999: 16
Year with most hurricanes between 1995 and 2001: 1995
2022-03-16 22:12:55,499 - climada.hazard.centroids.centr - INFO - Convert centroids to GeoSeries of Point shapes.
2022-03-16 22:12:57,355 - climada.util.coordinates - INFO - dist_to_coast: UTM 32617 (1/2)
2022-03-16 22:12:59,926 - climada.util.coordinates - INFO - dist_to_coast: UTM 32618 (2/2)
Number of centroids close to coast:  41

Part 5: Visualize Hazards#

There are three different plot functions: plot_intensity(), plot_fraction()and plot_rp_intensity(). Depending on the inputs, different properties can be visualized. Check the documentation of the functions:

help(haz_tc_fl.plot_intensity)
help(haz_tc_fl.plot_rp_intensity)
Help on method plot_intensity in module climada.hazard.base:

plot_intensity(event=None, centr=None, smooth=True, axis=None, adapt_fontsize=True, **kwargs) method of climada.hazard.base.Hazard instance
    Plot intensity values for a selected event or centroid.
    
    Parameters
    ----------
    event: int or str, optional
        If event > 0, plot intensities of
        event with id = event. If event = 0, plot maximum intensity in
        each centroid. If event < 0, plot abs(event)-largest event. If
        event is string, plot events with that name.
    centr: int or tuple, optional
        If centr > 0, plot intensity
        of all events at centroid with id = centr. If centr = 0,
        plot maximum intensity of each event. If centr < 0,
        plot abs(centr)-largest centroid where higher intensities
        are reached. If tuple with (lat, lon) plot intensity of nearest
        centroid.
    smooth: bool, optional
        Rescale data to RESOLUTIONxRESOLUTION pixels (see constant
        in module `climada.util.plot`)
    axis: matplotlib.axes._subplots.AxesSubplot, optional
        axis to use
    kwargs: optional
        arguments for pcolormesh matplotlib function
        used in event plots or for plot function used in centroids plots
    
    Returns
    -------
        matplotlib.axes._subplots.AxesSubplot
    
    Raises
    ------
        ValueError

Help on method plot_rp_intensity in module climada.hazard.base:

plot_rp_intensity(return_periods=(25, 50, 100, 250), smooth=True, axis=None, figsize=(9, 13), adapt_fontsize=True, **kwargs) method of climada.hazard.base.Hazard instance
    Compute and plot hazard exceedance intensity maps for different
    return periods. Calls local_exceedance_inten.
    
    Parameters
    ----------
    return_periods: tuple(int), optional
        return periods to consider
    smooth: bool, optional
        smooth plot to plot.RESOLUTIONxplot.RESOLUTION
    axis: matplotlib.axes._subplots.AxesSubplot, optional
        axis to use
    figsize: tuple, optional
        figure size for plt.subplots
    kwargs: optional
        arguments for pcolormesh matplotlib function used in event plots
    
    Returns
    -------
    axis, inten_stats:  matplotlib.axes._subplots.AxesSubplot, np.ndarray
        intenstats is return_periods.size x num_centroids
# 1. intensities of the largest event (defined as greater sum of intensities):
# all events:
haz_tc_fl.plot_intensity(event=-1) # largest historical event: 1992230N11325 hurricane ANDREW

# 2. maximum intensities at each centroid:
haz_tc_fl.plot_intensity(event=0)

# 3. intensities of hurricane 1998295N12284:
haz_tc_fl.plot_intensity(event='1998295N12284', cmap='BuGn') # setting color map

# 4. tropical cyclone intensities maps for the return periods [10, 50, 75, 100]
_, res = haz_tc_fl.plot_rp_intensity([10, 50, 75, 100])

# 5. intensities of all the events in centroid with id 50
haz_tc_fl.plot_intensity(centr=50)

# 6. intensities of all the events in centroid closest to lat, lon = (26.5, -81)
haz_tc_fl.plot_intensity(centr=(26.5, -81));
2022-03-16 22:13:50,822 - climada.hazard.base - INFO - Computing exceedance intenstiy map for return periods: [ 10  50  75 100]
2022-03-16 22:13:51,449 - climada.hazard.base - WARNING - Exceedance intenstiy values below 0 are set to 0.                    Reason: no negative intensity values were found in hazard.
../_images/0345488c2aa1f38ecddee8ebd3d0baac2da28f269de8b6b810eb1f4bd61419f5.png ../_images/c78fea5e6acda7a975cb588c7b0b2546bc980a1ba48616b651f39a85a16515de.png ../_images/60c0f4de49b2c7b1b99a5df8a11d7a5c4c8496090686376bc334de0ef1fab912.png ../_images/16b107f0f7b81013ac54e8b1d384d70eaf7e4c311fc2be4d2a43f74f3af93e6e.png ../_images/a658f83bb961f22ad9cc52f5f61423ce59d8a44df49ad2d1285a1a99e207b5e4.png ../_images/ccdc1257f3745c7086d3cf46357976bbdcac71d3069c80e698b315088a364c38.png
# 7. one figure with two plots: maximum intensities and selected centroid with all intensities:
from climada.util.plot import make_map
import matplotlib.pyplot as plt
plt.ioff()
fig, ax1, fontsize = make_map(1)  # map
ax2 = fig.add_subplot(2, 1, 2) # add regular axes
haz_tc_fl.plot_intensity(axis=ax1, event=0) # plot original resolution
ax1.plot(-80, 26, 'or', mfc='none', markersize=12)
haz_tc_fl.plot_intensity(axis=ax2, centr=(26, -80))
fig.subplots_adjust(hspace=6.5)
../_images/19fec1bd4b156d7d9683b590bee672807d501411d218bc9d739c35ee20457d7d.png

Part 6: Write (=save) hazards#

Hazards can be written and read in hdf5 format as follows:

haz_tc_fl.write_hdf5('results/haz_tc_fl.h5')

haz = Hazard.from_hdf5('results/haz_tc_fl.h5')
haz.check()
2022-03-16 22:17:30,344 - climada.hazard.base - INFO - Writing results/haz_tc_fl.h5
2022-03-16 22:17:30,572 - climada.hazard.base - INFO - Reading results/haz_tc_fl.h5

GeoTiff data is generated using write_raster():

haz_ven.write_raster('results/haz_ven.tif') # each event is a band of the tif file
2022-03-16 22:17:41,283 - climada.util.coordinates - INFO - Writting results/haz_ven.tif

Pickle will work as well:

from climada.util.save import save
# this generates a results folder in the current path and stores the output there
save('tutorial_haz_tc_fl.p', haz_tc_fl)
2022-03-16 22:17:55,836 - climada.util.save - INFO - Written file C:\Users\yyljy\Documents\climada_main\doc\tutorial\results\tutorial_haz_tc_fl.p