Source code for openquake.hmtk.faults.mfd.characteristic
# -*- coding: utf-8 -*-
# vim: tabstop=4 shiftwidth=4 softtabstop=4
#
# LICENSE
#
# Copyright (C) 2010-2023 GEM Foundation, G. Weatherill, M. Pagani,
# D. Monelli.
#
# The Hazard Modeller's Toolkit 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.
#
# You should have received a copy of the GNU Affero General Public License
# along with OpenQuake. If not, see <http://www.gnu.org/licenses/>
#
# DISCLAIMER
#
# The software Hazard Modeller's Toolkit (openquake.hmtk) provided herein
# is released as a prototype implementation on behalf of
# scientists and engineers working within the GEM Foundation (Global
# Earthquake Model).
#
# It is distributed for the purpose of open collaboration and in the
# hope that it will be useful to the scientific, engineering, disaster
# risk and software design communities.
#
# The software is NOT distributed as part of GEM's OpenQuake suite
# (https://www.globalquakemodel.org/tools-products) and must be considered as a
# separate entity. The software provided herein is designed and implemented
# by scientific staff. It is not developed to the design standards, nor
# subject to same level of critical review by professional software
# developers, as GEM's OpenQuake software suite.
#
# Feedback and contribution to the software is welcome, and can be
# directed to the hazard scientific staff of the GEM Model Facility
# (hazard@globalquakemodel.org).
#
# The Hazard Modeller's Toolkit (openquake.hmtk) is therefore distributed
# WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
# FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
# for more details.
#
# The GEM Foundation, and the authors of the software, assume no
# liability for use of the software.
"""
Module :mod: openquake.hmtk.faults.mfd.characteristic implements
:class:Characteristic the simple characteristic earthquake calculator
of recurrence.
"""
import numpy as np
from scipy.stats import truncnorm
from math import fabs
from openquake.hmtk.faults.mfd.base import _scale_moment, BaseMFDfromSlip
[docs]class Characteristic(BaseMFDfromSlip):
"""
Class to implement the characteristic earthquake model assuming a truncated
Gaussian distribution
:param str mfd_model:
Type of magnitude frequency distribution
:param float mfd_weight:
Weight of the mfd distribution (for subsequent logic tree processing)
:param float bin_width:
Width of the magnitude bin (rates are given for the centre point)
:param float mmin:
Minimum magnitude
:param float mmax:
Maximum magnitude
:param float mmax_sigma:
Uncertainty on maximum magnitude
:param float lower_bound:
Lower bound of Gaussian distribution (as number of standard deviations)
:param float upper_bound:
Upper bound of Gaussian distribution (as number of standard deviations)
:param float sigma:
Standard deviation (in magnitude units) of the Gaussian distribution
:param numpy.ndarray occurrence_rate:
Activity rates for magnitude in the range mmin to mmax in steps of
bin_width
"""
[docs] def setUp(self, mfd_conf):
"""
Input core configuration parameters as specified in the
configuration file
:param dict mfd_conf:
Configuration file containing the following attributes:
* 'Model_Weight' - Logic tree weight of model type (float)
* 'MFD_spacing' - Width of MFD bin (float)
* 'Minimum_Magnitude' - Minimum magnitude of activity rates (float)
* 'Maximum_Magnitude' - Characteristic magnituded (float)
(if not defined will use scaling relation)
* 'Maximum_Magnitude_Uncertainty' - Uncertainty on
maximum magnitude
(If not defined and the MSR has a sigma term then this will be
taken from sigma)
* 'Lower_Bound' - Lower bound in terms of number of sigma (float)
* 'Upper_Bound' - Upper bound in terms of number of sigma (float)
* 'Sigma' - Standard deviation (in magnitude units) of distribution
"""
self.mfd_model = "Characteristic"
self.mfd_weight = mfd_conf["Model_Weight"]
self.bin_width = mfd_conf["MFD_spacing"]
self.mmin = None
self.mmax = None
self.mmax_sigma = None
self.lower_bound = mfd_conf["Lower_Bound"]
self.upper_bound = mfd_conf["Upper_Bound"]
self.sigma = mfd_conf["Sigma"]
self.occurrence_rate = None
[docs] def get_mmax(self, mfd_conf, msr, rake, area):
"""
Gets the mmax for the fault - reading directly from the config file
or using the msr otherwise
:param dict mfd_config:
Configuration file (see setUp for paramters)
:param msr:
Instance of :class: nhlib.scalerel
:param float rake:
Rake of the fault (in range -180 to 180)
:param float area:
Area of the fault surface (km^2)
"""
if mfd_conf["Maximum_Magnitude"]:
self.mmax = mfd_conf["Maximum_Magnitude"]
else:
self.mmax = msr.get_median_mag(area, rake)
self.mmax_sigma = mfd_conf.get(
"Maximum_Magnitude_Uncertainty", None
) or msr.get_std_dev_mag(None, rake)
[docs] def get_mfd(self, slip, area, shear_modulus=30.0):
"""
Calculates activity rate on the fault
:param float slip:
Slip rate in mm/yr
:param fault_width:
Width of the fault (km)
:param float disp_length_ratio:
Displacement to length ratio (dimensionless)
:param float shear_modulus:
Shear modulus of the fault (GPa)
:returns:
* Minimum Magnitude (float)
* Bin width (float)
* Occurrence Rates (numpy.ndarray)
"""
# Working in Nm so convert: shear_modulus - GPa -> Nm
# area - km ** 2. -> m ** 2.
# slip - mm/yr -> m/yr
moment_rate = (
(shear_modulus * 1.0e9) * (area * 1.0e6) * (slip / 1000.0)
)
moment_mag = _scale_moment(self.mmax, in_nm=True)
characteristic_rate = moment_rate / moment_mag
if self.sigma and (fabs(self.sigma) > 1e-5):
self.mmin = self.mmax + (self.lower_bound * self.sigma)
mag_upper = self.mmax + (self.upper_bound * self.sigma)
mag_range = np.arange(
self.mmin, mag_upper + self.bin_width, self.bin_width
)
self.occurrence_rate = characteristic_rate * (
truncnorm.cdf(
mag_range + (self.bin_width / 2.0),
self.lower_bound,
self.upper_bound,
loc=self.mmax,
scale=self.sigma,
)
- truncnorm.cdf(
mag_range - (self.bin_width / 2.0),
self.lower_bound,
self.upper_bound,
loc=self.mmax,
scale=self.sigma,
)
)
else:
# Returns only a single rate
self.mmin = self.mmax
self.occurrence_rate = np.array([characteristic_rate], dtype=float)
return self.mmin, self.bin_width, self.occurrence_rate
