# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4
#
# Copyright (C) 2015-2023 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:`AbrahamsonGulerce2020SInter`,
:class:`AbrahamsonGulerce2020SSlab`
"""
import numpy as np
from openquake.hazardlib.gsim.base import GMPE, CoeffsTable, add_alias
from openquake.hazardlib import const
from openquake.hazardlib.imt import PGA, SA
# The regions for which the model is supported. If not listed then the
# global model (GLO) should be applied
SUPPORTED_REGIONS = ["GLO", "USA-AK", "CAS", "CAM", "JPN", "NZL", "SAM", "TWN"]
# Region-specific constants or references to the corresponding column
# of the coefficients table
REGIONAL_TERMS = {
"GLO": {
"C1s": 7.5,
"a1": "a1",
},
"USA-AK": {
"C1s": 7.9, # In-slab magnitude scaling breakpoint
"a1": "a31", # Region-specific base constant
"a2-adj": None, # Adjustment to the geometric spreading term
"a6-adj": "a24", # Adjustment to the anelastic attenuation term
"a12-adj": "a17", # Adjustment to the linear Vs30 scaling term
"Z25-adj": None, # Adjustment to the basin depth scaling term
},
"CAS": {
"C1s": 7.1,
"a1": "a32",
"a2-adj": None,
"a6-adj": "a25",
"a12-adj": "a18",
"Z25-adj": "a39",
},
"CAM": {
"C1s": 7.4,
"a1": "a33",
"a2-adj": None,
"a6-adj": "a26",
"a12-adj": "a19",
"Z25-adj": None,
},
"JPN": {
"C1s": 7.6,
"a1": "a34",
"a2-adj": None,
"a6-adj": "a27",
"a12-adj": "a20",
"Z25-adj": "a41",
},
"NZL": {
"C1s": 8.0,
"a1": "a35",
"a2-adj": None,
"a6-adj": "a28",
"a12-adj": "a21",
"Z25-adj": None,
},
"SAM": {
"C1s": 7.5,
"a1": "a36",
"a2-adj": None,
"a6-adj": "a29",
"a12-adj": "a22",
"Z25-adj": None,
},
"TWN": {
"C1s": 7.7,
"a1": "a37",
"a2-adj": "a16",
"a6-adj": "a30",
"a12-adj": "a23",
"Z25-adj": None,
},
}
# Region- and period-independent constant values
CONSTS = {
"C4": 10.0,
"a3": 0.1,
"a4": 0.73,
"a5": 0.0,
"a9": 0.4,
"a15": -0.1,
"a17": 0.0,
"a19": 0.0,
"a45": 0.34,
"d0": 0.47,
"n": 1.18,
"c": 1.88,
"tau_lin": 0.47,
"phi_amp": 0.3,
"alpha_phi3": 0.42,
"T1_phi2": 0.03,
"T2_phi2": 0.075,
"T3_phi2": 0.20,
"T4_phi2": 1.0,
"T1_phi3": 0.03,
"T2_phi3": 0.075,
"T3_phi3": 0.10,
"T4_phi3": 0.3,
"d3_phi2": 0.109,
"d4_phi2": 0.062,
"d5_phi2": 0.470,
"d3_phi3": 0.242,
"d4_phi3": 0.000,
"d5_phi3": 0.000,
}
[docs]def get_base_term(C, region, apply_adjust):
"""
Returns the region-specific base term (a1 - Equation 3.1)
:param C:
Coefficient dictionary for the specfic IMT
:param str region:
Region identifier
:param bool apply_adjust:
For Alaska and Cascadia apply the modeller-defined adjustment factors
to the base term (True) or else leave regionalised but unadjusted
(False)
"""
if region in ("USA-AK", "CAS") and apply_adjust:
# For Alaska and Cascdia, apply the adjustments to the coefficient a1
return C[REGIONAL_TERMS[region]["a1"]] + C["{:s}_Adj".format(region)]
else:
return C[REGIONAL_TERMS[region]["a1"]]
[docs]def get_magnitude_scaling_term(C, trt, region, mag):
"""
Returns the magnitude scaling term (defined in Eq 3.3) and regional
constant
:param str trt:
Tectonic region type
:param np.ndarray mag:
Earthquake magnitude
"""
f_mag_quad = C["a13"] * ((10.0 - mag) ** 2.0)
a4 = CONSTS["a4"]
if trt == const.TRT.SUBDUCTION_INTERFACE:
c1 = C["c1i"]
else:
c1 = REGIONAL_TERMS[region]["C1s"]
a4 += CONSTS["a45"]
f_mag = np.where(mag <= c1,
f_mag_quad + a4 * (mag - c1),
f_mag_quad + CONSTS["a5"] * (mag - c1))
return f_mag
[docs]def get_geometric_spreading_term(C, region, mag, rrup):
"""
Returns the geometric spreading term defined in equation 3.1
:param numpy.ndarray rrup:
Rupture distances (km)
"""
a_2 = C["a2"]
if region == "TWN":
a_2 += C[REGIONAL_TERMS[region]["a2-adj"]]
f_r = a_2 + CONSTS["a3"] * (mag - 7.0)
hff = CONSTS["C4"] * np.exp(CONSTS["a9"] * (mag - 6.0))
f_r *= np.log(rrup + hff)
return f_r
[docs]def get_anelastic_attenuation_term(C, region, rrup):
"""
Returns the regionally-adjusted anelastic attenuation term
"""
a_6 = C["a6"]
if region != "GLO":
a_6 += C[REGIONAL_TERMS[region]["a6-adj"]]
return a_6 * rrup
[docs]def get_inslab_scaling_term(C, trt, region, mag, rrup):
"""
For inslab events, returns the inslab scaling term defined in equation 3.5
and corrected in the Erratum
"""
if trt == const.TRT.SUBDUCTION_INTERFACE:
# Doesn't apply to interface events
return 0.0
hff = 10.0 * np.exp(CONSTS["a9"] * (mag - 6.0))
c1s = REGIONAL_TERMS[region]["C1s"]
return C["a10"] + (CONSTS["a4"] + CONSTS["a45"]) * (c1s - 7.5) +\
C["a14"] * np.log(rrup + hff)
[docs]def get_rupture_depth_scaling_term(C, trt, ctx):
"""
Returns the rupture depth scaling described in Equation 3.6, which takes
the value 0 for interface events
:param numpy.ndarray ztor:
Top of rupture depths (km)
"""
if trt == const.TRT.SUBDUCTION_INTERFACE:
# Not defined for interface events
return 0.0
if isinstance(ctx.ztor, np.ndarray):
f_dep = C["a11"] * (ctx.ztor - 50.0)
f_dep[ctx.ztor > 200.0] = C["a11"] * 150.0
idx = ctx.ztor <= 50.0
f_dep[idx] = C["a8"] * (ctx.ztor[idx] - 50.0)
elif ctx.ztor > 200.0:
f_dep = C["a11"] * 150.0
elif ctx.ztor <= 50.0:
f_dep = C["a8"] * (ctx.ztor - 50.0)
else:
f_dep = C["a11"] * (ctx.ztor - 50.0)
return f_dep
[docs]def get_site_amplification_term(C, region, vs30, pga1000):
"""
Returns the shallow site amplification term as descrbied in Equation 3.7,
and corrected in the Erratum
:param numpy.ndarray vs30:
30-m averaged shearwave velocity (m/s)
:param numpy.ndarray pga1000:
Peak Ground Acceleration (PGA), g, on a reference bedrock of 1000 m/s
"""
a12 = C["a12"]
if region != "GLO":
# Apply regional adjustment
a12 += C[REGIONAL_TERMS[region]["a12-adj"]]
# V* is defined according to Equation 3.8 (i.e. Vs30 clipped at 1000 m/s
# in the Erratum)
vstar = np.clip(vs30, -np.inf, 1000.0)
vnorm = vstar / C["vlin"]
fsite = a12 * np.log(vnorm)
idx = vstar >= C["vlin"]
# Linear site term
fsite[idx] += (C["b"] * CONSTS["n"] * np.log(vnorm[idx]))
idx = np.logical_not(idx)
# Nonlinear site term
fsite[idx] += C["b"] * (
np.log(pga1000[idx] + CONSTS["c"] * (vnorm[idx] ** CONSTS["n"])) -
np.log(pga1000[idx] + CONSTS["c"])
)
return fsite
[docs]def get_reference_basin_depth(region, vs30):
"""
For the Cascadia and Japan regions a reference basin depth, dependent
on the Vs30, is returned according to equations 2.1 and 2.2
"""
if region == "CAS":
ln_zref = np.clip(8.52 - 0.88 * np.log(vs30 / 200.0), 7.6, 8.52)
elif region == "JPN":
ln_zref = np.clip(7.3 - 2.066 * np.log(vs30 / 170.0), 4.1, 7.3)
else:
raise ValueError("No reference basin depth term defined for region %s"
% region)
return np.exp(ln_zref)
[docs]def get_basin_depth_scaling(C, region, vs30, z25):
"""
Returns the basin depth scaling term, applicable only for the Cascadia
and Japan regions, defined in equations 3.9 - 3.11 and corrected in the
Erratum
:param numpy.ndarray z25:
Depth to 2.5 m/s shearwave velocity layer (km)
"""
if region not in ("CAS", "JPN"):
# Basin depth defined only for Cascadia and Japan, so return 0
return 0.0
# Define the reference basin depth from Vs30
z25_ref = get_reference_basin_depth(region, vs30)
# Normalise the basin depth term (Equation 3.9)
ln_z25_prime = np.log((z25 + 50.0) / (z25_ref + 50.0))
f_basin = np.zeros(vs30.shape)
if region == "JPN":
# Japan Basin (Equation 3.10)
f_basin = C["a41"] * np.clip(ln_z25_prime, -2.0, np.inf)
else:
# Cascadia Basin (Equation 3.11)
idx = ln_z25_prime > 0.0
f_basin[idx] = C["a39"] * ln_z25_prime[idx]
return f_basin
[docs]def get_acceleration_on_reference_rock(C, trt, region, ctx, apply_adjustment):
"""
Returns acceleration on reference rock - intended for use primarily with
PGA. Overrides the Vs30 values and removes any basin depth terms
"""
# Set all Vs30 to 1000 m/s
vs30 = np.full_like(ctx.vs30, 1000.)
# On rock the amplification is only linear, so PGA1000 is not used
# set to a null array
null_pga1000 = np.zeros(vs30.shape)
# No basin depth is calculated here
return (get_base_term(C, region, apply_adjustment) +
get_magnitude_scaling_term(C, trt, region, ctx.mag) +
get_geometric_spreading_term(C, region, ctx.mag, ctx.rrup) +
get_anelastic_attenuation_term(C, region, ctx.rrup) +
get_rupture_depth_scaling_term(C, trt, ctx) +
get_inslab_scaling_term(C, trt, region, ctx.mag, ctx.rrup) +
get_site_amplification_term(C, region, vs30, null_pga1000))
[docs]def get_mean_acceleration(C, trt, region, ctx, pga1000, apply_adjustment):
"""
Returns the mean acceleration on soil
"""
if region in ("CAS", "JPN"):
# Convert z25 from km to m
z25 = 1000.0 * ctx.z2pt5
else:
# Basin depths will be ignored, so set zeros
z25 = np.zeros(ctx.vs30.shape)
return (get_base_term(C, region, apply_adjustment) +
get_magnitude_scaling_term(C, trt, region, ctx.mag) +
get_geometric_spreading_term(C, region, ctx.mag, ctx.rrup) +
get_anelastic_attenuation_term(C, region, ctx.rrup) +
get_rupture_depth_scaling_term(C, trt, ctx) +
get_inslab_scaling_term(C, trt, region, ctx.mag, ctx.rrup) +
get_site_amplification_term(C, region, ctx.vs30, pga1000) +
get_basin_depth_scaling(C, region, ctx.vs30, z25))
def _get_f2(t1, t2, t3, t4, alpha, period):
"""
Returns the linear phi short-period scaling term (f2) defined in
equation 5.3 and corrected in the Erratum
"""
if period <= t1:
f2 = 1.0 - alpha
elif period <= t2:
f2 = 1.0 - alpha * (np.log(period / t2) / np.log(t1 / t2))
elif period <= t3:
f2 = 1.0
elif period < t4:
f2 = np.log(period / t4) / np.log(t3 / t4)
else:
f2 = 0.0
return f2
def _get_a_phi2(rrup, d3, d4, d5):
"""
Returns Aphi, the short period linear phi normalising factor to be
added to the base linear within-event term, as defined in Equation 5.4
"""
rrup_norm = (rrup - 225.0) / 225.0
a_phi2 = d3 + d4 * rrup_norm + d5 * (rrup_norm ** 2.)
a_phi2[rrup < 225] = d3
a_phi2[rrup > 450.0] = d3 + d4 + d5
return a_phi2
[docs]def get_phi_lin_model(C, C_PGA, region, period, rrup):
"""
Returns the distance-dependent linear phi term for both PGA and the
required spectral period. The term is regionally dependent with additional
factors added on for Central America, Japan and South America
Several equations are used here, described fully in section 5.3
:param float period:
Spectral period of ground motion
"""
# Define phi1 term according to equation 5.2
rrup_norm = (rrup - 150.0) / 300.0
phi1_sq = np.clip(C["d1"] + C["d2"] * rrup_norm,
C["d1"],
C["d1"] + C["d2"])
phi1_sq_pga = np.clip(C_PGA["d1"] + C_PGA["d2"] * rrup_norm,
C_PGA["d1"],
C_PGA["d1"] + C_PGA["d2"])
if region in ("USA-AK", "CAS", "NZL", "TWN", "GLO"):
# For Alaska, Cascadia, New Zealand and Taiwan then phi_lin corresponds
# only to phi1
return np.sqrt(phi1_sq), np.sqrt(phi1_sq_pga)
# Phi 2 model
# Alpha according to equation 5.6
alpha = np.clip(1.0 - 0.0036 * (rrup - 250.0), 0.28, 1.0)
# Retrieve the specific corner periods needed for phi2
t1, t2, t3, t4 = CONSTS["T1_phi2"], CONSTS["T2_phi2"], CONSTS["T3_phi2"],\
CONSTS["T4_phi2"]
# Retrieve the function describing the trapezoidal shape of the increase
# in linear phi at short periods
f2_phi2 = _get_f2(t1, t2, t3, t4, alpha, period)
f2_phi2_pga = 1.0 - alpha
# Retrieve constants needed for the amplitude scaling factor for
# the phi2 increase
d3, d4, d5 = CONSTS["d3_phi2"], CONSTS["d4_phi2"], CONSTS["d5_phi2"]
a_phi2 = _get_a_phi2(rrup, d3, d4, d5)
phi2_sq = a_phi2 * f2_phi2
phi2_sq_pga = a_phi2 * f2_phi2_pga
if region == "CAM":
# For Central America and Mexico only the phi1 and phi2 terms are used
return np.sqrt(phi1_sq + phi2_sq), np.sqrt(phi1_sq_pga + phi2_sq_pga)
# Phi 3 model - applies to Japan and South America
# Alpha term is constant
alpha = CONSTS["alpha_phi3"]
# Retrieve normalised phi_lin scaling function
t1, t2, t3, t4 = CONSTS["T1_phi3"], CONSTS["T2_phi3"], CONSTS["T3_phi3"],\
CONSTS["T4_phi3"]
f2_phi3 = _get_f2(t1, t2, t3, t4, alpha, period)
f2_phi3_pga = 1.0 - alpha
# In this case the amplitude of additional variance is constant
a_phi3 = CONSTS["d3_phi3"]
phi3_sq = a_phi3 * f2_phi3
phi3_sq_pga = a_phi3 * f2_phi3_pga
return np.sqrt(phi1_sq + phi2_sq + phi3_sq),\
np.sqrt(phi1_sq_pga + phi2_sq_pga + phi3_sq_pga)
[docs]def get_partial_derivative_site_pga(C, vs30, pga1000):
"""
Defines the partial derivative of the site term with respect to the PGA
on reference rock, described in equation 5.9 (corrected in Erratum)
"""
dfsite_dlnpga = np.zeros(vs30.shape)
vstar = np.clip(vs30, -np.inf, 1000.0)
idx = vstar <= C["vlin"]
vnorm = vstar[idx] / C["vlin"]
dfsite_dlnpga[idx] = C["b"] * pga1000 * (
(-1.0 / (pga1000 + CONSTS["c"])) +
(1.0 / (pga1000 + CONSTS["c"] * (vnorm ** CONSTS["n"])))
)
return dfsite_dlnpga
[docs]def get_tau_phi(C, C_PGA, region, period, rrup, vs30, pga1000, ergodic):
"""
Get the heteroskedastic within-event and between-event standard deviation
"""
# Get the site-to-site variability
if region == "CAM":
phi_s2s = C["phi_s2s_g2"]
phi_s2s_pga = C_PGA["phi_s2s_g2"]
elif region in ("JPN", "SAM"):
phi_s2s = C["phi_s2s_g3"]
phi_s2s_pga = C_PGA["phi_s2s_g3"]
else:
phi_s2s = C["phi_s2s_g1"]
phi_s2s_pga = C["phi_s2s_g1"]
# Get linear tau and phi
# linear tau is period independent
tau = CONSTS["tau_lin"] * np.ones(vs30.shape)
phi_lin, phi_lin_pga = get_phi_lin_model(C, C_PGA, region, period, rrup)
# Find the sites where nonlinear site terms apply
idx = np.clip(vs30, -np.inf, 1000) < C["vlin"]
if not np.any(idx):
# Only linear term
if ergodic:
return tau, phi_lin
else:
# Remove the site-to-site variability from phi
phi = np.sqrt(phi_lin ** 2.0 - phi_s2s ** 2.0)
return tau, phi
# Process the nonlinear site terms
phi = phi_lin.copy()
partial_f_pga = get_partial_derivative_site_pga(C, vs30[idx], pga1000[idx])
phi_b = np.sqrt(phi_lin[idx] ** 2.0 - CONSTS["phi_amp"] ** 2.0)
phi_b_pga = np.sqrt(phi_lin_pga[idx] ** 2.0 - CONSTS["phi_amp"] ** 2.0)
# Get nonlinear tau and phi terms
tau_sq = tau[idx] ** 2.0
tau_sq = tau_sq + (partial_f_pga ** 2.0) * tau_sq +\
2.0 * partial_f_pga * tau_sq * C["rhoB"]
phi_nl_sq = (phi_lin[idx] ** 2.0) +\
(partial_f_pga ** 2.0) * (phi_b_pga ** 2.0) +\
(2.0 * partial_f_pga * phi_b_pga * phi_b * C["rhoW"])
tau[idx] = np.sqrt(tau_sq)
phi[idx] = np.sqrt(phi_nl_sq)
if not ergodic:
# Need to update the nonlinear within-event variability term to remove
# the site-to-site variability
phi_ss = np.sqrt(phi_lin ** 2.0 - phi_s2s ** 2.0)
phi_ss_pga = np.sqrt(phi_lin_pga ** 2.0 - phi_s2s_pga ** 2.0)
phi_ss_b = np.sqrt(phi_ss[idx] ** 2.0 - CONSTS["phi_amp"] ** 2.0)
phi_ss_b_pga = np.sqrt(phi_ss_pga[idx] ** 2.0 -
CONSTS["phi_amp"] ** 2.0)
phi_ss_nl_sq = (phi_ss[idx] ** 2.0) +\
(partial_f_pga ** 2.0) * (phi_ss_b_pga ** 2.0) +\
(2.0 * partial_f_pga * phi_ss_b_pga * phi_ss_b * C["rhoW"])
phi_ss[idx] = np.sqrt(phi_ss_nl_sq)
return tau, phi_ss
return tau, phi
[docs]def get_epistemic_adjustment(C, rrup):
"""
Returns the distance-dependent epistemic adjustment factor defined in
equation 6.1. In theory, this should only be applied to the global
model, but we do not enforce this constraint here.
"""
rrup_norm = np.clip(rrup, 50.0, 500.0) / 100.0
return C["e1"] + C["e2"] * rrup_norm + C["e3"] * (rrup_norm ** 2.0)
[docs]class AbrahamsonGulerce2020SInter(GMPE):
"""
Implements the 2020 Subduction ground motion model of Abrahamson &
Gulerce (2020):
Abrahamson N. and Gulurce Z. (2020) "Regionalized Ground-Motion Models
for Subduction Earthquakes based on the NGA-SUB Database", Pacific
Earthquake Engineering Research Center (PEER) Technical Report,
PEER 2020/25
The model is regionalised, defining specific adjustment factors for
(invoking region term in parenthesis):
- Global ("GLO" - for application to any subduction region for which
no region-specific adjustment is defined)
- Alaska ("USA-AK")
- Cascadia ("CAS")
- Central America & Mexico ("CAM")
- Japan ("JPN")
- New Zealand ("NZL")
- South America ("SAM")
- Taiwan ("TWN")
The region-specific adjustments primarily affect the constant term,
the anelastic attenuation term and the linear Vs30 scaling term. In
addition, however, further period-specific adjustment factors can be
applied for the Alaska and Cascadia regions using the boolean
input `apply_adjustment`. These adjustments scale the resulting
ground motion values to appropriate levels accounting for limited data
and the Alaska and Cascadia region, based on analysis undertaken by
the authors.
A general epistemic uncertainty median adjustment factor is also defined
based on the standard deviation of the median ground motion from five
regions with estimated regional terms. This term should be applied only
to the global model (though this is not strictly enforced), and it is
controlled via the use of `sigma_mu_epsilon`, the number of standard
deviations by which the adjustment will be multiplied (default = 0)
A non-ergodic aleatory uncertainty model can be returned by setting
`ergodic=False`.
The code implementation and test tables have been verified using
Fortran code supplied by Professor N. Abrahamson, and cross-checked against
an independent implementation from Feng Li, Jason Motha and James Paterson
from University of Canterbury (New Zealand).
Attributes:
region (str): Choice of region among the supported regions ("GLO",
"USA-AK", "CAS", "CAM", "JPN", "NZL", "SAM", "TWN")
ergodic (bool): Return the ergodic aleatory variability model (True)
or non-ergodic form (False)
apply_usa_adjustment (bool): Apply the modeller designated Alaska or
Cascadia adjustments (available only for
the regions "USA-AK" or "CAS")
sigma_mu_epsilon (float): Number of standard deviations to multiply
sigma mu (the standard deviation of the
median) for the epistemic uncertainty model
"""
#: Supported tectonic region type is subduction interface
DEFINED_FOR_TECTONIC_REGION_TYPE = const.TRT.SUBDUCTION_INTERFACE
#: Supported intensity measure types are spectral acceleration,
#: and peak ground acceleration
DEFINED_FOR_INTENSITY_MEASURE_TYPES = {PGA, SA}
#: Supported intensity measure component is RotD50
DEFINED_FOR_INTENSITY_MEASURE_COMPONENT = const.IMC.RotD50
#: Supported standard deviation types are inter-event, intra-event
#: and total, see section 4.5
DEFINED_FOR_STANDARD_DEVIATION_TYPES = {
const.StdDev.TOTAL, const.StdDev.INTER_EVENT, const.StdDev.INTRA_EVENT}
#: Site amplification is dependent only upon Vs30 for the majority of cases
#: but Z2.5 is added for the JPN and CAS regions
REQUIRES_SITES_PARAMETERS = {'vs30', }
#: Required rupture parameters are only magnitude for the interface model
REQUIRES_RUPTURE_PARAMETERS = {'mag', }
#: Required distance measure is closest distance to rupture, for
#: interface events
REQUIRES_DISTANCES = {'rrup', }
#: Defined for a reference velocity of 1000 m/s
DEFINED_FOR_REFERENCE_VELOCITY = 1000.0
def __init__(self, region="GLO", ergodic=True, apply_usa_adjustment=False,
sigma_mu_epsilon=0.0):
assert region in SUPPORTED_REGIONS, "Region %s not supported by %s" \
% (region, self.__class__.__name__)
self.region = region
self.ergodic = ergodic
self.apply_usa_adjustment = apply_usa_adjustment
self.sigma_mu_epsilon = sigma_mu_epsilon
# If running for Cascadia or Japan then z2.5 is needed
if region in ("CAS", "JPN"):
self.REQUIRES_SITES_PARAMETERS = \
self.REQUIRES_SITES_PARAMETERS.union({"z2pt5", })
[docs] def compute(self, ctx: np.recarray, imts, mean, sig, tau, phi):
"""
See :meth:`superclass method
<.base.GroundShakingIntensityModel.compute>`
for spec of input and result values.
"""
trt = self.DEFINED_FOR_TECTONIC_REGION_TYPE
C_PGA = self.COEFFS[PGA()]
pga1000 = get_acceleration_on_reference_rock(C_PGA, trt,
self.region, ctx,
self.apply_usa_adjustment)
pga1000 = np.exp(pga1000)
for m, imt in enumerate(imts):
C = self.COEFFS[imt]
mean[m] = get_mean_acceleration(C, trt, self.region, ctx, pga1000,
self.apply_usa_adjustment)
if self.sigma_mu_epsilon:
# Apply an epistmic adjustment factor
mean[m] += (self.sigma_mu_epsilon *
get_epistemic_adjustment(C, ctx.rrup))
# Get the standard deviations
tau_m, phi_m = get_tau_phi(C, C_PGA, self.region, imt.period,
ctx.rrup, ctx.vs30, pga1000,
self.ergodic)
tau[m] = tau_m
phi[m] = phi_m
sig += np.sqrt(tau ** 2.0 + phi ** 2.0)
# Coefficients taken from digital files supplied by Norm Abrahamson
COEFFS = CoeffsTable(sa_damping=5, table="""\
imt c1i vlin b a1 a2 a6 a7 a8 a10 a11 a12 a13 a14 a16 a17 a18 a19 a20 a21 a22 a23 a24 a25 a26 a27 a28 a29 a30 a31 a32 a33 a34 a35 a36 a37 a39 a41 USA-AK_Adj CAS_Adj d1 d2 rhoW rhoB phi_s2s_g1 phi_s2s_g2 phi_s2s_g3 e1 e2 e3
pga 8.20 865.1 -1.186 4.5960 -1.4500 -0.0043 3.2100 0.0440 3.210 0.0070 0.9000 0.0000 -0.4600 0.0900 0.000 -0.2000 0.000 0.0000 0.0400 0.0400 0.0000 0.0015 0.0007 0.0036 -0.0004 0.0025 0.0006 0.0033 3.7783 3.3468 3.8025 5.0361 4.6272 4.8044 3.5669 0.000 -0.029 0.487 0.828 0.325 0.137 1.000 1.000 0.396 0.396 0.545 0.550 -0.270 0.050
0.010 8.20 865.1 -1.186 4.5960 -1.4500 -0.0043 3.2100 0.0440 3.210 0.0070 0.9000 0.0000 -0.4600 0.0900 0.000 -0.2000 0.000 0.0000 0.0400 0.0400 0.0000 0.0015 0.0007 0.0036 -0.0004 0.0025 0.0006 0.0033 3.7783 3.3468 3.8025 5.0361 4.6272 4.8044 3.5669 0.000 -0.029 0.487 0.828 0.325 0.137 1.000 1.000 0.396 0.396 0.545 0.550 -0.270 0.050
0.020 8.20 865.1 -1.219 4.6780 -1.4500 -0.0043 3.2100 0.0440 3.210 0.0070 1.0080 0.0000 -0.4600 0.0900 0.000 -0.2000 0.000 0.0000 0.0400 0.0400 0.0000 0.0015 0.0006 0.0036 -0.0005 0.0025 0.0005 0.0033 3.8281 3.4401 3.9053 5.1375 4.6958 4.8943 3.6425 0.000 -0.024 0.519 0.825 0.325 0.137 0.990 0.990 0.396 0.396 0.545 0.550 -0.270 0.050
0.030 8.20 907.8 -1.273 4.7730 -1.4500 -0.0044 3.2100 0.0440 3.210 0.0070 1.1270 0.0000 -0.4600 0.0900 0.000 -0.2000 0.000 0.0000 0.0400 0.0400 0.0000 0.0015 0.0006 0.0037 -0.0007 0.0025 0.0005 0.0034 3.8933 3.5087 4.0189 5.2699 4.7809 5.0028 3.7063 0.000 -0.034 0.543 0.834 0.325 0.137 0.990 0.990 0.396 0.396 0.545 0.550 -0.270 0.050
0.050 8.20 1053.5 -1.346 5.0290 -1.4500 -0.0046 3.2100 0.0440 3.210 0.0070 1.3330 0.0000 -0.4600 0.0900 0.000 -0.2000 0.000 0.0000 0.0400 0.0400 0.0000 0.0011 0.0006 0.0039 -0.0009 0.0026 0.0004 0.0036 4.2867 3.6553 4.2952 5.6157 5.0211 5.2819 3.9184 0.000 -0.061 0.435 0.895 0.325 0.137 0.970 0.985 0.396 0.467 0.644 0.560 -0.270 0.050
0.075 8.20 1085.7 -1.471 5.3340 -1.4500 -0.0047 3.2100 0.0440 3.210 0.0070 1.5650 0.0000 -0.4600 0.0900 0.000 -0.2000 0.000 0.0000 0.0600 0.0600 0.0000 0.0011 0.0004 0.0039 -0.0009 0.0026 0.0003 0.0037 4.5940 3.9799 4.5464 6.0204 5.3474 5.6123 4.2207 0.000 -0.076 0.410 0.863 0.325 0.137 0.950 0.980 0.396 0.516 0.713 0.580 -0.270 0.050
0.100 8.20 1032.5 -1.624 5.4550 -1.4500 -0.0048 3.2100 0.0440 3.210 0.0070 1.6790 0.0000 -0.4600 0.0900 0.000 -0.2000 0.000 0.0000 0.1000 0.1000 0.0000 0.0012 0.0003 0.0039 -0.0008 0.0026 0.0003 0.0038 4.7077 4.1312 4.6138 6.1625 5.5065 5.7668 4.3536 0.000 -0.049 0.397 0.842 0.325 0.137 0.920 0.970 0.396 0.516 0.713 0.590 -0.270 0.050
0.150 8.20 877.6 -1.931 5.3760 -1.4250 -0.0047 3.2100 0.0440 3.210 0.0070 1.8530 0.0000 -0.4600 0.0900 0.000 -0.1860 0.000 -0.0550 0.1350 0.1350 0.0690 0.0013 -0.0002 0.0037 -0.0009 0.0022 0.0001 0.0037 4.6065 4.2737 4.5290 5.9614 5.5180 5.7313 4.3664 0.000 -0.026 0.428 0.737 0.325 0.137 0.900 0.960 0.396 0.516 0.647 0.590 -0.270 0.050
0.200 8.20 748.2 -2.188 4.9360 -1.3350 -0.0045 3.2100 0.0430 3.210 0.0062 2.0220 0.0000 -0.4600 0.0840 0.000 -0.1500 0.000 -0.1050 0.1700 0.1700 0.1400 0.0013 -0.0007 0.0031 -0.0010 0.0018 -0.0001 0.0035 4.1866 3.9650 4.1656 5.3920 5.1668 5.2943 4.0169 0.000 -0.011 0.442 0.746 0.325 0.137 0.870 0.940 0.396 0.516 0.596 0.570 -0.270 0.050
0.250 8.20 654.3 -2.381 4.6360 -1.2750 -0.0043 3.2100 0.0420 3.210 0.0056 2.1810 0.0000 -0.4600 0.0800 0.000 -0.1400 0.000 -0.1340 0.1700 0.1700 0.1640 0.0013 -0.0009 0.0027 -0.0011 0.0016 -0.0003 0.0033 3.8515 3.6821 3.9147 5.0117 4.8744 5.0058 3.7590 0.101 -0.009 0.494 0.796 0.325 0.137 0.840 0.930 0.396 0.501 0.539 0.530 -0.224 0.043
0.300 8.20 587.1 -2.518 4.4230 -1.2310 -0.0042 3.2100 0.0410 3.210 0.0051 2.2810 -0.0020 -0.4600 0.0780 0.000 -0.1200 0.000 -0.1500 0.1700 0.1700 0.1900 0.0014 -0.0010 0.0020 -0.0009 0.0014 -0.0002 0.0032 3.5783 3.5415 3.7846 4.7057 4.6544 4.7588 3.5914 0.184 0.005 0.565 0.782 0.325 0.137 0.820 0.910 0.396 0.488 0.488 0.490 -0.186 0.037
0.400 8.20 503.0 -2.657 4.1240 -1.1650 -0.0040 3.2100 0.0400 3.210 0.0043 2.3790 -0.0070 -0.4700 0.0750 0.000 -0.1000 0.000 -0.1500 0.1700 0.1700 0.2060 0.0015 -0.0010 0.0013 -0.0007 0.0011 0.0000 0.0030 3.2493 3.3256 3.5702 4.2896 4.3660 4.3789 3.3704 0.315 0.040 0.625 0.768 0.325 0.137 0.740 0.860 0.396 0.468 0.468 0.425 -0.126 0.028
0.500 8.20 456.6 -2.669 3.8380 -1.1150 -0.0037 3.2100 0.0390 3.210 0.0037 2.3390 -0.0110 -0.4800 0.0720 0.000 -0.0800 0.000 -0.1500 0.1700 0.1700 0.2200 0.0015 -0.0011 0.0009 -0.0007 0.0008 0.0002 0.0027 2.9818 3.1334 3.3552 3.9322 4.0779 4.0394 3.1564 0.416 0.097 0.634 0.728 0.325 0.137 0.660 0.800 0.396 0.451 0.451 0.375 -0.079 0.022
0.600 8.20 430.3 -2.599 3.5620 -1.0710 -0.0035 3.2100 0.0380 3.210 0.0033 2.2170 -0.0150 -0.4900 0.0700 0.000 -0.0600 0.000 -0.1500 0.1700 0.1700 0.2250 0.0015 -0.0012 0.0006 -0.0007 0.0006 0.0002 0.0025 2.7784 2.9215 3.0922 3.6149 3.8146 3.7366 2.9584 0.499 0.145 0.581 0.701 0.325 0.137 0.590 0.780 0.396 0.438 0.438 0.345 -0.041 0.016
0.750 8.15 410.5 -2.401 3.1520 -1.0200 -0.0032 3.2100 0.0370 3.210 0.0027 1.9460 -0.0210 -0.5000 0.0670 0.000 -0.0470 0.000 -0.1500 0.1700 0.1700 0.2170 0.0014 -0.0011 0.0003 -0.0007 0.0004 0.0002 0.0022 2.4780 2.5380 2.6572 3.1785 3.4391 3.2930 2.6556 0.600 0.197 0.497 0.685 0.325 0.137 0.500 0.730 0.396 0.420 0.420 0.300 0.005 0.009
1.000 8.10 400.0 -1.955 2.5440 -0.9500 -0.0029 3.2100 0.0350 3.210 0.0019 1.4160 -0.0280 -0.5100 0.0630 0.000 -0.0350 0.000 -0.1500 0.1700 0.1700 0.1850 0.0013 -0.0008 0.0001 -0.0008 0.0002 0.0001 0.0019 1.9252 1.9626 2.1459 2.5722 2.8056 2.6475 2.0667 0.731 0.269 0.469 0.642 0.325 0.137 0.410 0.690 0.396 0.396 0.396 0.240 0.065 0.000
1.500 8.05 400.0 -1.025 1.6360 -0.8600 -0.0026 3.2100 0.0340 3.210 0.0008 0.3940 -0.0410 -0.5200 0.0590 0.000 -0.0180 0.000 -0.1300 0.1700 0.1700 0.0830 0.0014 -0.0004 -0.0001 -0.0008 0.0001 0.0000 0.0016 0.9924 1.3568 1.3499 1.6499 1.8546 1.6842 1.3316 0.748 0.347 0.509 0.325 0.312 0.113 0.330 0.620 0.379 0.379 0.379 0.230 0.065 0.000
2.000 8.00 400.0 -0.299 1.0760 -0.8200 -0.0024 3.2100 0.0320 3.210 0.0000 -0.4170 -0.0500 -0.5300 0.0590 0.000 -0.0100 0.000 -0.1100 0.1700 0.1700 0.0450 0.0015 0.0002 0.0000 -0.0007 0.0002 0.0000 0.0014 0.4676 0.8180 0.8148 1.0658 1.3020 1.1002 0.7607 0.761 0.384 0.478 0.257 0.302 0.096 0.300 0.560 0.366 0.366 0.366 0.230 0.065 0.000
2.500 7.95 400.0 0.000 0.6580 -0.7980 -0.0022 3.2100 0.0310 3.210 0.0000 -0.7250 -0.0570 -0.5400 0.0600 0.000 -0.0050 0.000 -0.0950 0.1700 0.1700 0.0260 0.0014 0.0004 0.0000 -0.0007 0.0002 -0.0002 0.0012 0.0579 0.4389 0.3979 0.6310 0.8017 0.6737 0.3648 0.770 0.397 0.492 0.211 0.295 0.082 0.270 0.520 0.356 0.356 0.356 0.230 0.065 0.000
3.000 7.90 400.0 0.000 0.4240 -0.7930 -0.0021 3.1300 0.0300 3.130 0.0000 -0.6950 -0.0650 -0.5400 0.0590 0.000 0.0000 0.000 -0.0850 0.1700 0.1700 0.0350 0.0014 0.0007 0.0003 -0.0007 0.0004 -0.0002 0.0011 -0.1391 0.1046 0.1046 0.3882 0.5958 0.4126 0.1688 0.778 0.404 0.470 0.296 0.289 0.072 0.250 0.495 0.348 0.348 0.348 0.240 0.065 0.000
4.000 7.85 400.0 0.000 0.0930 -0.7930 -0.0020 2.9850 0.0290 2.985 0.0000 -0.6380 -0.0770 -0.5400 0.0500 0.000 0.0000 0.000 -0.0730 0.1700 0.1700 0.0530 0.0014 0.0010 0.0007 -0.0006 0.0006 -0.0002 0.0010 -0.3030 -0.1597 -0.2324 0.0164 0.3522 0.0097 -0.0323 0.790 0.397 0.336 0.232 0.280 0.055 0.220 0.430 0.335 0.335 0.335 0.270 0.065 0.000
5.000 7.80 400.0 0.000 -0.1450 -0.7930 -0.0020 2.8180 0.0280 2.818 0.0000 -0.5970 -0.0880 -0.5400 0.0430 0.000 0.0000 0.000 -0.0650 0.1700 0.1700 0.0720 0.0014 0.0013 0.0014 -0.0004 0.0008 -0.0001 0.0010 -0.4094 -0.2063 -0.5722 -0.2802 0.1874 -0.2715 -0.1516 0.799 0.378 0.228 0.034 0.273 0.041 0.190 0.400 0.324 0.324 0.324 0.300 0.065 0.000
6.000 7.80 400.0 0.000 -0.3200 -0.7930 -0.0020 2.6820 0.0270 2.682 0.0000 -0.5610 -0.0980 -0.5400 0.0380 0.000 0.0000 0.000 -0.0600 0.1700 0.1700 0.0860 0.0014 0.0015 0.0015 -0.0003 0.0011 0.0000 0.0010 -0.5010 -0.3223 -0.8631 -0.4822 -0.1243 -0.4591 -0.2217 0.807 0.358 0.151 -0.037 0.267 0.030 0.170 0.370 0.314 0.314 0.314 0.320 0.065 0.000
7.500 7.80 400.0 0.000 -0.5560 -0.7930 -0.0020 2.5150 0.0260 2.515 0.0000 -0.5300 -0.1100 -0.5400 0.0320 0.000 0.0000 0.000 -0.0550 0.1700 0.1700 0.1150 0.0014 0.0017 0.0015 -0.0002 0.0014 0.0001 0.0010 -0.6209 -0.4223 -1.1773 -0.7566 -0.3316 -0.6822 -0.3338 0.817 0.333 0.051 -0.178 0.259 0.017 0.140 0.320 0.301 0.301 0.301 0.350 0.065 0.000
10.00 7.80 400.0 0.000 -0.8600 -0.7930 -0.0020 2.3000 0.0250 2.300 0.0000 -0.4860 -0.1270 -0.5400 0.0240 0.000 0.0000 0.000 -0.0450 0.1700 0.1700 0.1510 0.0014 0.0017 0.0015 -0.0001 0.0017 0.0002 0.0010 -0.6221 -0.5909 -1.4070 -1.0870 -0.6783 -0.9173 -0.5441 0.829 0.281 -0.251 -0.313 0.250 0.000 0.100 0.280 0.286 0.286 0.286 0.350 0.065 0.000
""")
[docs]class AbrahamsonGulerce2020SSlab(AbrahamsonGulerce2020SInter):
"""
Implements the 2020 Subduction ground motion model of Abrahamson &
Gulerce (2020) for subduction in-slab earthquakes
Abrahamson N. and Gulurce Z. (2020) "Regionalized Ground-Motion Models
for Subduction Earthquakes based on the NGA-SUB Database", Pacific
Earthquake Engineering Research Center (PEER) Technical Report,
PEER 2020/25
"""
#: Required rupture parameters are magnitude and top-of-rupture depth
REQUIRES_RUPTURE_PARAMETERS = {'mag', 'ztor'}
#: Supported tectonic region type is subduction inslab
DEFINED_FOR_TECTONIC_REGION_TYPE = const.TRT.SUBDUCTION_INTRASLAB
# Long form regional aliases
REGION_ALIASES = {
"GLO": "",
"USA-AK": "Alaska",
"CAS": "Cascadia",
"CAM": "CentralAmericaMexico",
"JPN": "Japan",
"NZL": "NewZealand",
"SAM": "SouthAmerica",
"TWN": "Taiwan",
}
for region in SUPPORTED_REGIONS[1:]:
add_alias("AbrahamsonGulerce2020SInter" + REGION_ALIASES[region],
AbrahamsonGulerce2020SInter,
region=region)
add_alias("AbrahamsonGulerce2020SSlab" + REGION_ALIASES[region],
AbrahamsonGulerce2020SSlab,
region=region)
