Source code for openquake.hazardlib.gsim.eshm20_craton
# -*- coding: utf-8 -*-# vim: tabstop=4 shiftwidth=4 softtabstop=4## Copyright (C) 2014-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:`ESHM20Craton`"""importnumpyasnpfromopenquake.hazardlib.gsim.baseimportGMPE,CoeffsTable,add_aliasfromopenquake.hazardlib.imtimportPGA,SAfromopenquake.hazardlibimportconstfromopenquake.hazardlib.gsim.nga_eastimport(get_tau_at_quantile,get_phi_ss_at_quantile,TAU_EXECUTION,TAU_SETUP,PHI_SETUP,get_phi_ss,NGAEastGMPE,_get_f760,get_nonlinear_stddev,get_linear_stddev,_get_fv,get_fnl)fromopenquake.hazardlib.gsim.usgs_ceus_2019importget_stewart_2019_phis2sfromopenquake.hazardlib.gsim.kotha_2020importKothaEtAl2020ESHM20CONSTANTS={"Mref":4.5,"Rref":1.,"Mh":6.2,"h":5.0}
[docs]defget_distance_term(C,mag,rrup):""" Returns the distance attenuation factor """rval=np.sqrt(rrup**2.+CONSTANTS["h"]**2.)rref_val=np.sqrt(CONSTANTS["Rref"]**2.+CONSTANTS["h"]**2.)f_r=(C["c1"]+C["c2"]*(mag-CONSTANTS["Mref"]))*\
np.log(rval/rref_val)+(C["c3"]*(rval-rref_val)/100.)returnf_r
[docs]defget_hard_rock_mean(C,ctx):""" Returns the mean and standard deviations for the reference very hard rock condition (Vs30 = 3000 m/s) """returnget_magnitude_scaling(C,ctx.mag)+get_distance_term(C,ctx.mag,ctx.rrup)
[docs]defget_magnitude_scaling(C,mag):""" Returns the magnitude scaling term """d_m=mag-CONSTANTS["Mh"]returnnp.where(mag<=CONSTANTS["Mh"],C["e1"]+C["b1"]*d_m+C["b2"]*d_m**2.0,C["e1"]+C["b3"]*d_m)
[docs]defget_site_amplification(site_epsilon,imt,pga_r,ctx):""" Returns the sum of the linear (Stewart et al., 2019) and non-linear (Hashash et al., 2019) amplification terms """# Get the coefficients for the IMTC_LIN=NGAEastGMPE.COEFFS_LINEAR[imt]C_F760=NGAEastGMPE.COEFFS_F760[imt]C_NL=NGAEastGMPE.COEFFS_NONLINEAR[imt]ifstr(imt).startswith("PGA"):period=0.01elifstr(imt).startswith("PGV"):period=0.5else:period=imt.period# Get f760f760=_get_f760(C_F760,ctx.vs30,NGAEastGMPE.CONSTANTS)# Get the linear amplification factorf_lin=_get_fv(C_LIN,ctx,f760,NGAEastGMPE.CONSTANTS)# Get the nonlinear amplification from Hashash et al., (2017)f_nl,f_rk=get_fnl(C_NL,pga_r,ctx.vs30,period)# Mean amplificationampl=f_lin+f_nl# If an epistemic uncertainty is required then retrieve the epistemic# sigma of both models and multiply by the input epsilonifsite_epsilon:# In the case of the linear model sigma_f760 and sigma_fv are# assumed independent and the resulting sigma_flin is the root# sum of squares (SRSS)f760_stddev=_get_f760(C_F760,ctx.vs30,NGAEastGMPE.CONSTANTS,is_stddev=True)f_lin_stddev=np.sqrt(f760_stddev**2.+get_linear_stddev(C_LIN,ctx.vs30,NGAEastGMPE.CONSTANTS)**2)# Likewise, the epistemic uncertainty on the linear and nonlinear# model are assumed independent and the SRSS is takenf_nl_stddev=get_nonlinear_stddev(C_NL,ctx.vs30)*f_rksite_epistemic=np.sqrt(f_lin_stddev**2.+f_nl_stddev**2.)ampl+=(site_epsilon*site_epistemic)returnampl
[docs]defget_stddevs(ergodic,tau_model,TAU,PHI_SS,imt,ctx):""" Returns the standard deviations for either the ergodic or non-ergodic models """phi=get_phi_ss(imt,ctx.mag,PHI_SS)ifergodic:phi_s2s=get_stewart_2019_phis2s(imt,ctx.vs30)phi=np.sqrt(phi**2.+phi_s2s**2.)tau=TAU_EXECUTION[tau_model](imt,ctx.mag,TAU)sigma=np.sqrt(tau**2.+phi**2.)return[sigma,tau,phi]
[docs]classESHM20Craton(GMPE):""" Implements a scalable backbone GMPE for application to stable cratonic regions (primarily intended for cratonic Europe). The median ground motion is determined by fitting a parametric model to an extensive set of ground motion scenarios from the suite of NGA East ground motion models for 800 m/s site class. The form of the parametric model is based on that of :class:`openquake.hazardlib.gsim.kotha_2019.KothaEtAl2019`, and the scaling in terms of the number of standard deviations of the epistemic uncertainty (sigma). The aleatory uncertainty model is that of Al Atik (2015), which is common to all NGA East ground motion models and configurable by the user. :param float epsilon: Number of standard deviations above or below the median to be applied to the epistemic uncertainty sigma :param str tau_model: Choice of model for the inter-event standard deviation (tau), selecting from "global" {default}, "cena" or "cena_constant" :param str phi_model: Choice of model for the single-station intra-event standard deviation (phi_ss), selecting from "global" {default}, "cena" or "cena_constant" :param TAU: Inter-event standard deviation model :param PHI_SS: Single-station standard deviation model :param PHI_S2SS: Station term for ergodic standard deviation model :param bool ergodic: True if an ergodic model is selected, False otherwise :param float tau_quantile: Epistemic uncertainty quantile for the inter-event standard deviation models. Float in the range 0 to 1, or None (mean value used) :param float phi_ss_quantile: Epistemic uncertainty quantile for the intra-event standard deviation models. Float in the range 0 to 1, or None (mean value used) :param float phi_s2ss_quantile: Epistemic uncertainty quantile for the site-to-site standard deviation models. Float in the range 0 to 1, or None (mean value used) :param float site_epsilon: Number of standard deviations above or below median for the uncertainty in the site amplification model """#: Supported tectonic region type is 'active shallow crust'DEFINED_FOR_TECTONIC_REGION_TYPE=const.TRT.STABLE_CONTINENTAL#: The GMPE is defined only for PGA and SADEFINED_FOR_INTENSITY_MEASURE_TYPES={PGA,SA}#: Supported intensity measure component is the geometric mean of two#: horizontal componentsDEFINED_FOR_INTENSITY_MEASURE_COMPONENT=const.IMC.RotD50#: Supported standard deviation types are inter-event, intra-event#: and totalDEFINED_FOR_STANDARD_DEVIATION_TYPES={const.StdDev.TOTAL,const.StdDev.INTER_EVENT,const.StdDev.INTRA_EVENT}#: Median calibrated for Vs30 3000 m/s Vs30, no site term required Vs30REQUIRES_SITES_PARAMETERS={'vs30'}#: Requires only magnitudeREQUIRES_RUPTURE_PARAMETERS={'mag'}#: Required distance measure is RrupREQUIRES_DISTANCES={'rrup'}#: Defined for a reference velocity of 3000 m/sDEFINED_FOR_REFERENCE_VELOCITY=3000.0def__init__(self,**kwargs):""" Instantiates the class with additional terms controlling both the epistemic uncertainty in the median and the preferred aleatory uncertainty model ('global', 'cena_constant', 'cena'), and the quantile of the epistemic uncertainty model (float in the range 0 to 1, or None) """super().__init__(**kwargs)self.epsilon=kwargs.get("epsilon",0.0)self.tau_model=kwargs.get("tau_model","global")self.phi_model=kwargs.get("phi_model","global")self.ergodic=kwargs.get("ergodic",True)self.tau_quantile=kwargs.get("tau_quantile",None)self.phi_ss_quantile=kwargs.get("phi_ss_quantile",None)self.site_epsilon=kwargs.get("site_epsilon",0.0)self.PHI_S2SS=None# define the standard deviation model from the NGA East aleatory# uncertainty model according to the calibrations specified by the user# setup tauself.TAU=get_tau_at_quantile(TAU_SETUP[self.tau_model]["MEAN"],TAU_SETUP[self.tau_model]["STD"],self.tau_quantile)# setup phiself.PHI_SS=get_phi_ss_at_quantile(PHI_SETUP[self.phi_model],self.phi_ss_quantile)
[docs]defcompute(self,ctx:np.recarray,imts,mean,sig,tau,phi):""" Returns the mean and standard deviations """C_ROCK=self.COEFFS[PGA()]pga_r=get_hard_rock_mean(C_ROCK,ctx)form,imtinenumerate(imts):C=self.COEFFS[imt]# Get the desired spectral acceleration on rockifimt.string!="PGA":# Calculate the ground motion at required spectral period for# the reference rockmean[m]=get_hard_rock_mean(C,ctx)else:# Avoid re-calculating PGA if that was already done!mean[m]=np.copy(pga_r)mean[m]+=get_site_amplification(self.site_epsilon,imt,np.exp(pga_r),ctx)# Get standard deviation modelsig[m],tau[m],phi[m]=get_stddevs(self.ergodic,self.tau_model,self.TAU,self.PHI_SS,imt,ctx)ifself.epsilon:# If requested, apply epistemic uncertaintymean[m]+=(self.epsilon*C["sigma_mu"])
# Add aliases for the ESHM20 adjustments to the Craton Model# Define the adjustment factors for the 3- and 5-pnt Gaussian approximation# according to Miller & Rice (1983)MILLER_RICE_GAUSS_3PNT=[-1.732051,0.0,1.732051]MILLER_RICE_GAUSS_5PNT=[-2.856970,-1.355630,0.0,1.355630,2.856970]STRESS_BRANCHES=["VLow","Low","Mid","High","VHigh"]SITE_BRANCHES=["Low","Mid","High"]# Get the 15 branch set of aliasesforstress,eps1inzip(STRESS_BRANCHES,MILLER_RICE_GAUSS_5PNT):forsite,eps2inzip(SITE_BRANCHES,MILLER_RICE_GAUSS_3PNT):alias="ESHM20Craton{:s}Stress{:s}Site".format(stress,site)add_alias(alias,ESHM20Craton,epsilon=eps1,site_epsilon=eps2)# Add on the four branches of the KothaEtAl2020ESHM20 adjustmentsadd_alias("ESHM20CratonShallowHighStressMidAtten",KothaEtAl2020ESHM20,sigma_mu_epsilon=1.732051)add_alias("ESHM20CratonShallowHighStressSlowAtten",KothaEtAl2020ESHM20,sigma_mu_epsilon=1.732051,c3_epsilon=1.732051)add_alias("ESHM20CratonShallowMidStressMidAtten",KothaEtAl2020ESHM20,)add_alias("ESHM20CratonShallowMidStressSlowAtten",KothaEtAl2020ESHM20,c3_epsilon=1.732051)