Source code for openquake.calculators.scenario_damage

# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4
#
# Copyright (C) 2014-2016 GEM Foundation
#
# OpenQuake is free software: you can redistribute it and/or modify it
# under the terms of the GNU Affero General Public License as published
# by the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.
#
# OpenQuake is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with OpenQuake. If not, see <http://www.gnu.org/licenses/>.

import itertools
import numpy

from openquake.commonlib import parallel, riskmodels, calc
from openquake.risklib import scientific
from openquake.calculators import base

F32 = numpy.float32
F64 = numpy.float64


[docs]def dist_by_asset(data, multi_stat_dt): """ :param data: array of shape (N, L, R, 2, ...) :param multi_stat_dt: numpy dtype for statistical outputs :returns: array of shape (N, R) with records of type multi_stat_dt """ N, L, R = data.shape[:3] out = numpy.zeros((N, R), multi_stat_dt) for l, lt in enumerate(multi_stat_dt.names): out_lt = out[lt] for n, r in itertools.product(range(N), range(R)): out_lt[n, r] = data[n, l, r] return out
[docs]def dist_by_taxon(data, multi_stat_dt): """ :param data: array of shape (T, L, R, ...) :param multi_stat_dt: numpy dtype for statistical outputs :returns: array of shape (T, R) with records of type multi_stat_dt """ T, L, R = data.shape[:3] out = numpy.zeros((T, R), multi_stat_dt) for l, lt in enumerate(multi_stat_dt.names): out_lt = out[lt] for t, r in itertools.product(range(T), range(R)): out_lt[t, r] = scientific.mean_std(data[t, l, r]) return out
[docs]def dist_total(data, multi_stat_dt): """ :param data: array of shape (T, L, R, ...) :param multi_stat_dt: numpy dtype for statistical outputs :returns: array of shape (R,) with records of type multi_stat_dt """ T, L, R = data.shape[:3] total = data.sum(axis=0) out = numpy.zeros(R, multi_stat_dt) for l, lt in enumerate(multi_stat_dt.names): out_lt = out[lt] for r in range(R): out_lt[r] = scientific.mean_std(total[l, r]) return out
@parallel.litetask
[docs]def scenario_damage(riskinput, riskmodel, rlzs_assoc, monitor): """ Core function for a damage computation. :param riskinput: a :class:`openquake.risklib.riskinput.RiskInput` object :param riskmodel: a :class:`openquake.risklib.riskinput.CompositeRiskModel` instance :param rlzs_assoc: a class:`openquake.commonlib.source.RlzsAssoc` instance :param monitor: :class:`openquake.baselib.performance.Monitor` instance :returns: a dictionary {'d_asset': [(l, r, a, mean-stddev), ...], 'd_taxonomy': damage array of shape T, L, R, E, D, 'c_asset': [(l, r, a, mean-stddev), ...], 'c_taxonomy': damage array of shape T, L, R, E} `d_asset` and `d_taxonomy` are related to the damage distributions whereas `c_asset` and `c_taxonomy` are the consequence distributions. If there is no consequence model `c_asset` is an empty list and `c_taxonomy` is a zero-value array. """ c_models = monitor.consequence_models L = len(riskmodel.loss_types) R = len(rlzs_assoc.realizations) D = len(riskmodel.damage_states) E = monitor.oqparam.number_of_ground_motion_fields T = len(monitor.taxonomies) taxo2idx = {taxo: i for i, taxo in enumerate(monitor.taxonomies)} result = dict(d_asset=[], d_taxon=numpy.zeros((T, L, R, E, D), F64), c_asset=[], c_taxon=numpy.zeros((T, L, R, E), F64)) for out_by_lr in riskmodel.gen_outputs( riskinput, rlzs_assoc, monitor): for (l, r), out in sorted(out_by_lr.items()): c_model = c_models.get(out.loss_type) for asset, fraction in zip(out.assets, out.damages): t = taxo2idx[asset.taxonomy] damages = fraction * asset.number if c_model: # compute consequences means = [par[0] for par in c_model[asset.taxonomy].params] # NB: we add a 0 in front for nodamage state c_ratio = numpy.dot(fraction, [0] + means) consequences = c_ratio * asset.value(out.loss_type) result['c_asset'].append( (l, r, asset.ordinal, scientific.mean_std(consequences))) result['c_taxon'][t, l, r, :] += consequences # TODO: consequences for the occupants result['d_asset'].append( (l, r, asset.ordinal, scientific.mean_std(damages))) result['d_taxon'][t, l, r, :] += damages return result
@base.calculators.add('scenario_damage')
[docs]class ScenarioDamageCalculator(base.RiskCalculator): """ Scenario damage calculator """ pre_calculator = 'scenario' core_task = scenario_damage is_stochastic = True
[docs] def pre_execute(self): if 'gmfs' in self.oqparam.inputs: self.pre_calculator = None base.RiskCalculator.pre_execute(self) self.monitor.consequence_models = riskmodels.get_risk_models( self.oqparam, 'consequence') self.etags, gmfs = calc.get_gmfs(self.datastore) self.riskinputs = self.build_riskinputs(gmfs) self.monitor.taxonomies = sorted(self.taxonomies)
[docs] def post_execute(self, result): """ Compute stats for the aggregated distributions and save the results on the datastore. """ dstates = self.riskmodel.damage_states ltypes = self.riskmodel.loss_types L = len(ltypes) R = len(self.rlzs_assoc.realizations) D = len(dstates) N = len(self.assetcol) # damage distributions dt_list = [] for ltype in ltypes: dt_list.append((ltype, numpy.dtype([('mean', (F32, D)), ('stddev', (F32, D))]))) multi_stat_dt = numpy.dtype(dt_list) d_asset = numpy.zeros((N, L, R, 2, D), F32) for (l, r, a, stat) in result['d_asset']: d_asset[a, l, r] = stat self.datastore['dmg_by_asset'] = dist_by_asset( d_asset, multi_stat_dt) self.datastore['dmg_by_taxon'] = dist_by_taxon( result['d_taxon'], multi_stat_dt) self.datastore['dmg_total'] = dist_total( result['d_taxon'], multi_stat_dt) # consequence distributions if result['c_asset']: c_asset = numpy.zeros((N, L, R, 2), F32) for (l, r, a, stat) in result['c_asset']: c_asset[a, l, r] = stat multi_stat_dt = numpy.dtype( [(lt, [('mean', F32), ('stddev', F32)]) for lt in ltypes]) self.datastore['csq_by_asset'] = dist_by_asset( c_asset, multi_stat_dt) self.datastore['csq_by_taxon'] = dist_by_taxon( result['c_taxon'], multi_stat_dt) self.datastore['csq_total'] = dist_total( result['c_taxon'], multi_stat_dt)