Source code for openquake.hazardlib.gsim.toro_1997

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

"""
Module exports :class:`ToroEtAl1997MblgNSHMP2008`,
:class:`ToroEtAl1997MwNSHMP2008`
"""
import numpy as np

from openquake.hazardlib.gsim.base import CoeffsTable, GMPE
from openquake.hazardlib.gsim.utils import (
    mblg_to_mw_johnston_96, mblg_to_mw_atkinson_boore_87, clip_mean)
from openquake.hazardlib import const
from openquake.hazardlib.imt import PGA, SA
from openquake.baselib.general import CallableDict

_compute_finite_fault_correction = CallableDict()


@_compute_finite_fault_correction.add("Mblg")
def _compute_finite_fault_correction_Mblg(kind, mag):
    """
    Compute finite fault correction term as geometric mean of correction
    terms obtained from Mw values calculated with Johnston 1996 and
    Atkinson and Boore 1987 conversion equations.

    Implement equations as in lines 1653 - 1658 in hazgridXnga2.f
    """
    mw_j96 = mblg_to_mw_johnston_96(mag)
    mw_ab87 = mblg_to_mw_atkinson_boore_87(mag)

    t1 = np.exp(-1.25 + 0.227 * mw_j96)
    t2 = np.exp(-1.25 + 0.227 * mw_ab87)

    return np.sqrt(t1 * t2)


@_compute_finite_fault_correction.add("Mw")
def _compute_finite_fault_correction_Mw(kind, mag):
    """
    Compute finite fault correction term.
    """
    return np.exp(-1.25 + 0.227 * mag)


def _compute_mean(kind, C, mag, rjb):
    """
    Compute ground motion mean value.
    """
    # line 1686 in hazgridXnga2.f
    ffc = _compute_finite_fault_correction(kind, mag)
    d = np.sqrt(rjb ** 2 + (C['c7'] ** 2) * (ffc ** 2))

    # lines 1663, 1694-1696 in hazgridXnga2.f
    mean = (C['c1'] + C['c2'] * (mag - 6.) +
            C['c3'] * ((mag - 6.) ** 2) -
            C['c4'] * np.log(d) - C['c6'] * d)

    factor = np.log(rjb / 100.)
    idx = factor > 0
    mean[idx] -= (C['c5'] - C['c4']) * factor[idx]

    return mean


[docs]class ToroEtAl1997MblgNSHMP2008(GMPE): """ Implements GMPE developed by G. R. Toro, N. A. Abrahamson, J. F. Schneider and published in "Model of Strong Ground Motions from Earthquakes in Central and Eastern North America: Best Estimates and Uncertainties" (Seismological Research Letters, Volume 68, Number 1, 1997) as utilized by the National Seismic Hazard Mapping Project (NSHMP) for the 2008 US hazard model. This class replicates the algorithm for the Toro et. al. 1997 GMPE as coded in the subroutine ``getToro`` in the ``hazgridXnga2.f`` Fortran code available at: http://earthquake.usgs.gov/hazards/products/conterminous/2008/software/ The class assumes rupture magnitude to be in Mblg scale (given that MFDs for central and eastern US are given in this scale). The equation implements also the finite-fault correction as given in "Modification of the Toro et al. 1997 Attenuation Equations for Large Magnitudes and Short Distances" (available at: http://www.riskeng.com/downloads/attenuation_equations). The correction uses Mw. Therefore Mblg is converted to Mw using both the Atkinson & Boore 1987 and Johnston 1996 conversion equations and an average correction term is computed. Coefficients are given for the B/C site conditions. """ kind = "Mblg" #: Supported tectonic region type is stable continental crust, #: given that the equations have been derived for central and eastern #: north America DEFINED_FOR_TECTONIC_REGION_TYPE = const.TRT.STABLE_CONTINENTAL #: Supported intensity measure types are spectral acceleration, #: and peak ground acceleration DEFINED_FOR_INTENSITY_MEASURE_TYPES = {PGA, SA} #: Supported intensity measure component is the geometric mean of #: two horizontal components #: :attr:`~openquake.hazardlib.const.IMC.AVERAGE_HORIZONTAL`, DEFINED_FOR_INTENSITY_MEASURE_COMPONENT = const.IMC.AVERAGE_HORIZONTAL #: Supported standard deviation type is only total. DEFINED_FOR_STANDARD_DEVIATION_TYPES = {const.StdDev.TOTAL} #: No site parameters required REQUIRES_SITES_PARAMETERS = set() #: Required rupture parameter is only magnitude (Mblg). REQUIRES_RUPTURE_PARAMETERS = {'mag'} #: Required distance measure is rjb REQUIRES_DISTANCES = {'rjb'} #: Shear-wave velocity for reference soil conditions in [m s-1] DEFINED_FOR_REFERENCE_VELOCITY = 760.
[docs] def compute(self, ctx, imts, mean, sig, tau, phi): """ See :meth:`superclass method <.base.GroundShakingIntensityModel.compute>` for spec of input and result values. """ for m, imt in enumerate(imts): C = self.COEFFS[imt] mean[m] = clip_mean( imt, _compute_mean(self.kind, C, ctx.mag, ctx.rjb)) sig[m] = C['sigma']
#: Coefficient table obtained from coefficient arrays (tb1, tb2, tb3, tb4, #: tb5, tb6, tbh) defined from line 1596 - 1614 in hazgridXnga2.f COEFFS = CoeffsTable(sa_damping=5, table="""\ IMT c1 c2 c3 c4 c5 c6 c7 sigma pga 2.489 1.20 0.0 1.28 1.23 0.0018 9.3 0.7506 0.1 2.91 1.23 0.0 1.12 1.05 0.0043 8.5 0.7506 0.2 2.165 1.24 0.0 0.98 0.74 0.0039 7.5 0.7506 0.3 1.7323 1.51 -0.11 0.96 0.6881 0.0034 7.35 0.7506 0.5 1.109 1.785 -0.2795 0.93 0.6354 0.002732 7.05 0.7506 1.0 0.173 2.05 -0.34 0.90 0.59 0.0019 6.8 0.799 2.0 -0.788 2.52 -0.47 0.93 0.6 0.0012 7.0 0.799 """)
[docs]class ToroEtAl1997MwNSHMP2008(ToroEtAl1997MblgNSHMP2008): """ Extend :class:`ToroEtAl1997MblgNSHMP2008` but assumes magnitude to be in Mw scale. Coefficients are Mw-specific and no magnitude conversion is considered to take into account finite-fault correction. """ kind = "Mw" #: Coefficient table obtained from coefficient arrays (tc1, tc2, tc3, tc4, #: tc5, tc6, th) defined in subroutine getToro in hazgridXnga2.f COEFFS = CoeffsTable(sa_damping=5, table="""\ IMT c1 c2 c3 c4 c5 c6 c7 sigma pga 2.619 0.81 0.0 1.27 1.16 0.0021 9.3 0.7506 0.1 2.92 0.81 0.0 1.1 1.02 0.004 8.3 0.7506 0.2 2.295 0.84 0.0 0.98 0.66 0.0042 7.5 0.7506 0.3 1.8823 0.964 -0.059 0.951 0.601 0.00367 7.26 0.7506 0.5 1.2887 1.14 -0.1244 0.9227 0.5429 0.00306 7.027 0.7506 1.0 0.383 1.42 -0.2 0.90 0.49 0.0023 6.8 0.799 2.0 -0.558 1.86 -0.31 0.92 0.46 0.0017 6.9 0.799 """)