# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4
# Copyright (c) 2010-2014, 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 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/>.
# Disable:
# - 'Maximum number of public methods for a class'
# - 'Missing docstring' (because of all of the model Meta)
# pylint: disable=R0904,C0111
'''
Model representations of the OpenQuake DB tables.
'''
import os
import re
import collections
import operator
from datetime import datetime
import numpy
from scipy import interpolate
from django.db import connections, transaction
from django.core.exceptions import ObjectDoesNotExist
from django.contrib.gis.db import models as djm
from openquake.hazardlib.imt import from_string
from openquake.hazardlib import geo, correlation
from openquake.commonlib.riskmodels import loss_type_to_cost_type
from openquake.commonlib.readinput import get_mesh
from openquake.commonlib.oqvalidation import OqParam
from openquake.commonlib import logictree
from openquake.engine.db import fields
from openquake.engine import writer, logs, utils
#: Kind of supported curve statistics
STAT_CHOICES = (
(u'mean', u'Mean'),
(u'quantile', u'Quantile'))
#: System Reference ID used for geometry objects
DEFAULT_SRID = 4326
VS30_TYPE_CHOICES = (
(u"measured", u"Value obtained from on-site measurements"),
(u"inferred", u"Estimated value"),
)
IMT_CHOICES = (
(u'PGA', u'Peak Ground Acceleration'),
(u'PGV', u'Peak Ground Velocity'),
(u'PGD', u'Peak Ground Displacement'),
(u'SA', u'Spectral Acceleration'),
(u'IA', u'Arias Intensity'),
(u'RSD', u'Relative Significant Duration'),
(u'MMI', u'Modified Mercalli Intensity'),
)
#: Default Loss Curve Resolution used for probabilistic risk calculators
DEFAULT_LOSS_CURVE_RESOLUTION = 50
#: Minimum value for a seed number
MIN_SINT_32 = -(2 ** 31)
#: Maximum value for a seed number
MAX_SINT_32 = (2 ** 31) - 1
#: Kind of supported type of loss outputs
LOSS_TYPES = ["structural", "nonstructural", "fatalities", "contents"]
#: relative tolerance to consider two risk outputs (almost) equal
RISK_RTOL = 0.05
#: absolute tolerance to consider two risk outputs (almost) equal
RISK_ATOL = 0.01
# TODO: these want to be dictionaries
INPUT_TYPE_CHOICES = (
(u'unknown', u'Unknown'),
(u'source', u'Source Model'),
(u'source_model_logic_tree', u'Source Model Logic Tree'),
(u'gsim_logic_tree', u'Ground Shaking Intensity Model Logic Tree'),
(u'exposure', u'Exposure'),
(u'fragility', u'Fragility'),
(u'site_model', u'Site Model'),
(u'rupture_model', u'Rupture Model'),
# vulnerability models
(u'structural_vulnerability', u'Structural Vulnerability'),
(u'nonstructural_vulnerability', u'Non Structural Vulnerability'),
(u'contents_vulnerability', u'Contents Vulnerability'),
(u'business_interruption_vulnerability',
u'Business Interruption Vulnerability'),
(u'occupants_vulnerability', u'Occupants Vulnerability'),
(u'structural_vulnerability_retrofitted',
u'Structural Vulnerability Retrofitted'))
RAISE_EXC = object() # sentinel used in OqJob.get_param
[docs]class MissingParameter(KeyError):
"""Raised by OqJob.get_param when a parameter is missing in the database"""
############## Fix FloatField underflow error ##################
# http://stackoverflow.com/questions/9556586/floating-point-numbers-of-python-float-and-postgresql-double-precision
def _get_prep_value(self, value):
if value is None:
return None
val = float(value)
if abs(val) < 1E-300:
return 0.
return val
djm.FloatField.get_prep_value = _get_prep_value
[docs]def risk_almost_equal(o1, o2, key=lambda x: x, rtol=RISK_RTOL, atol=RISK_ATOL):
return numpy.testing.assert_allclose(
numpy.array(key(o1)), numpy.array(key(o2)), rtol=rtol, atol=atol)
[docs]def loss_curve_almost_equal(curve, expected_curve):
if getattr(curve, 'asset_value', None) == 0.0 and getattr(
expected_curve, 'asset_value', None) == 0.0:
return risk_almost_equal(curve.loss_ratios, expected_curve.loss_ratios)
elif curve.losses[curve.losses > 0].any():
poes = interpolate.interp1d(
curve.losses, curve.poes,
bounds_error=False, fill_value=0)(expected_curve.losses)
else:
poes = numpy.zeros(len(expected_curve.poes))
return risk_almost_equal(poes, expected_curve.poes)
[docs]def getcursor(route):
"""Return a cursor from a Django route"""
return connections[route].cursor()
[docs]def order_by_location(queryset):
"""
Utility function to order a queryset by location. This works even if
the location is of Geography object (a simple order_by('location') only
works for Geometry objects).
"""
return queryset.extra(
select={'x': 'ST_X(geometry(location))',
'y': 'ST_Y(geometry(location))'},
order_by=["x", "y"])
[docs]def queryset_iter(queryset, chunk_size):
"""
Given a QuerySet, split it into smaller queries and yield the result of
each.
:param queryset:
A :class:`django.db.models.query.QuerySet` to iterate over, in chunks
of ``chunk_size``.
:param int chunksize:
Chunk size for iteration over query results. For an unexecuted
QuerySet, this will result in splitting a (potentially large) query
into smaller queries.
"""
offset = 0
while True:
chunk = list(queryset[offset:offset + chunk_size].iterator())
if len(chunk) == 0:
raise StopIteration
else:
yield chunk
offset += chunk_size
[docs]def build_curves(rlz, curves_by_trt_model_gsim):
"""
Build on the fly the hazard curves for the current realization
by looking at the associations stored in the database table
`hzrdr.assoc_lt_rlz_trt_model`.
"""
# fixed a realization, there are T associations where T is the
# number of TrtModels
curves = 0
for art in AssocLtRlzTrtModel.objects.filter(rlz=rlz):
pnes = 1. - curves_by_trt_model_gsim[art.trt_model_id, art.gsim]
curves = 1. - (1. - curves) * pnes
return curves
## Tables in the 'hzrdi' (Hazard Input) schema.
[docs]class SiteModel(djm.Model):
'''
A model for site-specific parameters, used in hazard calculations.
'''
job = djm.ForeignKey('OqJob')
# Average shear wave velocity for top 30 m. Units m/s.
vs30 = djm.FloatField()
# 'measured' or 'inferred'
vs30_type = djm.TextField(choices=VS30_TYPE_CHOICES)
# Depth to shear wave velocity of 1.0 km/s. Units m.
z1pt0 = djm.FloatField()
# Depth to shear wave velocity of 2.5 km/s. Units km.
z2pt5 = djm.FloatField()
location = djm.PointField(srid=DEFAULT_SRID)
@property
def measured(self):
"""True or False depending on the field vs30_type"""
return self.vs30_type == 'measured'
def __repr__(self):
return (
'SiteModel(location="%s", vs30=%s, vs30_type=%s, z1pt0=%s, '
'z2pt5=%s)'
% (self.location.wkt, self.vs30, self.vs30_type, self.z1pt0,
self.z2pt5))
class Meta:
db_table = 'hzrdi\".\"site_model'
## Tables in the 'uiapi' schema.
[docs]class OqJob(djm.Model):
'''
An OpenQuake engine run started by the user
'''
user_name = djm.TextField()
hazard_calculation = djm.ForeignKey('OqJob', null=True)
LOG_LEVEL_CHOICES = (
(u'debug', u'Debug'),
(u'info', u'Info'),
(u'progress', u'Progress'),
(u'warn', u'Warn'),
(u'error', u'Error'),
(u'critical', u'Critical'),
)
log_level = djm.TextField(choices=LOG_LEVEL_CHOICES, default='progress')
STATUS_CHOICES = (
(u'pre_executing', u'Pre-Executing'),
(u'executing', u'Executing'),
(u'post_executing', u'Post-Executing'),
(u'post_processing', u'Post-Processing'),
(u'export', u'Exporting results'),
(u'clean_up', u'Cleaning up'),
(u'complete', u'Complete'),
)
status = djm.TextField(choices=STATUS_CHOICES, default='pre_executing')
oq_version = djm.TextField(null=True, blank=True)
hazardlib_version = djm.TextField(null=True, blank=True)
commonlib_version = djm.TextField(null=True, blank=True)
risklib_version = djm.TextField(null=True, blank=True)
is_running = djm.BooleanField(default=False)
duration = djm.IntegerField(default=0)
job_pid = djm.IntegerField(default=0)
supervisor_pid = djm.IntegerField(default=0)
last_update = djm.DateTimeField(editable=False, default=datetime.utcnow)
class Meta:
db_table = 'uiapi\".\"oq_job'
@property
def risk_calculation(self):
return RiskCalculation(self)
@property
def job_type(self):
"""
'hazard' or 'risk'
"""
# only if the job is of kind 'risk' the field hazard_calculation_id
# is not null and contains a reference to the previous hazard job
return 'hazard' if self.hazard_calculation is None else 'risk'
[docs] def get_param(self, name, missing=RAISE_EXC):
"""
`job.get_param(name)` returns the value of the requested parameter
or raise a MissingParameter exception if the parameter does not
exist in the database.
`job.get_param(name, missing)` returns the value of the requested
parameter or the `missing` value if the parameter does not
exist in the database.
:param name: the name of the parameter
:param missing: value returned if the parameter is missing
NB: since job_param.value is NOT NULL, `.get_param(name)`
can return None only if the parameter is missing.
"""
try:
return JobParam.objects.get(job=self, name=name).value
except ObjectDoesNotExist:
if missing is RAISE_EXC:
raise MissingParameter(name)
return missing
[docs] def get_oqparam(self):
"""
Return an OqParam object as read from the database
"""
oqparam = object.__new__(OqParam)
for row in JobParam.objects.filter(job=self):
setattr(oqparam, row.name, row.value)
return oqparam
[docs] def save_params(self, params):
"""
Save on the database table job_params the given parameters.
:param job: an :class:`OqJob` instance
:param params: a dictionary {name: string} of parameters
"""
for name, value in params.iteritems():
JobParam.objects.create(job=self, name=name, value=repr(value))
[docs] def save_hazard_sites(self):
"""
Populate the table HazardSite by inferring the points from
the sites, region, or exposure.
"""
assert self.job_type == 'hazard', self.job_type
oqparam = self.get_oqparam()
if 'exposure' in oqparam.inputs:
assets = self.exposuremodel.exposuredata_set.all()
# the coords here must be sorted; the issue is that the
# disaggregation calculator has a for loop of kind
# for site in sites:
# bin_edge, disagg_matrix = disaggregation(site, ...)
# the generated ruptures are random if the order of the sites
# is random, even if the seed is fixed; in particular for some
# ordering no ruptures are generated and the test
# qa_tests/hazard/disagg/case_1/test.py fails with a bad
# error message
coords = sorted(
set((asset.site.x, asset.site.y) for asset in assets))
lons, lats = zip(*coords)
mesh = geo.Mesh(numpy.array(lons), numpy.array(lats), None)
else:
mesh = get_mesh(oqparam)
sids = save_sites(self, ((p.longitude, p.latitude) for p in mesh))
return mesh, sids
def __repr__(self):
return '<%s %d, %s>' % (self.__class__.__name__,
self.id, self.job_type)
[docs]def oqparam(job_id):
"""
:param job_id: ID of :class:`openquake.engine.db.models.OqJob`
:returns: instance of :class:`openquake.commonlib.oqvalidation.OqParam`
"""
return OqJob.objects.get(pk=job_id).get_oqparam()
[docs]class JobStats(djm.Model):
'''
Capture various statistics about a job.
'''
oq_job = djm.OneToOneField('OqJob')
start_time = djm.DateTimeField(editable=False, default=datetime.utcnow)
stop_time = djm.DateTimeField(editable=False)
# The disk space occupation in bytes
disk_space = djm.IntegerField(null=True)
class Meta:
db_table = 'uiapi\".\"job_stats'
# created at the end of the pre_execute phase
[docs]class JobInfo(djm.Model):
'''
Store information about a job.
'''
oq_job = djm.OneToOneField('OqJob')
parent_job = djm.ForeignKey('OqJob')
num_sites = djm.IntegerField(null=False)
num_realizations = djm.IntegerField(null=False)
num_imts = djm.IntegerField(null=False)
num_levels = djm.FloatField(null=False)
input_weight = djm.FloatField(null=False)
output_weight = djm.FloatField(null=False)
class Meta:
db_table = 'uiapi\".\"job_info'
[docs]class JobParam(djm.Model):
'''
The parameters of a job
'''
job = djm.ForeignKey('OqJob')
name = djm.TextField(null=False)
value = fields.LiteralField(null=False)
class Meta:
db_table = 'uiapi\".\"job_param'
[docs]def save_sites(job, coords):
"""
Save all the gives sites on the hzrdi.hazard_site table.
:param coordinates: a sequence of (lon, lat) pairs
:returns: the ids of the inserted HazardSite instances
"""
sites = [HazardSite(hazard_calculation=job,
location='POINT(%s %s)' % (lon, lat))
for lon, lat in coords]
return writer.CacheInserter.saveall(sites)
[docs]class RiskCalculation(object):
#: Default maximum asset-hazard distance in km
DEFAULT_MAXIMUM_DISTANCE = 5
def __init__(self, job):
self.oqjob = job
vars(self).update(vars(job.get_oqparam()))
self.hazard_output = Output.objects.get(pk=self.hazard_output_id) \
if self.hazard_output_id else None
self.hazard_calculation = OqJob.objects.get(
pk=self.hazard_calculation_id)
if not hasattr(self, 'taxonomies_from_model'):
self.taxonomies_from_model = None
if not hasattr(self, 'time_event'):
self.time_event = None
if not hasattr(self, 'asset_correlation'):
self.asset_correlation = 0
if not hasattr(self, 'master_seed'):
self.master_seed = 42
if not hasattr(self, 'insured_losses'):
self.insured_losses = False
if not hasattr(self, 'risk_investigation_time'):
self.risk_investigation_time = None
if not hasattr(self, 'loss_curve_resolution'):
self.loss_curve_resolution = 50 # default
[docs] def get_hazard_param(self):
"""
Get the hazard parameters associated with the hazard job that generated
the output used as an input for the current risk calculation.
:returns:
:class:`openquake.commonlib.oqvalidation.OqParam` instance.
"""
return self.hazard_calculation.get_oqparam()
[docs] def hazard_outputs(self):
"""
Returns the list of hazard outputs to be considered. Apply
`filters` to the default queryset
"""
if self.hazard_output:
return [self.hazard_output]
elif self.hazard_calculation:
if self.calculation_mode in ["classical_risk", "classical_bcr"]:
filters = dict(output_type='hazard_curve_multi',
hazard_curve__lt_realization__isnull=False)
elif self.calculation_mode in [
"event_based_risk", "event_based_bcr"]:
filters = dict(
output_type='gmf', gmf__lt_realization__isnull=False)
elif self.calculation_mode == "event_based_fr":
filters = dict(output_type='ses')
elif self.calculation_mode in ['scenario_risk', 'scenario_damage']:
filters = dict(output_type='gmf_scenario')
else:
raise NotImplementedError
return self.hazard_calculation.output_set.filter(
**filters).order_by('id')
else:
raise RuntimeError("Neither hazard calculation "
"neither a hazard output has been provided")
@property
def best_maximum_distance(self):
"""
Get the asset-hazard maximum distance (in km)
:returns:
The minimum between the maximum distance provided by the user (if
not given, `DEFAULT_MAXIMUM_DISTANCE` is used as default) and the
step (if exists) used by the hazard calculation.
"""
dist = getattr(self, 'maximum_distance', self.DEFAULT_MAXIMUM_DISTANCE)
grid_spacing = self.hazard_calculation.get_param(
'region_grid_spacing', None)
if grid_spacing:
dist = min(dist, grid_spacing * numpy.sqrt(2) / 2)
return dist
@property
def preloaded_exposure_model(self):
pem_id = getattr(self, 'preloaded_exposure_model_id', None)
return ExposureModel.objects.get(pk=pem_id) if pem_id else None
@property
def exposure_model(self):
"""
Return the right exposure model by following rules in order:
1. if the `preloaded_exposure_model_id` is set in job_risk.ini, use it
2. if an exposure_file is defined in job_risk.ini, use it
3. if an exposure was used in the hazard job, use it
4. if no exposure is found, return None
"""
return (self.preloaded_exposure_model or
extract_from(
[self.oqjob, self.hazard_calculation], 'exposuremodel'))
@property
def investigation_time(self):
return (self.risk_investigation_time or
self.get_hazard_param().investigation_time)
[docs]class Imt(djm.Model):
"""
Table with the Intensity Measure Types
"""
imt_str = djm.TextField(null=False)
im_type = djm.TextField(choices=IMT_CHOICES)
sa_period = djm.FloatField(null=True)
sa_damping = djm.FloatField(null=True)
stored_imts = None
@classmethod
[docs] def get(cls, imt_str):
"""
:param imt_str: a string specifying the IMT
:returns: a :class:`openquake.engine.db.models.Imt` instance
"""
if cls.stored_imts is None: # the first time
cls.stored_imts = {imt.imt_str: imt
for imt in Imt.objects.filter()}
return cls.stored_imts[imt_str]
@classmethod
[docs] def save_new(cls, hazardlib_imts):
"""
Save the intensity measure types not already stored in the database.
:param hazardlib_imts: a list of hazardlib IMT tuples
"""
if cls.stored_imts is None: # the first time
cls.stored_imts = {imt.imt_str: imt
for imt in Imt.objects.filter()}
for x in hazardlib_imts:
imt_str = str(x)
if imt_str not in cls.stored_imts:
imt = cls.objects.create(
imt_str=imt_str,
im_type=x[0], sa_period=x[1], sa_damping=x[2])
cls.stored_imts[imt_str] = imt
class Meta:
db_table = 'hzrdi\".\"imt'
[docs]class ImtTaxonomy(djm.Model):
"""
Table with the associations IMT, taxonomy, as extracted from the risk
models.
"""
job = djm.ForeignKey('OqJob', null=False)
imt = djm.ForeignKey('Imt', null=False)
taxonomy = djm.TextField(null=True)
class Meta:
db_table = 'riski\".\"imt_taxonomy'
[docs]class OutputManager(djm.Manager):
"""
Manager class to filter and create Output objects
"""
[docs] def create_output(self, job, display_name, output_type):
"""
Create an output for the given `job`, `display_name` and
`output_type` (default to hazard_curve)
"""
return self.create(oq_job=job,
display_name=display_name,
output_type=output_type)
[docs]class Output(djm.Model):
'''
A single artifact which is a result of an OpenQuake job.
The data may reside in a file or in the database.
'''
#: Metadata of hazard outputs used by risk calculation. See
#: `hazard_metadata` property for more details
HazardMetadata = collections.namedtuple(
'hazard_metadata',
'investigation_time statistics quantile sm_path gsim_path')
#: Hold the full paths in the model trees of ground shaking
#: intensity models and of source models, respectively.
LogicTreePath = collections.namedtuple(
'logic_tree_path',
'gsim_path sm_path')
#: Hold the statistical params (statistics, quantile).
StatisticalParams = collections.namedtuple(
'statistical_params',
'statistics quantile')
oq_job = djm.ForeignKey('OqJob', null=False)
display_name = djm.TextField()
HAZARD_OUTPUT_TYPE_CHOICES = (
(u'disagg_matrix', u'Disaggregation Matrix'),
(u'gmf', u'Ground Motion Field'),
(u'gmf_scenario', u'Ground Motion Field'),
(u'hazard_curve', u'Hazard Curve'),
(u'hazard_curve_multi', u'Hazard Curve (multiple imts)'),
(u'hazard_map', u'Hazard Map'),
(u'ses', u'Stochastic Event Set'),
(u'uh_spectra', u'Uniform Hazard Spectra'),
)
RISK_OUTPUT_TYPE_CHOICES = (
(u'agg_loss_curve', u'Aggregate Loss Curve'),
(u'aggregate_loss', u'Aggregate Losses'),
(u'bcr_distribution', u'Benefit-cost ratio distribution'),
(u'collapse_map', u'Collapse Map Distribution'),
(u'dmg_dist_per_asset', u'Damage Distribution Per Asset'),
(u'dmg_dist_per_taxonomy', u'Damage Distribution Per Taxonomy'),
(u'dmg_dist_total', u'Total Damage Distribution'),
(u'event_loss', u'Event Loss Table'),
(u'event_loss_asset', u'Event Loss Asset'),
(u'loss_curve', u'Loss Curve'),
(u'event_loss_curve', u'Loss Curve'),
(u'loss_fraction', u'Loss fractions'),
(u'loss_map', u'Loss Map'),
)
output_type = djm.TextField(
choices=HAZARD_OUTPUT_TYPE_CHOICES + RISK_OUTPUT_TYPE_CHOICES)
last_update = djm.DateTimeField(editable=False, default=datetime.utcnow)
objects = OutputManager()
def __str__(self):
return "%d||%s||%s" % (self.id, self.output_type, self.display_name)
class Meta:
db_table = 'uiapi\".\"output'
ordering = ['id']
[docs] def is_hazard_curve(self):
return self.output_type in ['hazard_curve', 'hazard_curve_multi']
@property
def output_container(self):
"""
:returns: the output container associated with this output
"""
# FIXME(lp). Remove the following outstanding exceptions
if self.output_type in ['agg_loss_curve', 'event_loss_curve']:
return self.loss_curve
elif self.output_type == 'hazard_curve_multi':
return self.hazard_curve
elif self.output_type == 'gmf_scenario':
return self.gmf
elif self.output_type == 'event_loss_asset':
return self.event_loss
return getattr(self, self.output_type)
@property
def lt_realization_paths(self):
"""
:returns: an instance of `LogicTreePath` the output is
associated with. When the output is not associated with any
logic tree branch then it returns a LogicTreePath namedtuple
with a couple of None.
"""
hazard_output_types = [el[0] for el in self.HAZARD_OUTPUT_TYPE_CHOICES]
risk_output_types = [el[0] for el in self.RISK_OUTPUT_TYPE_CHOICES]
container = self.output_container
if self.output_type in hazard_output_types:
rlz = getattr(container, 'lt_realization_id', None)
if rlz is not None:
return self.LogicTreePath(
tuple(container.lt_realization.gsim_lt_path),
tuple(container.lt_realization.sm_lt_path))
else:
return self.LogicTreePath(None, None)
elif self.output_type in risk_output_types:
if getattr(container, 'hazard_output_id', None):
return container.hazard_output.lt_realization_paths
else:
return self.LogicTreePath(None, None)
raise RuntimeError("unexpected output type %s" % self.output_type)
@property
def statistical_params(self):
"""
:returns: an instance of `StatisticalParams` the output is
associated with
"""
if getattr(self.output_container, 'statistics', None) is not None:
return self.StatisticalParams(self.output_container.statistics,
self.output_container.quantile)
elif getattr(
self.output_container, 'hazard_output_id', None) is not None:
return self.output_container.hazard_output.statistical_params
else:
return self.StatisticalParams(None, None)
@property
def hazard_metadata(self):
"""
Given an Output produced by a risk calculation it returns the
corresponding hazard metadata.
:returns:
A `namedtuple` with the following attributes::
* investigation_time: the hazard investigation time (float)
* statistics: the kind of hazard statistics (None, "mean" or
"quantile")
* quantile: quantile value (when `statistics` is "quantile")
* sm_path: a list representing the source model path
* gsim_path: a list representing the gsim logic tree path
"""
investigation_time = self.oq_job\
.risk_calculation\
.hazard_calculation\
.get_param('investigation_time', None)
statistics, quantile = self.statistical_params
gsim_lt_path, sm_lt_path = self.lt_realization_paths
return self.HazardMetadata(investigation_time,
statistics, quantile,
sm_lt_path, gsim_lt_path)
## Tables in the 'hzrdr' schema.
[docs]class HazardMap(djm.Model):
'''
Hazard Map header (information which pertains to entire map)
'''
output = djm.OneToOneField('Output', related_name="hazard_map")
# FK only required for non-statistical results (i.e., mean or quantile
# curves).
lt_realization = djm.ForeignKey('LtRealization', null=True)
investigation_time = djm.FloatField()
imt = djm.TextField(choices=IMT_CHOICES)
statistics = djm.TextField(null=True, choices=STAT_CHOICES)
quantile = djm.FloatField(null=True)
sa_period = djm.FloatField(null=True)
sa_damping = djm.FloatField(null=True)
poe = djm.FloatField()
# lons, lats, and imls are stored as numpy arrays with a uniform size and
# shape
lons = fields.FloatArrayField()
lats = fields.FloatArrayField()
imls = fields.FloatArrayField()
class Meta:
db_table = 'hzrdr\".\"hazard_map'
def __str__(self):
return (
'HazardMap(poe=%(poe)s, imt=%(imt)s, sa_period=%(sa_period)s, '
'statistics=%(statistics)s, quantile=%(quantile)s)'
) % self.__dict__
def __repr__(self):
return self.__str__()
[docs]class HazardCurve(djm.Model):
'''
Hazard Curve header information
'''
output = djm.OneToOneField(
'Output', null=True, related_name="hazard_curve")
# FK only required for non-statistical results (i.e., mean or quantile
# curves).
lt_realization = djm.ForeignKey('LtRealization', null=True)
investigation_time = djm.FloatField()
imt = djm.TextField(choices=IMT_CHOICES, default=None, blank=True)
imls = fields.FloatArrayField()
STAT_CHOICES = (
(u'mean', u'Mean'),
(u'quantile', u'Quantile'),
)
statistics = djm.TextField(null=True, choices=STAT_CHOICES)
quantile = djm.FloatField(null=True)
sa_period = djm.FloatField(null=True)
sa_damping = djm.FloatField(null=True)
class Meta:
db_table = 'hzrdr\".\"hazard_curve'
@property
def imt_long(self):
"""
:returns: a string representing the imt associated with the
curve (if any) in the long form, e.g. SA(0.01)
"""
if self.imt:
if self.imt == "SA":
return "%s(%s)" % (self.imt, self.sa_damping)
else:
return self.imt
def __iter__(self):
assert self.output.output_type == 'hazard_curve_multi'
siblings = self.__class__.objects.filter(
output__oq_job=self.output.oq_job,
output__output_type='hazard_curve')
if not self.statistics:
return iter(siblings.filter(lt_realization__isnull=False))
elif self.quantile:
return iter(
siblings.filter(statistics="quantile", quantile=self.quantile))
else:
return iter(siblings.filter(statistics="mean"))
[docs]class HazardCurveDataManager(djm.GeoManager):
"""
Manager class to filter and create HazardCurveData objects
"""
[docs] def all_curves_for_imt(self, job, imt, sa_period, sa_damping):
"""
Helper function for creating a :class:`django.db.models.query.QuerySet`
for selecting all curves from all realizations for a given ``job_id``
and ``imt``.
:param job:
An :class:`openquake.engine.db.models.OqJob` instance.
:param str imt:
Intensity measure type.
:param sa_period:
Spectral Acceleration period value. Only relevant if the ``imt`` is
"SA".
:param sa_damping:
Spectrail Acceleration damping value. Only relevant if the ``imt``
is "SA".
"""
return self.filter(hazard_curve__output__oq_job=job,
hazard_curve__imt=imt,
hazard_curve__sa_period=sa_period,
hazard_curve__sa_damping=sa_damping,
# We only want curves associated with a logic tree
# realization (and not statistical aggregates):
hazard_curve__lt_realization__isnull=False)
[docs] def all_curves_simple(self, filter_args=None, order_by='id'):
"""
Get all :class:`HazardCurveData` records matching `filter_args` and
return the results in a simple, lean format: a sequence of (x, y, poes)
triples, where x and y are longitude and latitude of the `location`.
For querying large sets of hazard curve data, this is a rather lean
and efficient method for getting the results.
:param dict filter_args:
Optional. Dictionary of filter arguments to apply to the query.
:param str order_by:
Defaults to the primary key ('id'). Field by which to order
results. Currently, only one `ORDER BY` field is supported.
"""
if filter_args is None:
filter_args = dict()
return self\
.filter(**filter_args)\
.order_by(order_by)\
.extra(select={'x': 'ST_X(location)', 'y': 'ST_Y(location)'})\
.values_list('x', 'y', 'poes')\
.iterator()
[docs]class HazardCurveData(djm.Model):
'''
Hazard Curve data
Contains an list of PoE (Probability of Exceedance)
values and the geographical point associated with the curve
'''
hazard_curve = djm.ForeignKey('HazardCurve')
poes = fields.FloatArrayField()
location = djm.PointField(srid=DEFAULT_SRID)
# weight can be null/None if the weight is implicit:
weight = djm.DecimalField(decimal_places=100, max_digits=101, null=True)
objects = HazardCurveDataManager()
class Meta:
db_table = 'hzrdr\".\"hazard_curve_data'
[docs]class SESCollection(djm.Model):
"""
Stochastic Event Set Collection: A container for 1 or more Stochastic Event
Sets for a given logic tree realization.
See also :class:`SES` and :class:`SESRupture`.
"""
output = djm.OneToOneField('Output', related_name="ses")
trt_model = djm.OneToOneField(
'TrtModel', related_name='ses_collection', null=False)
ordinal = djm.IntegerField(null=False)
class Meta:
db_table = 'hzrdr\".\"ses_collection'
ordering = ['ordinal']
@classmethod
[docs] def create(cls, output):
"""
Create an LtSourceModel, a TrtModel and a SESCollection associated
to it and to the given output.
:param output: an output of type GMF
"""
lt_model = LtSourceModel.objects.create(
hazard_calculation=output.oq_job, ordinal=0,
sm_lt_path=[], sm_name='fake-from-rupture', weight=1)
# in order to save a ProbabilisticRupture, a TrtModel is needed;
# here we generate a fake one, corresponding to the tectonic
# region type NA i.e. Not Available
trt_model = TrtModel.objects.create(
lt_model=lt_model, tectonic_region_type='NA', num_sources=0,
num_ruptures=1, min_mag=0, max_mag=0, gsims=[])
return cls.objects.create(
output=output, trt_model=trt_model, ordinal=0)
@property
def sm_lt_path(self):
"""
The source model logic tree path corresponding to the collection
"""
return tuple(self.trt_model.lt_model.sm_lt_path)
[docs] def get_ruptures(self):
"""Return the SESRuptures associated to self"""
return SESRupture.objects.filter(rupture__ses_collection=self.id)
def __iter__(self):
"""
Iterator for walking through all child :class:`SES` objects.
"""
n = self.output.oq_job.get_param('ses_per_logic_tree_path', 1)
for ordinal in xrange(1, n + 1): # 1 for scenario
yield SES(self, ordinal)
def __len__(self):
"""
Return the ses_per_logic_tree_path parameter
"""
return self.output.oq_job.get_param('ses_per_logic_tree_path', 1)
def __repr__(self):
return '<%s=%d, trt_model=%s, ordinal=%d>' % (
self.__class__.__name__, self.id, self.trt_model.id, self.ordinal)
[docs]class SES(object):
"""
Stochastic Event Set: A container for 1 or more ruptures associated with a
specific investigation time span.
"""
# the ordinal must be > 0: the reason is that it appears in the
# exported XML file and the schema constraints the number to be
# nonzero
def __init__(self, ses_collection, ordinal=1):
self.ses_collection = ses_collection
self.ordinal = ordinal
self.investigation_time = self.ses_collection.output.oq_job.get_param(
'investigation_time', None)
def __cmp__(self, other):
return cmp(self.ordinal, other.ordinal)
def __iter__(self):
"""
Iterator for walking through all child :class:`SESRupture` objects.
"""
return SESRupture.objects.filter(
rupture__ses_collection=self.ses_collection.id,
ses_id=self.ordinal).order_by('tag').iterator()
[docs]def get_geom(surface, is_from_fault_source, is_multi_surface):
"""
The following fields can be interpreted different ways,
depending on the value of `is_from_fault_source`. If
`is_from_fault_source` is True, each of these fields should
contain a 2D numpy array (all of the same shape). Each triple
of (lon, lat, depth) for a given index represents the node of
a rectangular mesh. If `is_from_fault_source` is False, each
of these fields should contain a sequence (tuple, list, or
numpy array, for example) of 4 values. In order, the triples
of (lon, lat, depth) represent top left, top right, bottom
left, and bottom right corners of the the rupture's planar
surface. Update: There is now a third case. If the rupture
originated from a characteristic fault source with a
multi-planar-surface geometry, `lons`, `lats`, and `depths`
will contain one or more sets of 4 points, similar to how
planar surface geometry is stored (see above).
:param rupture: an instance of :class:
`openquake.hazardlib.source.rupture.BaseProbabilisticRupture`
:param is_from_fault_source: a boolean
:param is_multi_surface: a boolean
"""
if is_from_fault_source:
# for simple and complex fault sources,
# rupture surface geometry is represented by a mesh
surf_mesh = surface.get_mesh()
lons = surf_mesh.lons
lats = surf_mesh.lats
depths = surf_mesh.depths
else:
if is_multi_surface:
# `list` of
# openquake.hazardlib.geo.surface.planar.PlanarSurface
# objects:
surfaces = surface.surfaces
# lons, lats, and depths are arrays with len == 4*N,
# where N is the number of surfaces in the
# multisurface for each `corner_*`, the ordering is:
# - top left
# - top right
# - bottom left
# - bottom right
lons = numpy.concatenate([x.corner_lons for x in surfaces])
lats = numpy.concatenate([x.corner_lats for x in surfaces])
depths = numpy.concatenate([x.corner_depths for x in surfaces])
else:
# For area or point source,
# rupture geometry is represented by a planar surface,
# defined by 3D corner points
lons = numpy.zeros((4))
lats = numpy.zeros((4))
depths = numpy.zeros((4))
# NOTE: It is important to maintain the order of these
# corner points. TODO: check the ordering
for i, corner in enumerate((surface.top_left,
surface.top_right,
surface.bottom_left,
surface.bottom_right)):
lons[i] = corner.longitude
lats[i] = corner.latitude
depths[i] = corner.depth
return lons, lats, depths
[docs]class ProbabilisticRupture(djm.Model):
"""
A rupture as part of a Stochastic Event Set Collection.
"""
ses_collection = djm.ForeignKey('SESCollection')
magnitude = djm.FloatField(null=False)
_hypocenter = fields.FloatArrayField(null=False)
rake = djm.FloatField(null=False)
is_from_fault_source = djm.NullBooleanField(null=False)
is_multi_surface = djm.NullBooleanField(null=False)
surface = fields.PickleField(null=False)
site_indices = fields.IntArrayField(null=True)
# NB (MS): the proper solution would be to store the hypocenter as a 3D
# point, however I was unable to do so, due to a bug in Django 1.3
# (I was getting a GeometryProxy exception).
# The GEOS library we are using does not even support
# the WKT for 3D points; that's why I am storing the point as a
# 3D array, as a workaround; luckily, we never perform any
# geospatial query on the hypocenter.
# we will be able to do better when we will upgrade (Geo)Django
@property
def hypocenter(self):
"""Convert the 3D array into a hazardlib point"""
return geo.Point(*self._hypocenter)
class Meta:
db_table = 'hzrdr\".\"probabilistic_rupture'
@classmethod
[docs] def create(cls, rupture, ses_collection, site_indices=None):
"""
Create a ProbabilisticRupture row on the database.
:param rupture:
a hazardlib rupture
:param ses_collection:
a Stochastic Event Set Collection object
:param site_indices:
an array of indices for the site_collection
"""
iffs = isinstance(rupture.surface,
(geo.ComplexFaultSurface, geo.SimpleFaultSurface))
ims = isinstance(rupture.surface, geo.MultiSurface)
lons, lats, depths = get_geom(rupture.surface, iffs, ims)
hp = rupture.hypocenter
return cls.objects.create(
ses_collection=ses_collection,
magnitude=rupture.mag,
rake=rupture.rake,
is_from_fault_source=iffs,
is_multi_surface=ims,
surface=rupture.surface,
_hypocenter=[hp.longitude, hp.latitude, hp.depth],
site_indices=site_indices)
@property
def tectonic_region_type(self):
"""The TRT associated to the underlying trt_model"""
return self.ses_collection.trt_model.tectonic_region_type
_geom = None
@property
def geom(self):
"""
Extract the triple (lons, lats, depths) from the surface geometry
(cached).
"""
if self._geom is not None:
return self._geom
self._geom = get_geom(self.surface, self.is_from_fault_source,
self.is_multi_surface)
return self._geom
@property
def lons(self):
return self.geom[0]
@property
def lats(self):
return self.geom[1]
@property
def depths(self):
return self.geom[2]
@property
def strike(self):
return self.surface.get_strike()
@property
def dip(self):
return self.surface.get_dip()
@property
def mag(self):
return self.magnitude
def _validate_planar_surface(self):
"""
A rupture's planar surface (existing only in the case of ruptures from
area/point sources) may only consist of 4 points (top left, top right,
bottom right, and bottom left corners, in that order).
If the surface is not valid, a :exc:`ValueError` is raised.
This should only be used if `is_from_fault_source` is `False`.
"""
if not (4 == len(self.lons) == len(self.lats) == len(self.depths)):
raise ValueError(
"Incorrect number of points; there should be exactly 4")
@property
def top_left_corner(self):
if not (self.is_from_fault_source or self.is_multi_surface):
self._validate_planar_surface()
return self.lons[0], self.lats[0], self.depths[0]
return None
@property
def top_right_corner(self):
if not (self.is_from_fault_source or self.is_multi_surface):
self._validate_planar_surface()
return self.lons[1], self.lats[1], self.depths[1]
return None
@property
def bottom_left_corner(self):
if not (self.is_from_fault_source or self.is_multi_surface):
self._validate_planar_surface()
return self.lons[2], self.lats[2], self.depths[2]
return None
@property
def bottom_right_corner(self):
if not (self.is_from_fault_source or self.is_multi_surface):
self._validate_planar_surface()
return self.lons[3], self.lats[3], self.depths[3]
return None
[docs]class SESRupture(djm.Model):
"""
A rupture as part of a Stochastic Event Set.
"""
rupture = djm.ForeignKey('ProbabilisticRupture')
ses_id = djm.IntegerField(null=False)
tag = djm.TextField(null=False)
seed = djm.IntegerField(null=False)
class Meta:
db_table = 'hzrdr\".\"ses_rupture'
ordering = ['tag']
@classmethod
[docs] def create(cls, prob_rupture, ses_ordinal, source_id, rupt_no, rupt_occ,
seed):
"""
Create a SESRupture row in the database.
:param prob_rupture:
:class:`openquake.engine.db.models.ProbabilisticRupture` instance
:param int ses_ordinal:
ordinal for a :class:`openquake.engine.db.models.SES` instance
:param str source_id:
id of the source that generated the rupture
:param rupt_no:
the rupture number (an ordinal from source.iter_ruptures())
:param rupt_occ:
the occurrence number of the rupture in the given ses
:param int seed:
a seed that will be used when computing the GMF from the rupture
"""
tag = 'trt=%02d|ses=%04d|src=%s|rup=%03d-%02d' % (
prob_rupture.ses_collection.ordinal, ses_ordinal,
source_id, rupt_no, rupt_occ)
return cls.objects.create(
rupture=prob_rupture, ses_id=ses_ordinal, tag=tag, seed=seed)
@property
def surface(self):
"""The surface of the underlying rupture"""
return self.rupture.surface
@property
def hypocenter(self):
"""The hypocenter of the underlying rupture"""
return self.rupture.hypocenter
_Point = collections.namedtuple('_Point', 'x y')
[docs]def get_correl_model(job):
"""
Helper function for constructing the appropriate correlation model.
:returns:
A correlation object. See :mod:`openquake.hazardlib.correlation`
for more info.
"""
correl_model_name = job.get_param('ground_motion_correlation_model', None)
if correl_model_name is None:
# There's no correlation model for this calculation.
return None
correl_model_cls = getattr(
correlation, '%sCorrelationModel' % correl_model_name, None)
if correl_model_cls is None:
# There's no correlation model for this calculation.
return None
gmc_params = job.get_param('ground_motion_correlation_params', None)
return correl_model_cls(**gmc_params)
[docs]class Gmf(djm.Model):
"""
A collection of ground motion field (GMF) sets for a given logic tree
realization.
"""
output = djm.OneToOneField('Output', related_name="gmf")
lt_realization = djm.ForeignKey('LtRealization', null=False)
class Meta:
db_table = 'hzrdr\".\"gmf'
[docs] def check_export_size(self):
"""
Raise an error if the number of rows to export is bigger that
the configuration parameter `max_rows_export_gmfs`.
"""
max_rows = int(utils.config.get('hazard', 'max_rows_export_gmfs'))
if max_rows:
num_rows = GmfData.objects.filter(gmf=self).count()
if num_rows > max_rows:
raise RuntimeError(
'Cannot export the GMFs: the `max_rows_export_gmf limit` '
'is set to %d rows, but there are %d rows to export' % (
max_rows, num_rows))
def __iter__(self):
"""
Get the ground motion fields per SES ("GMF set") for
the XML export. Each "GMF set" should:
* have an `investigation_time` attribute
* have an `stochastic_event_set_id` attribute
* be iterable, yielding a sequence of "GMF" objects
Each "GMF" object should:
* have an `imt` attribute
* have an `sa_period` attribute (only if `imt` is 'SA')
* have an `sa_damping` attribute (only if `imt` is 'SA')
* have a `rupture_id` attribute (to indicate which rupture
contributed to this gmf)
* be iterable, yielding a sequence of "GMF node" objects
Each "GMF node" object should have:
* a `gmv` attribute (to indicate the ground motion value
* `lon` and `lat` attributes (to indicate the geographical location
of the ground motion field)
If a SES does not generate any GMF, it is ignored.
"""
self.check_export_size()
job = self.output.oq_job
curs = getcursor('job_init')
# extract all the gmf_data for the current realization
logs.LOG.info('reading gmf_data for gmf_id=%s', self.id)
curs.execute("""\
SELECT site_id, imt, sa_period, sa_damping, gmvs, rupture_ids
FROM hzrdr.gmf_data
WHERE gmf_id=%d
ORDER BY site_id, imt, sa_period, sa_damping""" % self.id)
gmfdata = curs.fetchall()
# find all the locations interested by the hazard calculation
logs.LOG.info('reading HazardSite for job_id=%d', job.id)
curs.execute("""\
SELECT id, ST_X(location::geometry), ST_Y(location::geometry)
FROM hzrdi.hazard_site WHERE hazard_calculation_id=%d""" % job.id)
loc = {site_id: (x, y) for site_id, x, y in curs}
for ses_coll in SESCollection.objects.filter(output__oq_job=job):
ruptag = dict(SESRupture.objects.filter(
rupture__ses_collection=ses_coll
).values_list('id', 'tag'))
gmf_dict = collections.defaultdict(list) # imt, tag -> [gmv-x-y]
for site_id, imt, sa_period, sa_damping, gmvs, rupture_ids \
in gmfdata:
x, y = loc[site_id]
for gmv, rupid in zip(gmvs, rupture_ids):
try:
tag = ruptag[rupid]
gmf_dict[imt, sa_period, sa_damping, tag].append(
(x, y, gmv))
except KeyError:
pass
logs.LOG.info('reordered GMF data for SES collection %d',
ses_coll.ordinal)
for gmfset in self.gmfsets(ses_coll, gmf_dict):
yield gmfset
[docs] def gmfsets(self, ses_coll, gmf_dict):
"""
:param ses_coll:
a SESCollection instance
:param gmf_dict:
a dictionary (imt, sa_period, sa_damping, rupture_tag)
-> [(x, y, gmv), ...]
:returns:
a list of :class:`openquake.engine.db.models.GmfSet` instances
"""
# a set of GMFs generate by the same SES, one per rupture
gmfset = collections.defaultdict(list) # ses_ordinal -> GMFs
for key in sorted(gmf_dict):
imt, sa_period, sa_damping, rupture_tag = key
# using a generator here saves a lot of memory
nodes = (_GroundMotionFieldNode(gmv, _Point(x, y))
for x, y, gmv in gmf_dict[key])
ses_ordinal = extract_ses_ordinal(rupture_tag)
gmfset[ses_ordinal].append(
_GroundMotionField(
imt, sa_period, sa_damping, rupture_tag, nodes))
return [GmfSet(ses, gmfset[ses.ordinal]) for ses in ses_coll
if ses.ordinal in gmfset]
[docs]class GmfSet(object):
"""
Small wrapper around the list of Gmf objects associated to the given SES.
"""
def __init__(self, ses, gmfset):
self.gmfset = gmfset
self.investigation_time = ses.investigation_time or 0
self.stochastic_event_set_id = ses.ordinal
def __iter__(self):
return iter(self.gmfset)
def __nonzero__(self):
return bool(self.gmfset)
def __str__(self):
return (
'GMFsPerSES(investigation_time=%f, '
'stochastic_event_set_id=%s,\n%s)' % (
self.investigation_time,
self.stochastic_event_set_id, '\n'.join(
sorted(str(g) for g in self.gmfset))))
class _GroundMotionField(object):
"""
The Ground Motion Field generated by the given rupture
"""
def __init__(self, imt, sa_period, sa_damping, rupture_id, gmf_nodes):
self.imt = imt
self.sa_period = sa_period
self.sa_damping = sa_damping
self.rupture_id = rupture_id
self.gmf_nodes = gmf_nodes
def __iter__(self):
return iter(self.gmf_nodes)
def __getitem__(self, key):
return self.gmf_nodes[key]
def __str__(self):
"""
String representation of a _GroundMotionField object showing the
content of the nodes (lon, lat an gmv). This is useful for debugging
and testing.
"""
mdata = ('imt=%(imt)s sa_period=%(sa_period)s '
'sa_damping=%(sa_damping)s rupture_id=%(rupture_id)s' %
vars(self))
nodes = sorted(map(str, self.gmf_nodes))
return 'GMF(%s\n%s)' % (mdata, '\n'.join(nodes))
class _GroundMotionFieldNode(object):
# the signature is not (gmv, x, y) because the XML writer expect a location
# object
def __init__(self, gmv, loc):
self.gmv = gmv
self.location = loc
def __lt__(self, other):
"""
A reproducible ordering by lon and lat; used in
:function:`openquake.commonlib.hazard_writers.gen_gmfs`
"""
return self.location < other.location
def __str__(self):
"""Return lon, lat and gmv of the node in a compact string form"""
return '<X=%9.5f, Y=%9.5f, GMV=%9.7f>' % (
self.location.x, self.location.y, self.gmv)
[docs]class GmfData(djm.Model):
"""
Ground Motion Field: A collection of ground motion values and their
respective geographical locations.
"""
gmf = djm.ForeignKey('Gmf')
task_no = djm.IntegerField(null=False)
imt = djm.TextField(choices=IMT_CHOICES)
sa_period = djm.FloatField(null=True)
sa_damping = djm.FloatField(null=True)
gmvs = fields.FloatArrayField()
rupture_ids = fields.IntArrayField(null=True)
site = djm.ForeignKey('HazardSite')
objects = djm.GeoManager()
class Meta:
db_table = 'hzrdr\".\"gmf_data'
ordering = ['gmf', 'task_no']
[docs]def get_gmvs_per_site(output, imt):
"""
Iterator for walking through all :class:`GmfData` objects associated
to a given output. Notice that values for the same site are
displayed together and ordered according to the rupture ids, so that
it is possible to get consistent outputs in the test cases.
:param output: instance of :class:`openquake.engine.db.models.Output`
:param string imt: a string with the IMT to extract
:returns: a list of ground motion values per each site
"""
imtype, sa_period, sa_damping = from_string(imt)
gmf_id = output.gmf.id
curs = getcursor('job_init')
query = '''\
SELECT site_id, array_concat(gmvs), array_concat(rupture_ids)
FROM hzrdr.gmf_data WHERE gmf_id=%s AND imt=%s {}
GROUP BY site_id ORDER BY site_id'''
if imtype == 'SA':
curs.execute(query.format('AND sa_period=%s AND sa_damping=%s'),
(gmf_id, imtype, sa_period, sa_damping))
else:
curs.execute(query.format(''), (gmf_id, imtype))
for site_id, gmvs, rup_ids in curs:
gmv = dict(zip(rup_ids, gmvs))
yield [gmv[r] for r in sorted(rup_ids)]
[docs]class DisaggResult(djm.Model):
"""
Storage for disaggregation historgrams. Each histogram is stored in
`matrix` as a 6-dimensional numpy array (pickled). The dimensions of the
matrix are as follows, in order:
* magnitude
* distance
* longitude
* latitude
* epsilon
* tectonic region type
Bin edges are defined for all of these dimensions (except tectonic region
type) as:
* `mag_bin_edges`
* `dist_bin_edges`
* `lat_bin_edges`
* `lon_bin_edges`
* `eps_bin_edges`
The size of the tectonic region type (TRT) dimension is simply determined
by the length of `trts`.
Additional metadata for the disaggregation histogram is stored, including
location (POINT geometry), disaggregation PoE (Probability of Exceedance)
and the corresponding IML (Intensity Measure Level) extracted from the
hazard curve, logic tree path information, and investigation time.
"""
output = djm.OneToOneField('Output', related_name="disagg_matrix")
lt_realization = djm.ForeignKey('LtRealization')
investigation_time = djm.FloatField()
imt = djm.TextField(choices=IMT_CHOICES)
iml = djm.FloatField()
poe = djm.FloatField()
sa_period = djm.FloatField(null=True)
sa_damping = djm.FloatField(null=True)
mag_bin_edges = fields.FloatArrayField()
dist_bin_edges = fields.FloatArrayField()
lon_bin_edges = fields.FloatArrayField()
lat_bin_edges = fields.FloatArrayField()
eps_bin_edges = fields.FloatArrayField()
trts = fields.CharArrayField()
location = djm.PointField(srid=DEFAULT_SRID)
matrix = fields.PickleField()
class Meta:
db_table = 'hzrdr\".\"disagg_result'
[docs]class UHS(djm.Model):
"""
UHS/Uniform Hazard Spectra:
* "Uniform" meaning "the same PoE"
* "Spectrum" because it covers a range/band of periods/frequencies
Records in this table contain metadata for a collection of UHS data.
"""
output = djm.OneToOneField('Output', null=True, related_name="uh_spectra")
# FK only required for non-statistical results (i.e., mean or quantile
# curves).
lt_realization = djm.ForeignKey('LtRealization', null=True)
investigation_time = djm.FloatField()
poe = djm.FloatField()
periods = fields.FloatArrayField()
STAT_CHOICES = (
(u'mean', u'Mean'),
(u'quantile', u'Quantile'),
)
statistics = djm.TextField(null=True, choices=STAT_CHOICES)
quantile = djm.FloatField(null=True)
class Meta:
db_table = 'hzrdr\".\"uhs'
def __iter__(self):
"""
Iterate over the :class:`UHSData` which belong this object.
"""
return self.uhsdata_set.iterator()
[docs]class UHSData(djm.Model):
"""
UHS curve for a given location.
"""
uhs = djm.ForeignKey('UHS')
imls = fields.FloatArrayField()
location = djm.PointField(srid=DEFAULT_SRID)
class Meta:
db_table = 'hzrdr\".\"uhs_data'
[docs]class LtSourceModel(djm.Model):
"""
Identify a logic tree source model.
"""
hazard_calculation = djm.ForeignKey('OqJob')
ordinal = djm.IntegerField()
sm_lt_path = fields.CharArrayField()
sm_name = djm.TextField(null=False)
weight = djm.DecimalField(decimal_places=100, max_digits=101, null=True)
class Meta:
db_table = 'hzrdr\".\"lt_source_model'
ordering = ['ordinal']
[docs] def get_num_sources(self):
"""
Return the number of sources in the model.
"""
return sum(info.num_sources for info in self.trtmodel_set.all())
[docs] def get_tectonic_region_types(self):
"""
Return the tectonic region types in the model,
ordered by number of sources.
"""
return self.trtmodel_set.filter(
lt_model=self, num_ruptures__gt=0).values_list(
'tectonic_region_type', flat=True)
[docs] def make_gsim_lt(self, trts=()):
"""
Helper to instantiate a GsimLogicTree object from the logic tree file.
"""
hc = self.hazard_calculation.get_oqparam()
trts = trts or self.get_tectonic_region_types()
fname = os.path.join(hc.base_path, hc.inputs['gsim_logic_tree'])
gsim_lt = logictree.GsimLogicTree(
fname, 'applyToTectonicRegionType', trts)
for trt in trts:
if not trt in gsim_lt.values:
raise ValueError(
"Found in %r a tectonic region type %r inconsistent with "
"the ones in %r" % (self.sm_name, trt, fname))
return gsim_lt
def __iter__(self):
"""
Yield the realizations corresponding to the given model
"""
return iter(self.ltrealization_set.all())
[docs]class TrtModel(djm.Model):
"""
Source submodel containing sources of the same tectonic region type.
"""
lt_model = djm.ForeignKey('LtSourceModel')
tectonic_region_type = djm.TextField(null=False)
num_sources = djm.IntegerField(null=False)
num_ruptures = djm.IntegerField(null=False)
min_mag = djm.FloatField(null=False)
max_mag = djm.FloatField(null=False)
gsims = fields.CharArrayField(null=True)
[docs] def get_realizations(self, gsim_name):
"""
Return the realizations associated to the current TrtModel and
the given GSIM.
:param str gsim_name: name of a GSIM class
"""
assert gsim_name in self.gsims, gsim_name
for art in AssocLtRlzTrtModel.objects.filter(
trt_model=self.id, gsim=gsim_name):
yield art.rlz
[docs] def get_rlzs_by_gsim(self):
"""
Return the realizations associated to the current TrtModel
as an ordered dictionary {gsim_name: [rlz, ...]}
"""
dic = collections.defaultdict(list)
for art in AssocLtRlzTrtModel.objects.filter(trt_model=self.id):
dic[art.gsim].append(art.rlz)
return collections.OrderedDict(
(gsim, dic[gsim]) for gsim in sorted(dic))
[docs] def get_gsim_instances(self):
"""
Return the GSIM instances associated to the current TrtModel
"""
return [logictree.GSIM[gsim]() for gsim in self.gsims]
class Meta:
db_table = 'hzrdr\".\"trt_model'
ordering = ['id']
# NB: the TrtModels are built in the right order, see
# BaseHazardCalculator.initialize_sources
[docs]class SourceInfo(djm.Model):
"""
Source specific infos
"""
trt_model = djm.ForeignKey('TrtModel')
source_id = djm.TextField(null=False)
source_class = djm.TextField(null=False)
num_sites = djm.IntegerField(null=False)
num_ruptures = djm.IntegerField(null=False)
occ_ruptures = djm.IntegerField(null=True)
uniq_ruptures = djm.IntegerField(null=True)
calc_time = djm.FloatField(null=False)
class Meta:
db_table = 'hzrdr\".\"source_info'
ordering = ['trt_model', 'source_id']
[docs]class AssocLtRlzTrtModel(djm.Model):
"""
Associations between logic tree realizations and TrtModels. Fixed
a realization and a TRT, the gsim is unique.
"""
rlz = djm.ForeignKey('LtRealization')
trt_model = djm.ForeignKey('TrtModel')
gsim = djm.TextField(null=False)
class Meta:
db_table = 'hzrdr\".\"assoc_lt_rlz_trt_model'
ordering = ['id']
[docs]class LtRealization(djm.Model):
"""
Identify a logic tree branch.
"""
lt_model = djm.ForeignKey('LtSourceModel')
ordinal = djm.IntegerField()
weight = djm.DecimalField(decimal_places=100, max_digits=101)
gsim_lt_path = fields.CharArrayField()
@property
def sm_lt_path(self):
"""
The source model logic tree path extracted from the underlying
source model
"""
return self.lt_model.sm_lt_path
class Meta:
db_table = 'hzrdr\".\"lt_realization'
ordering = ['ordinal']
[docs] def get_gsim_instances(self):
"""
Return the GSIM instances associated to the current realization
by looking at the association table.
"""
return [logictree.GSIM[art.gsim]() for art in
AssocLtRlzTrtModel.objects.filter(rlz=self)]
## Tables in the 'riskr' schema.
[docs]class LossFraction(djm.Model):
"""
Holds metadata for loss fraction data
"""
output = djm.OneToOneField("Output", related_name="loss_fraction")
variable = djm.TextField(choices=(("taxonomy", "taxonomy"),
("magnitude_distance",
"Magnitude Distance"),
("coordinate", "Coordinate")))
hazard_output = djm.ForeignKey(
"Output", related_name="risk_loss_fractions")
statistics = djm.TextField(null=True, choices=STAT_CHOICES)
quantile = djm.FloatField(null=True)
poe = djm.FloatField(null=True)
loss_type = djm.TextField(choices=zip(LOSS_TYPES, LOSS_TYPES))
class Meta:
db_table = 'riskr\".\"loss_fraction'
def __iter__(self):
return iter(self.lossfractiondata_set.all())
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
return (self.output.output_type,
self.output.hazard_metadata,
self.statistics, self.quantile,
self.variable, self.poe, self.loss_type)
[docs] def display_value(self, value, rc):
"""
Converts `value` in a form that is best suited to be
displayed.
:param rc:
A `RiskCalculation` object used to get the bin width
:returns: `value` if the attribute `variable` is equal to
taxonomy. if the attribute `variable` is equal to
`magnitude-distance`, then it extracts two integers (comma
separated) from `value` and convert them into ranges
encoded back as csv.
"""
if self.variable == "taxonomy":
return value
elif self.variable == "magnitude_distance":
magnitude, distance = map(float, value.split(","))
return "%.4f,%.4f|%.4f,%.4f" % (
magnitude * rc.mag_bin_width,
(magnitude + 1) * rc.mag_bin_width,
distance * rc.distance_bin_width,
(distance + 1) * rc.distance_bin_width)
elif self.variable == "coordinate":
lon, lat = map(float, value.split(","))
return "%.4f,%.4f|%.4f,%.4f" % (
lon * rc.coordinate_bin_width,
(lon + 1) * rc.coordinate_bin_width,
lat * rc.coordinate_bin_width,
(lat + 1) * rc.coordinate_bin_width)
else:
raise RuntimeError(
"disaggregation of type %s not supported" % self.variable)
[docs] def total_fractions(self):
"""
:returns: a dictionary mapping values of `variable` (e.g. a
taxonomy) to tuples yielding the associated absolute losses
(e.g. the absolute losses for assets of a taxonomy) and the
percentage (expressed in decimal format) over the total losses
"""
cursor = connections['job_init'].cursor()
total = self.lossfractiondata_set.aggregate(
djm.Sum('absolute_loss')).values()[0]
if not total:
return {}
query = """
SELECT value, sum(absolute_loss)
FROM riskr.loss_fraction_data
WHERE loss_fraction_id = %s
GROUP BY value
"""
cursor.execute(query, (self.id,))
rc = self.output.oq_job.risk_calculation
loss_fraction = collections.namedtuple('loss_fraction', 'bin loss')
return collections.OrderedDict(
sorted(
[loss_fraction(
self.display_value(value, rc),
(loss, loss / total))
for value, loss in cursor],
key=operator.attrgetter('loss'),
reverse=True))
[docs] def iteritems(self):
"""
Yields tuples with two elements. The first one is a location
(described by a lon/lat tuple), the second one is a dictionary
modeling the disaggregation of the losses on such location. In
this dictionary, each key is a value of `variable`, and each
corresponding value is a tuple holding the absolute losses and
the fraction of losses occurring in that location.
"""
rc = self.output.oq_job.risk_calculation
cursor = connections['job_init'].cursor()
# Partition by lon,lat because partitioning on geometry types
# seems not supported in postgis 1.5
query = """
SELECT lon, lat, value,
fraction_loss,
SUM(fraction_loss) OVER w,
COUNT(*) OVER w
FROM (SELECT ST_X(location) as lon,
ST_Y(location) as lat,
value, sum(absolute_loss) as fraction_loss
FROM riskr.loss_fraction_data
WHERE loss_fraction_id = %s
GROUP BY location, value) g
WINDOW w AS (PARTITION BY lon, lat)
"""
cursor.execute(query, (self.id, ))
def display_value_and_fractions(value, absolute_loss, total_loss):
display_value = self.display_value(value, rc)
if total_loss > 0:
fraction = absolute_loss / total_loss
else:
# When a rupture is associated with a positive ground
# shaking (gmv > 0) but with a loss = 0, we still
# store this information. In that case, total_loss =
# absolute_loss = 0
fraction = 0
return display_value, fraction
# We iterate on loss fraction data by location in two steps.
# First we fetch a loss fraction for a single location and a
# single value. In the same query we get the number `count` of
# bins stored for such location. Then, we fetch `count` - 1
# fractions to finalize the fractions on the current location.
while 1:
data = cursor.fetchone()
if data is None:
raise StopIteration
lon, lat, value, absolute_loss, total_loss, count = data
display_value, fraction = display_value_and_fractions(
value, absolute_loss, total_loss)
node = [(lon, lat),
{display_value: (absolute_loss, fraction)}]
data = cursor.fetchmany(count - 1)
for lon, lat, value, absolute_loss, total_loss, count in data:
display_value, fraction = display_value_and_fractions(
value, absolute_loss, total_loss)
node[1][display_value] = (absolute_loss, fraction)
node[1] = collections.OrderedDict(
sorted([(k, v) for k, v in node[1].items()],
key=lambda kv: kv[0]))
yield node
[docs] def to_array(self):
"""
:returns: the loss fractions as numpy array
:NOTE: (not memory efficient)
"""
def to_tuple():
for (lon, lat), data in self.iteritems():
for taxonomy, (absolute_loss, fraction) in data.items():
yield lon, lat, taxonomy, absolute_loss, fraction
return numpy.array(list(to_tuple()), dtype='f4, f4, S3, f4, f4')
[docs]class LossFractionData(djm.Model):
loss_fraction = djm.ForeignKey(LossFraction)
location = djm.PointField(srid=DEFAULT_SRID)
value = djm.TextField()
absolute_loss = djm.TextField()
class Meta:
db_table = 'riskr\".\"loss_fraction_data'
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return (self.loss_fraction.output_hash +
("%.5f" % self.location.x, "%.5f" % self.location.y,
self.value))
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(
self, data, operator.attrgetter('absolute_loss'))
def to_csv_str(self):
"""
Convert LossFractionData into a CSV string
"""
return '%.5f,%.5f,%s,%s' % (
self.location.x, self.location.y, self.value, self.absolute_loss)
[docs] def to_csv_str(self):
"""
Convert LossFraction into a CSV string
"""
return '%.5f,%.5f,%s,%s' % (
self.location.x, self.location.y, self.value, self.absolute_loss)
[docs]class LossMap(djm.Model):
'''
Holds metadata for loss maps
'''
output = djm.OneToOneField("Output", related_name="loss_map")
hazard_output = djm.ForeignKey("Output", related_name="risk_loss_maps")
insured = djm.BooleanField(default=False)
poe = djm.FloatField(null=True)
statistics = djm.TextField(null=True, choices=STAT_CHOICES)
quantile = djm.FloatField(null=True)
loss_type = djm.TextField(choices=zip(LOSS_TYPES, LOSS_TYPES))
class Meta:
db_table = 'riskr\".\"loss_map'
def __iter__(self):
return iter(self.lossmapdata_set.all())
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
return (self.output.output_type,
self.output.hazard_metadata,
self.statistics, self.quantile,
self.insured, self.poe, self.loss_type)
[docs]class LossMapData(djm.Model):
'''
Holds an asset, its position and a value plus (for
non-scenario maps) the standard deviation for its loss
'''
loss_map = djm.ForeignKey("LossMap")
asset_ref = djm.TextField()
value = djm.FloatField()
std_dev = djm.FloatField(default=0.0, null=True)
location = djm.PointField(srid=DEFAULT_SRID)
class Meta:
db_table = 'riskr\".\"loss_map_data'
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.loss_map.output_hash + (self.asset_ref,)
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(
self, data, operator.attrgetter('value', 'stddev'))
[docs]class AggregateLoss(djm.Model):
output = djm.OneToOneField("Output", related_name="aggregate_loss")
insured = djm.BooleanField(default=False)
mean = djm.FloatField()
std_dev = djm.FloatField()
loss_type = djm.TextField(choices=zip(LOSS_TYPES, LOSS_TYPES))
class Meta:
db_table = 'riskr\".\"aggregate_loss'
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
return (self.output.output_type,
self.output.hazard_metadata,
self.insured,
self.mean, self.std_dev,
self.loss_type)
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.output_hash
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(
self, data, lambda x: operator.attrgetter('mean', 'std_dev'))
[docs] def to_csv_str(self):
"""
Convert AggregateLoss into a CSV string
"""
return '\n'.join(data.to_csv_str('row-%d' % i)
for i, data in enumerate(self, 1))
[docs]class LossCurve(djm.Model):
'''
Holds the parameters common to a set of loss curves
'''
output = djm.OneToOneField("Output", related_name="loss_curve")
hazard_output = djm.ForeignKey("Output", related_name="risk_loss_curves")
aggregate = djm.BooleanField(default=False)
insured = djm.BooleanField(default=False)
# If the curve is a result of an aggregation over different
# hazard_output the following fields must be set
statistics = djm.TextField(null=True, choices=STAT_CHOICES)
quantile = djm.FloatField(null=True)
loss_type = djm.TextField(choices=zip(LOSS_TYPES, LOSS_TYPES))
class Meta:
db_table = 'riskr\".\"loss_curve'
def __iter__(self):
if self.aggregate:
return iter([self.aggregatelosscurvedata])
else:
return iter(self.losscurvedata_set.all())
[docs] def to_csv_str(self):
"""
Convert LossCurve into a CSV string
"""
return '\n'.join(data.to_csv_str('row-%d' % i)
for i, data in enumerate(self, 1))
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
return (self.output.output_type,
self.output.hazard_metadata,
self.statistics, self.quantile,
self.aggregate, self.insured, self.loss_type)
[docs]class LossCurveData(djm.Model):
'''
Holds the probabilities of exceedance for a given loss curve
'''
loss_curve = djm.ForeignKey("LossCurve")
asset_ref = djm.TextField()
asset_value = djm.FloatField()
loss_ratios = fields.FloatArrayField()
poes = fields.FloatArrayField()
location = djm.PointField(srid=DEFAULT_SRID)
average_loss_ratio = djm.FloatField()
stddev_loss_ratio = djm.FloatField(blank=True, null=True)
class Meta:
db_table = 'riskr\".\"loss_curve_data'
@property
def losses(self):
return numpy.array(self.loss_ratios) * self.asset_value
@property
def average_loss(self):
return self.average_loss_ratio * self.asset_value
@property
def stddev_loss(self):
if self.stddev_loss_ratio is not None:
return self.stddev_loss_ratio * self.asset_value
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.loss_curve.output_hash + (self.asset_ref,)
[docs] def assertAlmostEqual(self, data):
return loss_curve_almost_equal(self, data)
[docs] def to_csv_str(self, label):
"""
Convert LossCurveData into a CSV string.
:param str label:
an identifier for the curve (for instance the asset_ref)
"""
ratios = [label, 'Ratios'] + map(str, self.loss_ratios)
data = ','.join(ratios) + '\n'
data += ','.join(map(str, [self.asset_value, 'PoE'] + list(self.poes)))
return data
[docs]class AggregateLossCurveData(djm.Model):
'''
Holds the probabilities of exceedance for the whole exposure model
'''
loss_curve = djm.OneToOneField("LossCurve")
losses = fields.FloatArrayField()
poes = fields.FloatArrayField()
average_loss = djm.FloatField()
stddev_loss = djm.FloatField(blank=True, null=True)
class Meta:
db_table = 'riskr\".\"aggregate_loss_curve_data'
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.loss_curve.output_hash
[docs] def assertAlmostEqual(self, data):
return loss_curve_almost_equal(self, data)
[docs] def to_csv_str(self, label):
"""
Convert AggregateLossCurveData into a CSV string.
:param str label:
an identifier for the curve (for instance the cost type)
"""
data = ','.join(map(str, [label, 'Losses'] + list(self.losses))) + '\n'
data += ','.join(map(str, ['', 'PoE'] + list(self.poes)))
return data
[docs]class EventLoss(djm.Model):
"""
Holds the aggregate loss we have for each rupture
"""
#: Foreign key to an :class:`openquake.engine.db.models.Output`
#: object with output_type == event_loss
output = djm.OneToOneField('Output', related_name="event_loss")
hazard_output = djm.ForeignKey(
"Output", related_name="risk_event_loss_tables")
loss_type = djm.TextField(choices=zip(LOSS_TYPES, LOSS_TYPES))
class Meta:
db_table = 'riskr\".\"event_loss'
def __iter__(self):
return iter(self.eventlossdata_set.all().order_by('-aggregate_loss'))
[docs] def to_csv_str(self):
"""
Convert EventLoss into a CSV with fields rupture_tag, aggregate_loss
"""
return '\n'.join('%s,%s' % (self.rupture.tag, self.aggregate_loss)
for data in self)
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
# FIXME(lp) this is not db-sequence independent
return (self.output.output_type, self.output.hazard_metadata,
self.loss_type)
[docs]class EventLossData(djm.Model):
event_loss = djm.ForeignKey(EventLoss)
rupture = djm.ForeignKey('SESRupture')
aggregate_loss = djm.FloatField()
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.event_loss.output_hash + (self.rupture_id,)
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(
self, data, operator.attrgetter('aggregate_loss'))
class Meta:
db_table = 'riskr\".\"event_loss_data'
[docs] def to_csv_str(self):
"""
Convert EventLossData into a CSV string
"""
return '%s,%s,%s' % (self.rupture.tag, self.rupture.rupture.mag,
self.aggregate_loss)
[docs]class EventLossAsset(djm.Model):
event_loss = djm.ForeignKey('EventLoss', null=False)
rupture = djm.ForeignKey('SESRupture', null=False)
asset = djm.ForeignKey('ExposureData', null=False)
loss = djm.FloatField(null=False)
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.event_loss.output_hash + (self.rupture_id, self.asset_id)
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(self, data, operator.attrgetter('loss'))
class Meta:
db_table = 'riskr\".\"event_loss_asset'
[docs] def to_csv_str(self):
"""
Convert EventLossAsset into a CSV string
"""
return '%s,%s,%s,%s' % (self.rupture.tag, self.rupture.rupture.mag,
self.asset.asset_ref, self.loss)
[docs]class BCRDistribution(djm.Model):
'''
Holds metadata for the benefit-cost ratio distribution
'''
output = djm.OneToOneField("Output", related_name="bcr_distribution")
hazard_output = djm.ForeignKey(
"Output", related_name="risk_bcr_distribution")
loss_type = djm.TextField(choices=zip(LOSS_TYPES, LOSS_TYPES))
class Meta:
db_table = 'riskr\".\"bcr_distribution'
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
return (self.output.output_type,
self.output.hazard_metadata,
self.loss_type)
[docs]class BCRDistributionData(djm.Model):
'''
Holds the actual data for the benefit-cost ratio distribution
'''
bcr_distribution = djm.ForeignKey("BCRDistribution")
asset_ref = djm.TextField()
average_annual_loss_original = djm.FloatField()
average_annual_loss_retrofitted = djm.FloatField()
bcr = djm.FloatField()
location = djm.PointField(srid=DEFAULT_SRID)
class Meta:
db_table = 'riskr\".\"bcr_distribution_data'
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.bcr_distribution.output_hash + (self.asset_ref,)
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(
self, data,
operator.attrgetter('average_annual_loss_original',
'average_loss_retrofitted',
'bcr'))
[docs]class DmgState(djm.Model):
"""Holds the damage_states associated to a given output"""
# they actually come from the fragility model xml input
risk_calculation = djm.ForeignKey("OqJob")
dmg_state = djm.TextField(
help_text="The name of the damage state")
lsi = djm.PositiveSmallIntegerField(
help_text="limit state index, to order the limit states")
class Meta:
db_table = 'riskr\".\"dmg_state'
[docs]class DmgDistPerAsset(djm.Model):
"""Holds the actual data for damage distributions per asset."""
dmg_state = djm.ForeignKey("DmgState")
exposure_data = djm.ForeignKey("ExposureData")
mean = djm.FloatField()
stddev = djm.FloatField()
class Meta:
db_table = 'riskr\".\"dmg_dist_per_asset'
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
return (self.output.output_type,
self.dmg_state.dmg_state, self.exposure_data.asset_ref)
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.output_hash
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(
self, data, operator.attrgetter('mean', 'stddev'))
@property
def output(self):
return self.dmg_state.rc_calculation.oqjob.output_set.get(
output_type="dmg_dist_per_asset")
[docs]class DmgDistPerTaxonomy(djm.Model):
"""Holds the actual data for damage distributions per taxonomy."""
dmg_state = djm.ForeignKey("DmgState")
taxonomy = djm.TextField()
mean = djm.FloatField()
stddev = djm.FloatField()
class Meta:
db_table = 'riskr\".\"dmg_dist_per_taxonomy'
@property
def output(self):
return self.dmg_state.rc_calculation.oqjob.output_set.get(
output_type="dmg_dist_per_taxonomy")
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
return (self.output.output_type,
self.dmg_state.dmg_state, self.taxonomy)
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.output_hash
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(
self, data, operator.attrgetter('mean', 'stddev'))
[docs]class DmgDistTotal(djm.Model):
"""Holds the actual 'total damage distribution' values for for an entire
calculation. There should be one record per calculation per damage state.
"""
dmg_state = djm.ForeignKey("DmgState")
mean = djm.FloatField()
stddev = djm.FloatField()
class Meta:
db_table = 'riskr\".\"dmg_dist_total'
@property
def output(self):
return self.dmg_state.rc_calculation.oqjob.output_set.get(
output_type="dmg_dist_total")
@property
def output_hash(self):
"""
:returns:
a (db-sequence independent) tuple that identifies this output among
which the ones created in the same calculation
"""
return (self.output.output_type, self.dmg_state.dmg_state, )
@property
def data_hash(self):
"""
A db-sequence independent tuple that identifies this output
"""
return self.output_hash
[docs] def assertAlmostEqual(self, data):
return risk_almost_equal(
self, data, operator.attrgetter('mean', 'stddev'))
## Tables in the 'riski' schema.
[docs]class ExposureModel(djm.Model):
'''
A risk exposure model
'''
job = djm.OneToOneField("OqJob")
name = djm.TextField()
description = djm.TextField(null=True)
category = djm.TextField()
taxonomy_source = djm.TextField(
null=True, help_text="the taxonomy system used to classify the assets")
AREA_CHOICES = (
(u'aggregated', u'Aggregated area value'),
(u'per_asset', u'Per asset area value'),
)
area_type = djm.TextField(null=True, choices=AREA_CHOICES)
area_unit = djm.TextField(null=True)
deductible_absolute = djm.BooleanField(default=True)
insurance_limit_absolute = djm.BooleanField(default=True)
class Meta:
db_table = 'riski\".\"exposure_model'
[docs] def taxonomies_in(self, region_constraint):
"""
:param str region_constraint:
polygon in wkt format the assets must be contained into
:returns:
A dictionary mapping each taxonomy with the number of assets
contained in `region_constraint`
"""
return ExposureData.objects.taxonomies_contained_in(
self.id, region_constraint)
[docs] def unit(self, loss_type):
if loss_type == "fatalities":
return "people"
else:
return self.costtype_set.get(name=loss_type).unit
[docs] def has_insurance_bounds(self):
return not self.exposuredata_set.filter(
(djm.Q(cost__deductible_absolute__isnull=True) |
djm.Q(cost__insurance_limit_absolute__isnull=True))).exists()
[docs] def has_retrofitted_costs(self):
return not (
self.exposuredata_set.filter(
cost__converted_retrofitted_cost__isnull=True)).exists()
[docs] def has_time_event(self, time_event):
return (
self.exposuredata_set.filter(occupancy__period=time_event).count()
==
self.exposuredata_set.count())
[docs] def supports_loss_type(self, loss_type):
"""
:returns:
True if the exposure contains the asset data needed
for computing losses for `loss_type`
"""
if loss_type != "fatalities":
ct = loss_type_to_cost_type(loss_type)
return self.exposuredata_set.filter(
cost__cost_type__name=ct).exists()
else:
if self.category == "population":
return not self.exposuredata_set.filter(
number_of_units__isnull=True).exists()
else:
return not self.exposuredata_set.filter(
occupancy__isnull=True).exists()
[docs]class CostType(djm.Model):
COST_TYPE_CHOICES = (
("structuralCost", "structuralCost"),
("retrofittedStructuralCost", "retrofittedStructuralCost"),
("nonStructuralCost", "nonStructuralCost"),
("contentsCost", "contentsCost"),
("businessInterruptionCost", "businessInterruptionCost"))
CONVERSION_CHOICES = (
(u'aggregated', u'Aggregated economic value'),
(u'per_area', u'Per area economic value'),
(u'per_asset', u'Per asset economic value'),
)
exposure_model = djm.ForeignKey(ExposureModel)
name = djm.TextField(choices=COST_TYPE_CHOICES)
conversion = djm.TextField(choices=CONVERSION_CHOICES)
unit = djm.TextField(null=True)
retrofitted_conversion = djm.TextField(
null=True, choices=CONVERSION_CHOICES)
retrofitted_unit = djm.TextField(null=True)
class Meta:
db_table = 'riski\".\"cost_type'
[docs]class Occupancy(djm.Model):
'''
Asset occupancy data
'''
exposure_data = djm.ForeignKey("ExposureData")
period = djm.TextField()
occupants = djm.FloatField()
class Meta:
db_table = 'riski\".\"occupancy'
[docs]class AssetManager(djm.GeoManager):
"""
Asset manager
"""
[docs] def get_asset_chunk(self, rc, assocs):
"""
:param assocs:
a list of :class:`openquake.engine.db.models.AssetSite` objects
:returns:
a list of instances of
:class:`openquake.engine.db.models.ExposureData` (ordered
by location) associated with
the `openquake.engine.db.models.ExposureModel` associated
with `rc`.
It also add an annotation to each ExposureData object to provide the
occupants value for the risk calculation given in input and the cost
for each cost type considered in `rc`
"""
assocs = sorted(assocs, key=lambda assoc: assoc.asset.id)
asset_ids = tuple(assoc.asset.id for assoc in assocs)
query, args = self._get_asset_chunk_query_args(rc, asset_ids)
with transaction.commit_on_success('job_init'):
annotated_assets = list(self.raw(query, args))
# add asset_site_id attribute to each asset
for ass, assoc in zip(annotated_assets, assocs):
ass.asset_site_id = assoc.id
return annotated_assets
def _get_asset_chunk_query_args(self, rc, asset_ids):
"""
Build a parametric query string and the corresponding args for
#get_asset_chunk
"""
args = (rc.exposure_model.id, asset_ids)
people_field, occupants_cond, occupancy_join, occupants_args = (
self._get_people_query_helper(
rc.exposure_model.category, rc.time_event))
args += occupants_args
cost_type_fields, cost_type_joins = self._get_cost_types_query_helper(
rc.exposure_model.costtype_set.all())
query = """\
SELECT riski.exposure_data.*,
{people_field} AS people,
{costs}
FROM riski.exposure_data
{occupancy_join}
ON riski.exposure_data.id = riski.occupancy.exposure_data_id
{costs_join}
WHERE exposure_model_id = %s
AND riski.exposure_data.id IN %s
AND {occupants_cond}
GROUP BY riski.exposure_data.id
ORDER BY riski.exposure_data.id
""".format(people_field=people_field,
occupants_cond=occupants_cond,
costs=cost_type_fields,
costs_join=cost_type_joins,
occupancy_join=occupancy_join)
return query, args
def _get_people_query_helper(self, category, time_event):
"""
Support function for _get_asset_chunk_query_args
"""
args = ()
# if the exposure model is of type "population" we extract the
# data from the `number_of_units` field
if category == "population":
occupants_field = "number_of_units"
occupants_cond = "1 = 1"
occupancy_join = ""
else:
# otherwise we will "left join" the occupancy table
occupancy_join = "LEFT JOIN riski.occupancy"
occupants_field = "AVG(riski.occupancy.occupants)"
# and the time_event is not specified we compute the
# number of occupants by averaging the occupancy data for
# each asset, otherwise we get the unique proper occupants
# value.
if time_event is None:
occupants_cond = "1 = 1"
else:
args += (time_event,)
occupants_cond = "riski.occupancy.period = %s"
return occupants_field, occupants_cond, occupancy_join, args
def _get_cost_types_query_helper(self, cost_types):
"""
Support function for _get_asset_chunk_query_args
"""
# For each cost type associated with the exposure model we
# join the `cost` table to the current queryset in order to
# lookup for a cost value for each asset.
# Actually we extract 4 values: the cost, the retrofitted
# cost, the deductible and the insurance limit
costs = []
costs_join = ""
for cost_type in cost_types:
# here the max value is irrelevant as we are sureto join
# against one row
costs.append("max(%s.converted_cost) AS %s" % (cost_type.name,
cost_type.name))
costs.append(
"max(%s.converted_retrofitted_cost) AS retrofitted_%s" % (
cost_type.name, cost_type.name))
costs.append(
"max(%s.deductible_absolute) AS deductible_%s" % (
cost_type.name, cost_type.name))
costs.append(
"max(%s.insurance_limit_absolute) AS insurance_limit_%s" % (
cost_type.name, cost_type.name))
costs_join += """
LEFT JOIN riski.cost AS %(name)s
ON %(name)s.cost_type_id = '%(id)s' AND
%(name)s.exposure_data_id = riski.exposure_data.id""" % dict(
name=cost_type.name, id=cost_type.id)
return ", ".join(costs), costs_join
[docs] def taxonomies_contained_in(self, exposure_model_id, region_constraint):
"""
:param exposure_model_id:
ID of an ExposureModel
:param region_constraint:
A string describing a region in WKT format
:returns:
A dictionary which map each taxonomy associated with
`exposure_model` and contained in `region_constraint` with the
number of assets.
"""
cursor = connections['job_init'].cursor()
cursor.execute("""
SELECT riski.exposure_data.taxonomy, COUNT(*)
FROM riski.exposure_data WHERE
exposure_model_id = %s AND ST_COVERS(ST_GeographyFromText(%s), site)
group by riski.exposure_data.taxonomy
""", [exposure_model_id, "SRID=4326; %s" % region_constraint])
return dict(cursor)
[docs]class ExposureData(djm.Model):
'''
Per-asset risk exposure data
'''
NO_RETROFITTING_COST = "no retrofitting cost"
exposure_model = djm.ForeignKey("ExposureModel")
asset_ref = djm.TextField()
taxonomy = djm.TextField()
site = djm.PointField(geography=True)
number_of_units = djm.FloatField(
null=True, help_text="number of assets, people, etc.")
area = djm.FloatField(null=True)
objects = AssetManager()
class Meta:
db_table = 'riski\".\"exposure_data'
ordering = ['asset_ref']
def __str__(self):
return "%s (%s-%s @ %s)" % (
self.id, self.exposure_model_id, self.asset_ref, self.site)
@staticmethod
[docs] def per_asset_value(
cost, cost_type, area, area_type, number_of_units, category):
"""
Return per-asset value for the given exposure data set.
Calculate per asset value by considering the given exposure
data as follows:
case 1: cost type: aggregated:
cost = economic value
case 2: cost type: per asset:
cost * number (of assets) = economic value
case 3: cost type: per area and area type: aggregated:
cost * area = economic value
case 4: cost type: per area and area type: per asset:
cost * area * number = economic value
The same "formula" applies to contenst/retrofitting cost analogously.
:returns:
The per-asset value as a `float`.
:raises:
`ValueError` in case of a malformed (risk exposure data) input.
"""
if category is not None and category == "population":
return number_of_units
if cost_type == "aggregated":
return cost
elif cost_type == "per_asset":
return cost * number_of_units
elif cost_type == "per_area":
if area_type == "aggregated":
return cost * area
elif area_type == "per_asset":
return cost * area * number_of_units
raise ValueError("Invalid input")
# we expect several annotations depending on "loss_type" to be
# present. See `get_asset_chunk` for details.
[docs] def value(self, loss_type):
"""
Extract the value of the asset for the given `loss_type`.
Although the Django Model definition does not have a value for
each loss type, we rely on the fact that an annotation on the
asset named `loss_type` is present.
"""
if loss_type == "fatalities":
return getattr(self, "people")
return getattr(self, loss_type)
[docs] def retrofitted(self, loss_type):
"""
Extract the retrofitted cost of the asset for the given
`loss_type`. See the method `value` for details.
"""
return getattr(self, "retrofitted_%s" % loss_type)
[docs] def deductible(self, loss_type):
"""
Extract the deductible limit of the asset for the given
`loss_type`. See the method `value` for details.
"""
return (getattr(self, "deductible_%s" % loss_type) /
getattr(self, loss_type))
[docs] def insurance_limit(self, loss_type):
"""
Extract the insurance limit of the asset for the given
`loss_type`. See the method `value` for details.
"""
return (getattr(self, "insurance_limit_%s" % loss_type) /
getattr(self, loss_type))
[docs]def make_absolute(limit, value, is_absolute=None):
"""
:returns:
`limit` if `is_absolute` is True or `limit` is None,
else `limit` * `value`
"""
if limit is not None:
if not is_absolute:
return value * limit
else:
return limit
[docs]class Cost(djm.Model):
exposure_data = djm.ForeignKey(ExposureData)
cost_type = djm.ForeignKey(CostType)
converted_cost = djm.FloatField()
converted_retrofitted_cost = djm.FloatField(null=True)
deductible_absolute = djm.FloatField(null=True, blank=True)
insurance_limit_absolute = djm.FloatField(null=True, blank=True)
class Meta:
db_table = 'riski\".\"cost'
## Tables in the 'htemp' schema.
[docs]class HazardSite(djm.Model):
"""
Contains pre-computed site parameter matrices. ``lons`` and ``lats``
represent the calculation sites of interest. The associated site parameters
are from the closest point in a site model in relation to each calculation
point of interest.
Used only if a calculation defines a site model (otherwise, reference
parameters are use for all points of interest).
"""
hazard_calculation = djm.ForeignKey('OqJob')
location = djm.PointField(srid=DEFAULT_SRID)
class Meta:
db_table = 'hzrdi\".\"hazard_site'
[docs]class AssetSite(djm.Model):
"""
Contains the association exposure_data_id -> site_id, as generated
by the current risk job.
"""
job = djm.ForeignKey('OqJob', null=False)
asset = djm.ForeignKey('ExposureData', null=False)
site = djm.ForeignKey('HazardSite', null=False)
class Meta:
db_table = 'riskr\".\"asset_site'
ordering = ['asset']
def __repr__(self):
return '<%s=%d, (%d, %d)>' % (self.__class__.__name__,
self.id, self.asset.id, self.site.id)
[docs]class Epsilon(djm.Model):
"""
Contains the association (asset_site, ses_collection) -> epsilons,
as generated by the current risk job for event based and scenario
computations.
"""
asset_site = djm.ForeignKey('AssetSite', null=False)
ses_collection = djm.ForeignKey('SESCollection', null=False)
epsilons = fields.FloatArrayField(null=False)
@classmethod
[docs] def saveall(cls, ses_coll, asset_sites, epsilon_matrix):
"""
Insert the epsilon matrix associated to the given
SES collection for each asset_sites association.
:param ses_coll:
a :class:`openquake.engine.db.models.SESCollection` instance
:param asset_sites:
a list of :class:`openquake.engine.db.models.AssetSite` instances
:param epsilon_matrix:
a numpy matrix with NxE elements, where `N` is the number of assets
and `E` the number of events for the given SESCollection
"""
assert len(asset_sites) == len(epsilon_matrix), (
len(asset_sites), len(epsilon_matrix))
data = [cls(asset_site=asset_site,
ses_collection=ses_coll,
epsilons=list(epsilons))
for asset_site, epsilons in zip(asset_sites, epsilon_matrix)]
return writer.CacheInserter.saveall(data)
class Meta:
db_table = 'riskr\".\"epsilon'
ordering = ['asset_site']
def __repr__(self):
return '<%s %r %r>' % (self.__class__.__name__, self.asset_site,
self.ses_collection)