Merge pull request 'upload Lasocki2013 ground motion model' (!20) from lasocki2013_gmm into master

Reviewed-on: #20
This commit is contained in:
asia 2025-07-11 11:28:58 +02:00
commit 9e75e0e768

142
gsim/lasocki_2013.py Normal file
View File

@ -0,0 +1,142 @@
# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4
#
# Copyright (C) 2013-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:`Lasocki2013`.
"""
import numpy as np
from scipy.constants import g
from openquake.hazardlib.gsim.base import GMPE, CoeffsTable
from openquake.hazardlib import const
from openquake.hazardlib.imt import PGA, PGV, SA
def get_magnitude_energy(mag):
"""
Converts magnitude to energy-based magnitude term.
"""
return 1.15 + 1.96 * mag
def get_distance_term(coeffs, repi):
"""
Computes the distance term using the given GMPE coefficients.
"""
R_h = np.sqrt(repi ** 2 + coeffs["c7"] ** 2)
return np.log10(R_h)
def get_standard_deviation(coeffs, coeffs_cov, magE, repi):
"""
Computes the standard deviation term.
"""
Cb = np.array(list(coeffs_cov)).reshape(3,3)
R_h = np.sqrt(repi ** 2 + coeffs["c7"] ** 2)
X0 = np.array([1, magE[0], np.log10(R_h[0]**2)])
variance_term = np.sqrt(X0 @ Cb @ X0 + coeffs["sigma"]**2)
return variance_term
class Lasocki2013(GMPE):
"""
Implement equation developed by Lasocki in "REPORT ON THE ATTENUATION
RELATIONS OF PEAK GROUND ACCELERATION AND SPECTRAL ORDINATES OF GROUND
MOTION FOR MINING-INDUCED SEISMIC EVENTS IN THE REGION OF THE ŻELAZNY
MOST REPOSITORY", 2009.
Equation coefficients provided for the random horizontal component
"""
#: Supported tectonic region type is induced, given
#: that the equations have been derived for the LGCD mining area
DEFINED_FOR_TECTONIC_REGION_TYPE = const.TRT.INDUCED
#: Supported intensity measure types are spectral acceleration,
#: and peak ground acceleration
DEFINED_FOR_INTENSITY_MEASURE_TYPES = {PGA, SA}
#: Supported intensity measure component is random horizontal
#: :attr:`~openquake.hazardlib.const.IMC.RANDOM_HORIZONTAL`,
DEFINED_FOR_INTENSITY_MEASURE_COMPONENT = const.IMC.RANDOM_HORIZONTAL
# DEFINED_FOR_INTENSITY_MEASURE_COMPONENT = const.IMC.VERTICAL
#: Supported standard deviation type is total
DEFINED_FOR_STANDARD_DEVIATION_TYPES = {const.StdDev.TOTAL}
#: site params are not required
REQUIRES_SITES_PARAMETERS = set()
#: Required rupture parameter is magnitude
REQUIRES_RUPTURE_PARAMETERS = {'mag'}
#: Required distance measure is epicentral distance
#: see paragraph 'Predictor Variables', page 6.
REQUIRES_DISTANCES = {'repi'}
def compute(self, ctx: np.recarray, imts, mean, sig, tau, phi):
"""
Computes mean ground motion values and standard deviations for Lasocki (2013) GMPE.
Parameters:
repi (float): Epicentral distance (m)
mag (float): Earthquake magnitude
imts (list of str): List of intensity measure types (e.g., 'PHA', 'PVA')
mean (np.array): Array to store computed mean values
sig (np.array): Array to store computed standard deviations
"""
# Loop through each IMT and compute values
for i, imt in enumerate(imts):
C = self.COEFFS[imt]
C_cov = self.COEFFS_COV[imt]
mag = ctx.mag
repi = ctx.repi*1000.0 # Convert distance from km to m
# Compute magnitude energy term
magE = get_magnitude_energy(mag)
# Compute GMPE terms
mag_term = C['c1'] + C['c2'] * magE
dist_term = C['c5'] * get_distance_term(C, repi)
# Compute mean ground motion
imean = mag_term + dist_term
mean_value = np.log((10 ** imean) / g) # Convert to natural log scale and divide by g
mean[i] = mean_value
# Compute standard deviation
sigma = get_standard_deviation(C, C_cov, magE, repi)
sig[i] = sigma
#: coefficient table provided by report
COEFFS = CoeffsTable(sa_damping=5, table="""\
IMT c1 c2 c5 c7 sigma N
pga 1.25 0.31 -1.34 558 0.196 1196
0.6 -3.86 0.55 -0.65 183 0.242 1194
1.0 -3.94 0.62 -0.66 308 0.262 1197
2.0 -1.14 0.47 -1.02 741 0.235 1195
5.0 0.99 0.31 -1.15 690 0.234 1206
10.0 3.06 0.27 -1.65 906 0.203 1192
20.0 2.62 0.27 -1.60 435 0.196 1191
50.0 2.09 0.27 -1.48 375 0.204 1191
""")
COEFFS_COV = CoeffsTable(sa_damping=5, table="""\
IMT Cb00 Cb01 Cb02 Cb10 Cb11 Cb12 Cb20 Cb21 Cb22
pga 0.005586 -0.000376 -0.000752 -0.000376 0.000103 -0.000111 -0.000752 -0.000111 0.000440
0.6 0.007509 -0.000662 -0.000688 -0.000662 0.000161 -0.000154 -0.000688 -0.000154 0.000516
1.0 0.009119 -0.00075 -0.000948 -0.00075 0.000189 -0.000187 -0.000948 -0.000187 0.000657
2.0 0.008563 -0.000514 -0.001282 -0.000514 0.000147 -0.000164 -0.001282 -0.000164 0.000696
5.0 0.008283 -0.000516 -0.001202 -0.000516 0.000145 -0.000161 -0.001202 -0.000161 0.000668
10.0 0.006904 -0.00036 -0.001145 -0.00036 0.000108 -0.000126 -0.001145 -0.000126 0.000578
20.0 0.005389 -0.000396 -0.000658 -0.000396 0.000104 -0.000107 -0.000658 -0.000107 0.000408
50.0 0.005874 -0.000449 -0.000678 -0.000449 0.000114 -0.000114 -0.000678 -0.000114 0.000428
""")