Skip to content
This repository has been archived by the owner on Nov 10, 2017. It is now read-only.

Commit

Permalink
Merge pull request #400 from gem/better-filtering
Browse files Browse the repository at this point in the history
Move the rupture filtering inside the ContextMaker
  • Loading branch information
micheles committed Dec 4, 2015
2 parents 3bda4c5 + 9758b4d commit ade5f47
Show file tree
Hide file tree
Showing 5 changed files with 85 additions and 45 deletions.
1 change: 1 addition & 0 deletions debian/changelog
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
[Michele Simionato]
* Avoided doing the rupture filtering twice
* Optimized the case of multiple GSIMs by instantiating the
ruptures/sites/distances contexts only once

Expand Down
19 changes: 13 additions & 6 deletions openquake/hazardlib/calc/hazard_curve.py
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,7 @@
from openquake.baselib.python3compat import range, raise_
from openquake.baselib.performance import DummyMonitor
from openquake.hazardlib.calc import filters
from openquake.hazardlib.gsim.base import ContextMaker
from openquake.hazardlib.gsim.base import ContextMaker, FarAwayRupture
from openquake.hazardlib.imt import from_string
from openquake.baselib.general import deprecated

Expand Down Expand Up @@ -93,7 +93,8 @@ def hazard_curves(
def calc_hazard_curves(
sources, sites, imtls, gsim_by_trt, truncation_level=None,
source_site_filter=filters.source_site_noop_filter,
rupture_site_filter=filters.rupture_site_noop_filter):
rupture_site_filter=filters.rupture_site_noop_filter,
maximum_distance=None):
"""
Compute hazard curves on a list of sites, given a set of seismic sources
and a set of ground shaking intensity models (one per tectonic region type
Expand Down Expand Up @@ -158,11 +159,13 @@ def calc_hazard_curves(
return curves


# TODO: remove the rupture_site_filter, since its work is now done by the
# maximum_distance parameter; see what would break
def hazard_curves_per_trt(
sources, sites, imtls, gsims, truncation_level=None,
source_site_filter=filters.source_site_noop_filter,
rupture_site_filter=filters.rupture_site_noop_filter,
monitor=DummyMonitor()):
maximum_distance=None, monitor=DummyMonitor()):
"""
Compute the hazard curves for a set of sources belonging to the same
tectonic region type for all the GSIMs associated to that TRT.
Expand All @@ -175,7 +178,7 @@ def hazard_curves_per_trt(
by the intensity measure types; the size of each field is given by the
number of levels in ``imtls``.
"""
cmaker = ContextMaker(gsims)
cmaker = ContextMaker(gsims, maximum_distance)
gnames = list(map(str, gsims))
imt_dt = numpy.dtype([(imt, float, len(imtls[imt]))
for imt in sorted(imtls)])
Expand All @@ -194,15 +197,19 @@ def hazard_curves_per_trt(
(rupture, s_sites) for rupture in source.iter_ruptures()))
for rupture, r_sites in rupture_sites:
with ctx_mon:
sctx, rctx, dctx = cmaker.make_contexts(r_sites, rupture)
try:
sctx, rctx, dctx = cmaker.make_contexts(
r_sites, rupture)
except FarAwayRupture:
continue
for i, gsim in enumerate(gsims):
with pne_mon:
for imt in imts:
poes = gsim.get_poes(
sctx, rctx, dctx, imt, imts[imt],
truncation_level)
pno = rupture.get_probability_no_exceedance(poes)
expanded_pno = r_sites.expand(pno, placeholder=1)
expanded_pno = sctx.sites.expand(pno, 1.0)
curves[i][str(imt)] *= expanded_pno
except Exception as err:
etype, err, tb = sys.exc_info()
Expand Down
96 changes: 61 additions & 35 deletions openquake/hazardlib/gsim/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -108,21 +108,52 @@ def forbid_instantiation(cls):
cls.instantiable = True


def get_distances(rupture, mesh, param='rjb'):
"""
:param rupture: a rupture
:param mesh: a mesh of points
:param param: the kind of distance to compute (default rjb)
:returns: an array of distances from the given mesh
"""
if param == 'rrup':
dist = rupture.surface.get_min_distance(mesh)
elif param == 'rx':
dist = rupture.surface.get_rx_distance(mesh)
elif param == 'ry0':
dist = rupture.surface.get_ry0_distance(mesh)
elif param == 'rjb':
dist = rupture.surface.get_joyner_boore_distance(mesh)
elif param == 'rhypo':
dist = rupture.hypocenter.distance_to_mesh(mesh)
elif param == 'repi':
dist = rupture.hypocenter.distance_to_mesh(mesh, with_depths=False)
elif param == 'rcdpp':
dist = rupture.get_cdppvalue(mesh)
else:
raise ValueError('Unknown distance measure %r' % param)
return dist


class FarAwayRupture(Exception):
"""Raised if the rupture is outside the maximum distance for all sites"""


class ContextMaker(object):
"""
A class to manage the creation of contexts for distances, sites, rupture.
"""
REQUIRES = ['DISTANCES', 'SITES_PARAMETERS', 'RUPTURE_PARAMETERS']

def __init__(self, gsims):
def __init__(self, gsims, maximum_distance=None):
self.gsims = gsims
self.maximum_distance = maximum_distance
for req in self.REQUIRES:
reqset = set()
for gsim in gsims:
reqset.update(getattr(gsim, 'REQUIRES_' + req))
setattr(self, 'REQUIRES_' + req, reqset)

def make_distances_context(self, site_collection, rupture):
def make_distances_context(self, site_collection, rupture, dist_dict=()):
"""
Create distances context object for given site collection and rupture.
Expand All @@ -136,6 +167,9 @@ def make_distances_context(self, site_collection, rupture):
:class:
`~openquake.hazardlib.source.rupture.BaseProbabilisticRupture`).
:param dist_dict:
A dictionary of already computed distances, keyed by distance name
:returns:
Source to site distances as instance of :class:
`DistancesContext()`. Only those values that are required by GSIM
Expand All @@ -146,30 +180,11 @@ def make_distances_context(self, site_collection, rupture):
"""
dctx = DistancesContext()
for param in self.REQUIRES_DISTANCES:
if param == 'rrup':
dist = rupture.surface.get_min_distance(site_collection.mesh)
elif param == 'rx':
dist = rupture.surface.get_rx_distance(site_collection.mesh)
elif param == 'ry0':
dist = rupture.surface.get_ry0_distance(site_collection.mesh)
elif param == 'rjb':
dist = rupture.surface.get_joyner_boore_distance(
site_collection.mesh
)
elif param == 'rhypo':
dist = rupture.hypocenter.distance_to_mesh(
site_collection.mesh
)
elif param == 'repi':
dist = rupture.hypocenter.distance_to_mesh(
site_collection.mesh, with_depths=False
)
elif param == 'rcdpp':
dist = rupture.get_cdppvalue(site_collection.mesh)
if param in dist_dict: # already computed distances
distances = dist_dict[param]
else:
raise ValueError('%s requires unknown distance measure %r' %
(type(self).__name__, param))
setattr(dctx, param, dist)
distances = get_distances(rupture, site_collection.mesh, param)
setattr(dctx, param, distances)
return dctx

def make_sites_context(self, site_collection):
Expand All @@ -189,6 +204,7 @@ def make_sites_context(self, site_collection):
"""
sctx = SitesContext()
sctx.sites = site_collection
for param in self.REQUIRES_SITES_PARAMETERS:
try:
value = getattr(site_collection, param)
Expand All @@ -204,9 +220,8 @@ def make_rupture_context(self, rupture):
:param rupture:
Instance of
:class:`~openquake.hazardlib.source.rupture.Rupture` (or
subclass of
:class:`~openquake.hazardlib.source.rupture.BaseProbabilisticRupture`).
:class:`openquake.hazardlib.source.rupture.Rupture` or subclass of
:class:`openquake.hazardlib.source.rupture.BaseProbabilisticRupture`
:returns:
Rupture parameters as instance of :class:
Expand Down Expand Up @@ -244,16 +259,16 @@ def make_rupture_context(self, rupture):

def make_contexts(self, site_collection, rupture):
"""
Create context objects for given site collection and rupture.
Filter the site collection with respect to the rupture and
create context objects.
:param site_collection:
Instance of :class:`openquake.hazardlib.site.SiteCollection`.
:param rupture:
Instance of
:class:`~openquake.hazardlib.source.rupture.Rupture` (or
subclass of
:class:`~openquake.hazardlib.source.rupture.BaseProbabilisticRupture`).
:class:`openquake.hazardlib.source.rupture.Rupture` or subclass of
:class:`openquake.hazardlib.source.rupture.BaseProbabilisticRupture`
:returns:
Tuple of three items: sites context, rupture context and
Expand All @@ -267,9 +282,20 @@ def make_contexts(self, site_collection, rupture):
If any of declared required parameters (that includes site, rupture
and distance parameters) is unknown.
"""
return (self.make_sites_context(site_collection),
self.make_rupture_context(rupture),
self.make_distances_context(site_collection, rupture))
rctx = self.make_rupture_context(rupture)
distances = get_distances(rupture, site_collection.mesh, 'rjb')
sites = site_collection
if self.maximum_distance:
mask = distances <= self.maximum_distance
if mask.any():
sites = site_collection.filter(mask)
distances = distances[mask]
else:
raise FarAwayRupture

sctx = self.make_sites_context(sites)
dctx = self.make_distances_context(sites, rupture, {'rjb': distances})
return (sctx, rctx, dctx)


@functools.total_ordering
Expand Down
8 changes: 7 additions & 1 deletion openquake/hazardlib/tests/calc/hazard_curve_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,12 @@
from openquake.hazardlib.gsim.base import ContextMaker


class FakeSiteContext(object):
def __init__(self, sites):
self.sites = sites
self.mesh = sites.mesh


class HazardCurvesTestCase(unittest.TestCase):
class FakeRupture(object):
def __init__(self, probability, tectonic_region_type):
Expand Down Expand Up @@ -71,7 +77,7 @@ def __str__(self):
def setUp(self):
self.orig_make_contexts = ContextMaker.make_contexts
ContextMaker.make_contexts = lambda self, sites, rupture: (
sites, rupture, None)
FakeSiteContext(sites), rupture, None)
self.truncation_level = 3.4
self.imts = {'PGA': [1, 2, 3], 'PGD': [2, 4]}
self.time_span = 49.2
Expand Down
6 changes: 3 additions & 3 deletions openquake/hazardlib/tests/gsim/base_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -428,7 +428,7 @@ def make_contexts(self, site_collection, rupture):
def test_unknown_site_param_error(self):
self.gsim_class.REQUIRES_SITES_PARAMETERS.add('unknown!')
err = "ContextMaker requires unknown site parameter 'unknown!'"
sites = SiteCollection([self.site1])
sites = SiteCollection([self.site1, self.site2])
self._assert_value_error(self.make_contexts, err,
site_collection=sites, rupture=self.rupture)

Expand All @@ -441,8 +441,8 @@ def test_unknown_rupture_param_error(self):

def test_unknown_distance_error(self):
self.gsim_class.REQUIRES_DISTANCES.add('jump height')
err = "ContextMaker requires unknown distance measure 'jump height'"
sites = SiteCollection([self.site1])
err = "Unknown distance measure 'jump height'"
sites = SiteCollection([self.site1, self.site2])
self._assert_value_error(self.make_contexts, err,
site_collection=sites, rupture=self.rupture)

Expand Down

0 comments on commit ade5f47

Please sign in to comment.