Source code for openquake.hazardlib.gsim.bayless_abrahamson_2018
# -*- coding: utf-8 -*-# vim: tabstop=4 shiftwidth=4 softtabstop=4## Copyright (C) 2021 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:`BaylessAbrahamson2018`"""importosimportnumpyasnpfromopenquake.hazardlib.gsim.baseimportGMPE,CoeffsTablefromopenquake.hazardlibimportconstfromopenquake.hazardlib.imtimportEASBA_COEFFS=os.path.join(os.path.dirname(__file__),"bayless_abrahamson_2018.csv")def_magnitude_scaling(C,ctx):""" Compute the magnitude scaling term """t1=C['c1']t2=C['c2']*(ctx.mag-6.)tmp=np.log(1.0+np.exp(C['cn']*(C['cM']-ctx.mag)))t3=((C['c2']-C['c3'])/C['cn'])*tmpreturnt1+t2+t3def_path_scaling(C,ctx):""" Compute path scaling term """tmp=np.maximum(ctx.mag-C['chm'],0.)t1=C['c4']*np.log(ctx.rrup+C['c5']*np.cosh(C['c6']*tmp))tmp=(ctx.rrup**2+50.**2)**0.5t2=(-0.5-C['c4'])*np.log(tmp)t3=C['c7']*ctx.rrupreturnt1+t2+t3def_normal_fault_effect(C,ctx):""" Compute the correction coefficient for the normal faulting """idx=(ctx.rake>-150)&(ctx.rake<-30)fnm=np.zeros_like(ctx.rake)fnm[idx]=1.0returnC['c10']*fnmdef_ztor_scaling(C,ctx):""" Compute top of rupture term """returnC['c9']*np.minimum(ctx.ztor,20.)def_get_ln_ir_outcrop(self,ctx):""" Compute the natural logarithm of peak PGA at rock outcrop, Ir. See eq.10e page 2092. """# Set imt and cofficients. The value of frequency is consistent with# the one used in the matlab script - line 156im=EAS(5.0118725)tc=self.COEFFS[im]# Get meanmean=(_magnitude_scaling(tc,ctx)+_path_scaling(tc,ctx)+_ztor_scaling(tc,ctx)+_normal_fault_effect(tc,ctx))return1.238+0.846*meandef_linear_site_response(C,ctx):""" Compute the site response term """tmp=np.minimum(ctx.vs30,1000.)fsl=C['c8']*np.log(tmp/1000.)returnfsldef_get_basin_term(C,ctx,region=None):""" Compute the soil depth scaling term """# Set the c11 coefficient - See eq.13b at page 2093c11=np.ones_like(ctx.vs30)*C['c11a']c11[(ctx.vs30<=300)&(ctx.vs30>200)]=C['c11b']c11[(ctx.vs30<=500)&(ctx.vs30>300)]=C['c11c']c11[(ctx.vs30>500)]=C['c11d']# Compute the Z1ref parametertmp=(ctx.vs30**4+610**4)/(1360**4+610**4)z1ref=1/1000.*np.exp(-7.67/4*np.log(tmp))# Return the fz1 parameter. The z1pt0 is converted from m (standard in OQ)# to km as indicated in the papertmp=np.minimum(ctx.z1pt0/1000,np.ones_like(ctx.z1pt0)*2.0)+0.01returnc11*np.log(tmp/(z1ref+0.01))def_get_stddevs(C,ctx):""" Compute the standard deviations """# Set components of stdtau=np.ones_like(ctx.mag)*C['s1']phi_s2s=np.ones_like(ctx.mag)*C['s3']phi_ss=np.ones_like(ctx.mag)*C['s5']above=ctx.mag>6tau[above]=C['s2']phi_s2s[above]=C['s4']phi_ss[above]=C['s6']below=(ctx.mag>4)&(ctx.mag<=6)tau[below]=C['s1']+(C['s2']-C['s1'])/2.*(ctx.mag[below]-4)phi_s2s[below]=C['s3']+(C['s4']-C['s3'])/2.*(ctx.mag[below]-4)phi_ss[below]=C['s5']+(C['s6']-C['s5'])/2.*(ctx.mag[below]-4)# Collect the requested stdssigma=np.sqrt(tau**2+phi_s2s**2+phi_ss**2+C['c1a']**2)phi=np.sqrt(phi_s2s**2+phi_ss**2)sigma=np.ones_like(ctx.vs30)*sigmatau=np.ones_like(ctx.vs30)*tauphi=np.ones_like(ctx.vs30)*phireturnsigma,tau,phidef_get_nl_site_response(C,ctx,imt,ln_ir_outcrop,freq_nl,coeff_nl):""" :param ln_ir_outcrop: The peak ground acceleration (PGA [g]) at rock outcrop :param freq_nl: Frequencies for the coefficients f3, f4 and f5 used to compute the non-linear term :param coeff_nl: A :class:`numpy.ndarray` instance of shape [# freq, 3] with the values of the coefficients f3, f4 and f5 which are used to compute the non-linear term """vsref=760.0NSITES=len(ctx.vs30)NFREQS=coeff_nl.shape[0]# Updating the matrix with the coefficients. This has shape (number of# frequencies) x (number of coeffs) x (number of sites)coeff_nl=np.expand_dims(coeff_nl,2)coeff_nl=np.repeat(coeff_nl,NSITES,2)# Updating the matrix with Vs30 values. This has shape (number of# frequencies) x (number of sites)# vs30 = np.matlib.repmat(ctx.vs30, NFREQS, 1)# ln_ir_outcrop = np.matlib.repmat(ln_ir_outcrop, NFREQS, 1)vs30=np.tile(ctx.vs30,(NFREQS,1))ln_ir_outcrop=np.tile(ln_ir_outcrop,(NFREQS,1))# Computingt1=np.exp(coeff_nl[:,2]*(np.minimum(vs30,vsref)-360.))t2=np.exp(coeff_nl[:,2]*(vsref-360.))f2=coeff_nl[:,1]*(t1-t2)f3=coeff_nl[:,0]fnl_0=f2*np.log((np.exp(ln_ir_outcrop)+f3)/f3)idxs=np.argmin(fnl_0,axis=0)# Applying correction as described at page 2093 in BA18fnl=[]fori,jinenumerate(idxs):fnl_0[j:,i]=min(fnl_0[:,i])tmp=np.interp(imt.frequency,freq_nl,fnl_0[:,i])fnl.append(tmp)returnnp.array(fnl)def_get_dimunition_factor(ctx,imt):max_freq=23.988321kappa=np.exp(-0.4*np.log(ctx.vs30/760)-3.5)D=np.exp(-np.pi*kappa*(imt.frequency-max_freq))returnDdef_get_mean_linear(C,ctx):mean=(_magnitude_scaling(C,ctx)+_path_scaling(C,ctx)+_linear_site_response(C,ctx)+_ztor_scaling(C,ctx)+_normal_fault_effect(C,ctx)+_get_basin_term(C,ctx))returnmeandef_get_mean(self,C,ctx,imt):ifimt.frequency>=24.:im=EAS(23.988321)C=self.COEFFS[im]t1=np.exp(_get_mean_linear(C,ctx))mean=np.log(_get_dimunition_factor(ctx,imt)*t1)else:mean=_get_mean_linear(C,ctx)returnmean
[docs]classBaylessAbrahamson2018(GMPE):""" Implements the Bayless and Abrahamson (2018, 2019) model. References: - Bayless, J., and N. A. Abrahamson (2018b). An empirical model for Fourier amplitude spectra using the NGA-West2 database, PEER Rept. No. 2018/07, Pacific Earthquake Engineering Research Center, University of California, Berkeley, California. - Bayless, J. and N.A. Abrahamson (2019). Summary of the BA18 Ground-Motion Model for Fourier Amplitude Spectra for Crustal Earthquakes in California. Bull. Seism. Soc. Am., 109(5): 2088–2105 Disclaimer: The authors describe a smoothing technique that needs to be applied to the non linear component of the site response. We did not implement these smoothing functions in this initial versions since the match with the values in the verification tables is good even without it. """#: Supported tectonic region type is active shallow crust, see title!DEFINED_FOR_TECTONIC_REGION_TYPE=const.TRT.ACTIVE_SHALLOW_CRUST#: Supported intensity measure typesDEFINED_FOR_INTENSITY_MEASURE_TYPES={EAS}#: Supported intensity measure componentDEFINED_FOR_INTENSITY_MEASURE_COMPONENT=const.IMC.HORIZONTAL#: Supported standard deviation typesDEFINED_FOR_STANDARD_DEVIATION_TYPES={const.StdDev.TOTAL,const.StdDev.INTER_EVENT,const.StdDev.INTRA_EVENT}#: Required site parametersREQUIRES_SITES_PARAMETERS={'vs30','z1pt0'}#: Required rupture parametersREQUIRES_RUPTURE_PARAMETERS={'mag','rake','ztor'}#: Required distance measuresREQUIRES_DISTANCES={'rrup'}