Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP
Fetching contributors…

Cannot retrieve contributors at this time

170 lines (140 sloc) 7.219 kB
"""Impact function based on Padang 2009 post earthquake survey
This impact function estimates percentual damage to buildings as a
function of ground shaking measured in MMI.
Buildings are assumed to fall the 9 classes below as described in
the Geoscience Australia/ITB 2009 Padang earthquake
survey (http://trove.nla.gov.au/work/38470066).
Class Building Type Median (MMI) Beta (MMI)
-------------------------------------------------------------------------
1 URM with river rock walls 7.5 0.11
2 URM with Metal Roof 8.3 0.1
3 Timber frame with masonry in-fill 8.8 0.11
4 RC medium rise Frame with Masonry in-fill walls 8.4 0.05
5 Timber frame with stucco in-fill 9.2 0.11
6 Concrete Shear wall high rise* Hazus C2H 9.7 0.15
7 RC low rise Frame with Masonry in-fill walls 9 0.08
8 Confined Masonry 8.9 0.07
9 Timber frame residential 10.5 0.15
"""
from safe.impact_functions.core import (FunctionProvider,
get_hazard_layer,
get_exposure_layer,
get_question,
format_int)
from safe.storage.vector import Vector
from safe.common.utilities import ugettext as tr
from safe.common.numerics import lognormal_cdf
from safe.common.tables import Table, TableRow
from safe.impact_functions.mappings import osm2padang, sigab2padang
from safe.engine.interpolation import assign_hazard_values_to_exposure_data
# Damage curves for each of the nine classes derived from the Padang survey
damage_curves = {1: dict(median=7.5, beta=0.11),
2: dict(median=8.3, beta=0.1),
3: dict(median=8.8, beta=0.11),
4: dict(median=8.4, beta=0.05),
5: dict(median=9.2, beta=0.11),
6: dict(median=9.7, beta=0.15),
7: dict(median=9.0, beta=0.08),
8: dict(median=8.9, beta=0.07),
9: dict(median=10.5, beta=0.15)}
class PadangEarthquakeBuildingDamageFunction(FunctionProvider):
"""Risk plugin for Padang earthquake damage to buildings
:param requires category=='hazard' and \
subcategory=='earthquake' and \
layertype=='raster' and \
unit=='MMI'
:param requires category=='exposure' and \
subcategory=='structure' and \
layertype=='vector' and \
datatype in ['osm', 'itb', 'sigab']
"""
title = tr('Be damaged depending on building type')
def run(self, layers):
"""Risk plugin for Padang building survey
"""
# Extract data
H = get_hazard_layer(layers) # Ground shaking
E = get_exposure_layer(layers) # Building locations
question = get_question(H.get_name(),
E.get_name(),
self)
# Map from different kinds of datasets to Padang vulnerability classes
datatype = E.get_keywords()['datatype']
vclass_tag = 'VCLASS'
if datatype.lower() == 'osm':
# Map from OSM attributes
Emap = osm2padang(E)
elif datatype.lower() == 'sigab':
# Map from SIGAB attributes
Emap = sigab2padang(E)
else:
Emap = E
# Interpolate hazard level to building locations
I = assign_hazard_values_to_exposure_data(H, Emap,
attribute_name='MMI')
# Extract relevant numerical data
attributes = I.get_data()
N = len(I)
# Calculate building damage
count_high = count_medium = count_low = count_none = 0
for i in range(N):
mmi = float(attributes[i]['MMI'])
building_type = Emap.get_data(vclass_tag, i)
damage_params = damage_curves[building_type]
beta = damage_params['beta']
median = damage_params['median']
percent_damage = lognormal_cdf(mmi,
median=median,
sigma=beta) * 100
# Add calculated impact to existing attributes
attributes[i][self.target_field] = percent_damage
# Calculate statistics
if percent_damage < 10:
count_none += 1
if 10 <= percent_damage < 33:
count_low += 1
if 33 <= percent_damage < 66:
count_medium += 1
if 66 <= percent_damage:
count_high += 1
# Generate impact report
table_body = [question,
TableRow([tr('Buildings'), tr('Total')],
header=True),
TableRow([tr('All'), N]),
TableRow([tr('No damage'), format_int(count_none)]),
TableRow([tr('Low damage'), format_int(count_low)]),
TableRow([tr('Medium damage'),
format_int(count_medium)]),
TableRow([tr('High damage'), format_int(count_high)])]
table_body.append(TableRow(tr('Notes'), header=True))
table_body.append(tr('Levels of impact are defined by post 2009 '
'Padang earthquake survey conducted by Geoscience '
'Australia and Institute of Teknologi Bandung.'))
table_body.append(tr('Unreinforced masonry is assumed where no '
'structural information is available.'))
impact_summary = Table(table_body).toNewlineFreeString()
impact_table = impact_summary
map_title = tr('Earthquake damage to buildings')
# Create style
style_classes = [dict(label=tr('No damage'), min=0, max=10,
colour='#00ff00', transparency=0),
dict(label=tr('Low damage'), min=10, max=33,
colour='#ffff00', transparency=0),
dict(label=tr('Medium damage'), min=33, max=66,
colour='#ffaa00', transparency=0),
dict(label=tr('High damage'), min=66, max=100,
colour='#ff0000', transparency=0)]
style_info = dict(target_field=self.target_field,
style_classes=style_classes)
# Create vector layer and return
V = Vector(data=attributes,
projection=E.get_projection(),
geometry=E.get_geometry(),
name='Estimated pct damage',
keywords={'impact_summary': impact_summary,
'impact_table': impact_table,
'map_title': map_title,
'target_field': self.target_field},
style_info=style_info)
return V
Jump to Line
Something went wrong with that request. Please try again.