Source code for openquake.hazardlib.gsim.berge_thierry_2003

# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4
#
# Copyright (C) 2014-2018 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:`BergeThierryEtAl2003SIGMA`.
"""
import copy
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, SA


[docs]class BergeThierryEtAl2003SIGMA(GMPE): """ Implements GMPE developed by Catherine Berge-Thierry, Fabrice Cotton, Oona Scoti, Daphne-Anne Griot-Pommera, and Yoshimitsu Fukushima and published as "New Empirical Response Spectral Attenuation Laws For Moderate European Earthquakes" (2003, Journal of Earthquake Engineering, 193-222) The class implements also adjustment of the sigma value as required by the SIGMA project to make the GMPE usable with Mw (the GMPE was originally developed for Ms). """ #: Supported tectonic region type is active shallow crust, see #: `Introduction`, page 194. DEFINED_FOR_TECTONIC_REGION_TYPE = const.TRT.ACTIVE_SHALLOW_CRUST #: Supported intensity measure types are spectral acceleration, and peak #: ground acceleration. The original manuscript provide coefficients only #: SA. For PGA, coefficients are assumed equal to the ones of SA for the #: smallest period (0.03 s) DEFINED_FOR_INTENSITY_MEASURE_TYPES = set([ PGA, SA ]) #: Supported intensity measure component is horizontal, see page 196. DEFINED_FOR_INTENSITY_MEASURE_COMPONENT = const.IMC.HORIZONTAL #: Supported standard deviation type is total, see table 3, page 203 DEFINED_FOR_STANDARD_DEVIATION_TYPES = set([ const.StdDev.TOTAL ]) #: Required site parameters is Vs30, used to distinguish between rock sites #: (Vs30 >= 800) m/s and alluvium sites (300 < Vs < 800), see section 2.2.3 #: page 201 REQUIRES_SITES_PARAMETERS = set(('vs30', )) #: Required rupture parameters is magnitude, see equation 1 page 201 REQUIRES_RUPTURE_PARAMETERS = set(('mag', )) #: Required distance measure is hypocentral distance, see equation 1 page #: 201 REQUIRES_DISTANCES = set(('rhypo', ))
[docs] def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types): return self._get_mean_and_stddevs(sites, rup, dists, imt, stddev_types)
def _get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types): """ See :meth:`superclass method <.base.GroundShakingIntensityModel.get_mean_and_stddevs>` for spec of input and result values. """ # extract dictionaries of coefficients specific to required # intensity measure type C = self.COEFFS[imt] # clip distance at 4 km, minimum distance for which the equation is # valid (see section 2.2.4, page 201). This also avoids singularity # in the equation rhypo = dists.rhypo rhypo[rhypo < 4.] = 4. mean = C['a'] * rup.mag + C['b'] * rhypo - np.log10(rhypo) mean[sites.vs30 >= 800] += C['c1'] mean[sites.vs30 < 800] += C['c2'] # convert from log10 to ln, and from cm/s2 to g mean = mean * np.log(10) - 2 * np.log(10) - np.log(g) stddevs = self._get_stddevs(C, stddev_types, rhypo.shape[0]) return mean, stddevs def _get_stddevs(self, C, stddev_types, num_sites): """ Return total standard deviation. """ assert all(stddev_type in self.DEFINED_FOR_STANDARD_DEVIATION_TYPES for stddev_type in stddev_types) # adjustement for sigma, as explained in section 3.7.3 in # Delivrable_SIGMA-V02-2012-D4-18+reviews sigma = np.zeros(num_sites) + C['sigma'] * np.log(10) sigma = np.sqrt(sigma ** 2 + (C['a'] ** 2) * (0.2 ** 2)) stddevs = [sigma for _ in stddev_types] return stddevs #: Coefficient tables are constructed from the electronic suplements of #: the original paper. Original coefficients in function of frequency. COEFFS = CoeffsTable(sa_damping=5, table="""\ IMT a b c1 c2 sigma pga 0.3118000 -0.0009303 1.537000 1.57300 0.2923 0.029412 0.3118000 -0.0009303 1.537000 1.57300 0.2923 0.030000 0.3114000 -0.0009334 1.541000 1.57600 0.2924 0.032258 0.3097000 -0.0009422 1.558000 1.58900 0.2928 0.034000 0.3083000 -0.0009547 1.573000 1.60200 0.2935 0.035714 0.3068000 -0.0009822 1.593000 1.61800 0.2947 0.040000 0.3033000 -0.0011190 1.653000 1.66500 0.2982 0.045455 0.3016000 -0.0012820 1.701000 1.70000 0.3015 0.050000 0.2992000 -0.0013410 1.740000 1.72900 0.3009 0.053000 0.2981000 -0.0014290 1.766000 1.74900 0.3022 0.055556 0.2969000 -0.0014320 1.785000 1.76500 0.3027 0.058824 0.2960000 -0.0014600 1.809000 1.78400 0.3040 0.059999 0.2960000 -0.0014720 1.814000 1.78800 0.3041 0.062500 0.2965000 -0.0015070 1.826000 1.79600 0.3047 0.064998 0.2944000 -0.0015150 1.851000 1.81700 0.3059 0.066667 0.2933000 -0.0015130 1.865000 1.82900 0.3059 0.068966 0.2924000 -0.0015460 1.881000 1.84200 0.3048 0.069999 0.2920000 -0.0015700 1.888000 1.84900 0.3046 0.071429 0.2909000 -0.0015830 1.901000 1.86100 0.3047 0.074074 0.2879000 -0.0015550 1.926000 1.88700 0.3059 0.075002 0.2871000 -0.0015400 1.933000 1.89300 0.3066 0.076923 0.2857000 -0.0015200 1.947000 1.90600 0.3072 0.080000 0.2866000 -0.0015310 1.954000 1.91100 0.3061 0.083333 0.2858000 -0.0015520 1.968000 1.92700 0.3058 0.084998 0.2856000 -0.0015830 1.976000 1.93500 0.3053 0.086957 0.2848000 -0.0015710 1.987000 1.94600 0.3042 0.090001 0.2819000 -0.0015360 2.011000 1.97300 0.3020 0.090909 0.2809000 -0.0015280 2.020000 1.98100 0.3016 0.095238 0.2781000 -0.0015430 2.054000 2.01000 0.3017 0.100000 0.2786000 -0.0014740 2.059000 2.01600 0.3019 0.105260 0.2776000 -0.0014220 2.072000 2.03400 0.3043 0.110000 0.2783000 -0.0014420 2.075000 2.04500 0.3073 0.111110 0.2783000 -0.0014380 2.077000 2.04700 0.3079 0.117650 0.2793000 -0.0014380 2.081000 2.05600 0.3103 0.120000 0.2806000 -0.0014360 2.076000 2.05400 0.3106 0.125000 0.2831000 -0.0014150 2.066000 2.05100 0.3121 0.129030 0.2863000 -0.0014400 2.052000 2.04200 0.3137 0.130010 0.2867000 -0.0014250 2.050000 2.04000 0.3142 0.133330 0.2887000 -0.0013970 2.040000 2.03300 0.3156 0.137930 0.2903000 -0.0013380 2.032000 2.02800 0.3161 0.140000 0.2915000 -0.0013220 2.027000 2.02200 0.3166 0.142860 0.2933000 -0.0013070 2.018000 2.01400 0.3173 0.148150 0.2950000 -0.0012560 2.009000 2.01000 0.3191 0.149990 0.2955000 -0.0012180 2.004000 2.00900 0.3196 0.153850 0.2955000 -0.0010520 1.997000 2.00700 0.3201 0.160000 0.2939000 -0.0008056 1.996000 2.01300 0.3200 0.166670 0.2952000 -0.0007097 1.989000 2.00600 0.3215 0.170010 0.2974000 -0.0006986 1.978000 1.99400 0.3230 0.173910 0.3016000 -0.0007341 1.957000 1.97000 0.3247 0.179990 0.3089000 -0.0007793 1.915000 1.93000 0.3267 0.181820 0.3109000 -0.0007826 1.901000 1.91900 0.3271 0.190010 0.3147000 -0.0007369 1.868000 1.89400 0.3262 0.190480 0.3149000 -0.0007337 1.867000 1.89300 0.3262 0.200000 0.3167000 -0.0006889 1.843000 1.88100 0.3250 0.208330 0.3196000 -0.0006719 1.814000 1.86100 0.3261 0.217390 0.3254000 -0.0006750 1.770000 1.82500 0.3281 0.219780 0.3271000 -0.0006918 1.758000 1.81500 0.3292 0.227270 0.3303000 -0.0006678 1.726000 1.79200 0.3320 0.238100 0.3340000 -0.0006171 1.683000 1.76200 0.3367 0.239980 0.3344000 -0.0005988 1.677000 1.75800 0.3371 0.250000 0.3365000 -0.0005750 1.651000 1.73600 0.3394 0.259740 0.3430000 -0.0007075 1.609000 1.69700 0.3422 0.263160 0.3442000 -0.0007200 1.599000 1.68800 0.3429 0.277780 0.3501000 -0.0007520 1.550000 1.64500 0.3444 0.280030 0.3511000 -0.0007530 1.542000 1.63800 0.3447 0.290020 0.3555000 -0.0007836 1.506000 1.60500 0.3458 0.300030 0.3590000 -0.0008520 1.477000 1.58100 0.3477 0.303030 0.3602000 -0.0008737 1.466000 1.57300 0.3483 0.316960 0.3671000 -0.0009272 1.412000 1.52500 0.3491 0.320000 0.3690000 -0.0009468 1.397000 1.51200 0.3487 0.333330 0.3742000 -0.0010100 1.352000 1.47200 0.3474 0.340020 0.3752000 -0.0010060 1.337000 1.46100 0.3469 0.344830 0.3760000 -0.0009698 1.326000 1.45200 0.3471 0.357140 0.3807000 -0.0009114 1.286000 1.41500 0.3481 0.359970 0.3822000 -0.0009039 1.275000 1.40500 0.3484 0.370370 0.3867000 -0.0008635 1.237000 1.37200 0.3492 0.379940 0.3909000 -0.0008074 1.199000 1.33900 0.3500 0.384620 0.3931000 -0.0007955 1.179000 1.32100 0.3507 0.400000 0.3997000 -0.0007078 1.119000 1.26700 0.3517 0.416670 0.4028000 -0.0006613 1.078000 1.23200 0.3512 0.419990 0.4034000 -0.0006513 1.070000 1.22600 0.3513 0.434780 0.4070000 -0.0006167 1.029000 1.19100 0.3521 0.439950 0.4089000 -0.0006118 1.011000 1.17400 0.3527 0.454550 0.4148000 -0.0005931 0.960100 1.12500 0.3547 0.459980 0.4165000 -0.0005816 0.944000 1.10900 0.3549 0.476190 0.4222000 -0.0005404 0.893100 1.05800 0.3555 0.480080 0.4239000 -0.0005484 0.879500 1.04500 0.3556 0.500000 0.4323000 -0.0005680 0.815000 0.97970 0.3555 0.520020 0.4372000 -0.0005396 0.764200 0.93240 0.3568 0.526320 0.4379000 -0.0005050 0.752200 0.92080 0.3570 0.539960 0.4394000 -0.0004330 0.727100 0.89590 0.3574 0.555560 0.4418000 -0.0003601 0.694100 0.86610 0.3587 0.559910 0.4425000 -0.0003380 0.684400 0.85710 0.3592 0.580050 0.4472000 -0.0002702 0.635700 0.81170 0.3610 0.588240 0.4492000 -0.0002522 0.615700 0.79260 0.3609 0.599880 0.4516000 -0.0002175 0.589500 0.76820 0.3603 0.619960 0.4559000 -0.0001953 0.547800 0.72820 0.3604 0.625000 0.4569000 -0.0001995 0.538300 0.71870 0.3609 0.640200 0.4596000 -0.0001666 0.510600 0.69100 0.3625 0.660070 0.4637000 -0.0001549 0.472500 0.65200 0.3647 0.666670 0.4655000 -0.0001550 0.457300 0.63630 0.3649 0.679810 0.4688000 -0.0001668 0.428400 0.60810 0.3655 0.699790 0.4732000 -0.0001700 0.385700 0.56760 0.3667 0.714290 0.4771000 -0.0002019 0.354800 0.53540 0.3680 0.750190 0.4847000 -0.0003009 0.287100 0.46810 0.3700 0.769230 0.4875000 -0.0003122 0.254500 0.43800 0.3703 0.800000 0.4940000 -0.0002568 0.190600 0.37820 0.3714 0.833330 0.5010000 -0.0001932 0.126400 0.31450 0.3746 0.850340 0.5040000 -0.0001433 0.096150 0.28420 0.3758 0.900090 0.5098000 3.28E-05 0.020060 0.21300 0.3747 0.909090 0.5104000 8.39E-05 0.008195 0.20220 0.3744 1.000000 0.5199000 0.0002516 -0.116200 0.08290 0.3737 1.100100 0.5273000 0.0003908 -0.212300 -0.02900 0.3794 1.111100 0.5278000 0.0004074 -0.220700 -0.03875 0.3804 1.200500 0.5361000 0.0004479 -0.313300 -0.13380 0.3838 1.250000 0.5409000 0.0004860 -0.367900 -0.18910 0.3877 1.300400 0.5444000 0.0005329 -0.411300 -0.23770 0.3920 1.400600 0.5481000 0.0007676 -0.487000 -0.31550 0.3935 1.428600 0.5494000 0.0008272 -0.509500 -0.33850 0.3941 1.499300 0.5527000 0.0009124 -0.560400 -0.39350 0.3932 1.600000 0.5557000 0.0009844 -0.618600 -0.45830 0.3919 1.666700 0.5580000 0.0010850 -0.656400 -0.50260 0.3919 1.798600 0.5620000 0.0012450 -0.725800 -0.58340 0.3942 2.000000 0.5622000 0.0013750 -0.796300 -0.66600 0.4030 2.197800 0.5617000 0.0016520 -0.865600 -0.73950 0.4055 2.398100 0.5641000 0.0018290 -0.940600 -0.81580 0.4093 2.500000 0.5654000 0.0019210 -0.978700 -0.85420 0.4110 2.597400 0.5677000 0.0020060 -1.019000 -0.89800 0.4130 2.801100 0.5666000 0.0022770 -1.071000 -0.94950 0.4202 3.003000 0.5683000 0.0024490 -1.130000 -1.01400 0.4255 3.205100 0.5686000 0.0025360 -1.179000 -1.06900 0.4301 3.333300 0.5705000 0.0025330 -1.220000 -1.11100 0.4329 3.401400 0.5715000 0.0025410 -1.243000 -1.13500 0.4340 3.597100 0.5727000 0.0025730 -1.300000 -1.19400 0.4359 3.802300 0.5712000 0.0026620 -1.350000 -1.24200 0.4365 4.000000 0.5722000 0.0027110 -1.417000 -1.30300 0.4344 4.504500 0.5856000 0.0024490 -1.662000 -1.52000 0.4278 5.000000 0.5990000 0.0021050 -1.886000 -1.72900 0.4233 5.494500 0.6106000 0.0019410 -2.072000 -1.90400 0.4260 5.988000 0.6160000 0.0018800 -2.201000 -2.02800 0.4284 6.993000 0.6175000 0.0017660 -2.373000 -2.19200 0.4299 8.000000 0.6145000 0.0017210 -2.491000 -2.30800 0.4284 9.009000 0.6122000 0.0016370 -2.592000 -2.40800 0.4236 10.00000 0.6086000 0.0015630 -2.668000 -2.48500 0.4183 """)
[docs]class BergeThierryEtAl2003SIGMAMwW(BergeThierryEtAl2003SIGMA): """ Mw version of the Berge-Thierry et al. (2003) GMPE. For this conversion we use the Weatherill et al. (2016) conversion equation between Ms and Mw """
[docs] def get_mean_and_stddevs(self, sites, rup, dists, imt, stddev_types): newrup = copy.copy(rup) if rup.mag <= 6.064: newrup.mag = (rup.mag - 2.369) / 0.616 else: newrup.mag = (rup.mag - 0.100) / 0.994 return self._get_mean_and_stddevs(sites, newrup, dists, imt, stddev_types)