Implements a new set of methods: Ground Motion Intensity Conversion Equations (GMICEs)

[g-weatherill]

Ground Motion Intensity Conversion Equations (GMICEs) refer to a class of empirical function that allows for conversion between conventional ground motion measures (e.g. PGA, PGV, Sa) and macroseismic intensity. Whilst not necessarily (as yet) integral to conventional PSHA analysis they are an important component of Shakemap(TM) applications of the oq-hazardlib. GMICEs work in both directions: from ground motion to macroseismic intensity and vice-versa. Both are supported here. Furthermore, we anticipate two possible use cases: 1) observed (or reported) values of ground motion and/or intensity are to be input to the GMICEs, 2) ground motions and/or intensities are predicted at a site using a conventional GMPE or IPE. Both cases are supported here too.

The implementation here is intended for subsequent integration into the Shakemap software (https://github.com/usgs/shakemap)

In a manner similar to GSIMs, the GMICEs are tested using custormised verification tables. A check_gmice test suite has also been added in a manner similar to that used by the GSIMs.

This PR implements the following:

1. Base class for GMICEs (openquake.hazardlib.gmice.base.GMICE)
2. Two GMICEs from Worden et al. (2012) - one case without magnitude and distance dependence (WordenEtAl2012) and one case with magnitude and distance dependence (WordenEtAl2012MagDist)
3. New test suite for GMICE verifications

GMICE Usage

In the following we demonstrate the usage of the GMICEs. Here we generate a set of ground motions based on a scenario, using the Cauzzi et al (2014) GMPE, and a set of macroseismic intensities from the same scenario using the Allen et al. (2012) IPE. We adopt the Worden et al. (2012) GMICE with magnitude and distance dependence

import numpy as np
from openquake.hazardlib import const
from openquake.hazardlib.imt import PGA, PGV, SA, MMI
from openquake.hazardlib.gsim.cauzzi_2014 import CauzziEtAl2014
from openquake.hazardlib.gsim.allen_2012_ipe import AllenEtAl2012
from openquake.hazardlib.gsim.base import (RuptureContext, DistancesContext, SitesContext)
from openquake.hazardlib.gmice.worden_2012 import (WordenEtAl2012, WordenEtAl2012MagDist)

Use the context objects to set up a rupture scenario (here Mw = 6.5, rrup from 1.0 km to 100 km, Vs30 (random between 180 and 800 m/s)

# Build Rupture
rctx = RuptureContext()
rctx.mag = 6.5
rctx.rake = 0.0
# Build Distances
dctx = DistancesContext()
dctx.rrup = np.arange(1.0, 101., 1.)
# Sites
sctx = SitesContext()
sctx.vs30 = np.random.uniform(180.0, 800.0, len(dctx.rrup))

ipe = AllenEtAl2012()
gmpe = CauzziEtAl2014()

# Generate intensity data set
mmi, _ = ipe.get_mean_and_stddevs(sctx, rctx, dctx, MMI(),[const.StdDev.TOTAL])

# Generate PGA data set
pga, _ = gmpe.get_mean_and_stddevs(sctx, rctx, dctx, PGA(),[const.StdDev.TOTAL])
pga = np.exp(pga) # PGA returned in ln(PGA) - get actual values

The GMICE has some additional metadata that can be accessed:

gmice = WordenEtAl2012MagDist()
# The range of supported intensities
print(gmice.INTENSITY_RANGE)

(1.0, 10.0)

# A 'pretty' version of the GMPE name
print(gmice.SOURCE_NAME)

'Worden et al. (2012) [Magnitude & Distance Dependent]'

# A 'scale' attribute - for Shakemap only
print(gmice.SCALE)

'scale_wgrw12.ps' We can convert from the PGA to MMI as follows:

# Note the intensity measure type IMT (here PGA()) needs to indicate the input ground motion intensity measure type
mmi_from_pga, sigma_mmi_from_pga = gmice.get_mean_intensity_and_stddevs(pga, sctx, rctx, dctx, PGA(), [const.StdDev.TOTAL])

Or, from the GMPE directly:

mmi_from_pga, sigma_mmi_from_pga = gmice.get_mean_intensity_and_stddevs_from_gmpe(CauzziEtAl2014(), sctx, rctx, dctx, PGA(), [const.StdDev.TOTAL])

Similarly, we can convert from MMI to PGA via:

# Note the intensity measure type IMT (here PGA()) needs to indicate the expected output intensity measure
pga_from_mmi, sigma_pga_from_mmi = gmice.get_mean_gm_and_stddevs(mmi, sctx, rctx, dctx, PGA(), [const.StdDev.TOTAL]

and directly via the IPE

# Note the intensity measure type IMT (here PGA()) needs to indicate the expected output intensity measure
pga_from_mmi, sigma_pga_from_mmi = gmice.get_mean_gm_and_stddevs_from_ipe(AllenEtAl2012(), sctx, rctx, dctx, PGA(), [const.StdDev.TOTAL]

The test verification tables are generated directly from the corresponding Shakemap class (https://github.com/usgs/shakemap/blob/master/shakemap/gmice/wgrw12.py). Note that, contrary to the original Shakemap implementation, accelerations are input fractions of g (g in our case equal to 9.80665 m/s/s). PGV is input in cm/s.

In addition to the OpenQuake team we encourage feedback from the shakemap developers (@cbworden, @emthompson-usgs, @mhearne-usgs)

[cbworden]

Some comments:

• The basic functionality looks good to me, and I think we can integrate it into the new ShakeMap without too much difficulty.
• Some GMICEs use more than one ground motion type to determine intensity. For instance, Wald et al (1999) has equations for both PGA and PGV, but recommends using PGA for lower intensities and PGV for higher intensities. For ShakeMap we've always followed that recommendation, and then smoothed the transition with a linear ramp. I'm not sure how to accommodate that functionality within the current framework.
• The metadata functions are helpful, but I wonder if it would be better to just have a get_metadata() function that returns a dictionary of metadata. That might make future modifications easier.
• Worden et al (2012) was derived with the IMC of maximum of 2 horizontal, not average horizontal.

[micheles]

The tests are running here: https://ci.openquake.org/job/zdevel_oq-hazardlib/650

[micheles]

LGTM from the IT point of view except the two minor comments before.