rupture_mechanism.py 7.93 KB
Newer Older
1
2
3
4
5
6
"""
Class to manage the distribution of rupture mechanisms
"""
from __future__ import annotations
import numpy as np
from enum import Enum
7
from typing import Optional, Dict, Tuple, List
8
9
10
11
12
13
14
15
16
17
18
19
from openquake.hazardlib.geo import NodalPlane
from openquake.hazardlib.pmf import PMF
from shakyground2 import valid


class DIP_RANGES(Enum):
    """
    For the "global" default rupture mechanism distribution the possible ranges
    of dip values depend on the style-of-faulting. The possible range of values
    for each style-of-faulting is provided and they are assumed to be uniformly
    distributed.
    """
Marius Kriegerowski's avatar
linting    
Marius Kriegerowski committed
20

21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
    R = (20.0, 45.0)
    SS = (75.0, 91.0)
    N = (45.0, 75.0)


class RuptureMechanism(object):
    """
    General class to handle the mechanism properties of the rupture, which include the
    strike, dip and rake. The mechanism properties largely constrol the dimensions of the 3D
    rupture plane that is needed for calculating the finite-source to site distances (e.g.
    Joyner-Boore distance, Rupture Distance, Rx, Ry0 etc.)

    for the majority of reported earthquakes a single rupture mechanism is not known unless
    accompanied by a 3D finite fault rupture model. In the absence of any information other
    than the hypocentre, a "global" distribution of strike, dip and rake values is assumed.
    In the case that a focal mechanism is available then the RuptureMechanism is described by
    two equiprobable planes.

    Attributes:
        mechanism: Distribution of rupture mechanisms (probability mass function) as instance
                   of :class:`openquake.hazardlib.pmf.PMF`

    """
Marius Kriegerowski's avatar
linting    
Marius Kriegerowski committed
44

45
46
47
48
49
50
51
52
53
54
55
56
57
58
    def __init__(self, mechanism: Optional[PMF] = None):
        if mechanism:
            assert isinstance(mechanism, PMF)
            self.mechanism = mechanism
        else:
            # Build the default distribution
            self.mechanism = self.build_mechanism_distribution()

    def __iter__(self):
        # Iterate over the mechanisms in the data set yielding the probability
        # of the nodal plane and the nodal plane itself
        for prob, nodal_plane in self.mechanism.data:
            yield prob, nodal_plane

59
60
61
62
    def __len__(self):
        # Lenth corresponds to the number of mechanisms in the distribution
        return len(self.mechanism.data)

63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
    def sample(self, nsamples: int) -> Tuple[np.ndarray, np.ndarray, np.ndarray]:
        """
        Samples the focal mechanism distribution to return "nsamples" strike,
        dip and rake values

        Args:
            nsamples: Number of sample
        Returns:
            Tuple of vectors of samples strikes, dips and rakes
        """
        assert nsamples > 0
        sample_planes = self.mechanism.sample(nsamples)
        strikes = np.empty(nsamples, dtype=float)
        dips = np.empty(nsamples, dtype=float)
        rakes = np.empty(nsamples, dtype=float)
        for i, sample_plane in enumerate(sample_planes):
            strikes[i] = sample_plane.strike
            dips[i] = sample_plane.dip
            rakes[i] = sample_plane.rake
        return strikes, dips, rakes

    @classmethod
Marius Kriegerowski's avatar
linting    
Marius Kriegerowski committed
85
86
87
88
89
90
    def from_strike_dip_rake(
        cls,
        strike: Optional[float] = None,
        dip: Optional[float] = None,
        rake: Optional[float] = None,
    ) -> RuptureMechanism:
91
92
93
94
95
96
97
98
99
        """
        Constructs the rupture mechanism distribution from a simple strike, dip
        and rake combination, permitting None values if undefined

        Args:
            strike: Strike of fault (in decimal degrees between 0 and 360)
            dip: Dip of fault (in decimal degrees between 0 and 90)
            rake: Rake of fault (in decimal degrees between -180 and 180)
        """
Marius Kriegerowski's avatar
linting    
Marius Kriegerowski committed
100
101
102
103
104
        return cls(
            cls.build_mechanism_distribution(
                valid.strike(strike), valid.dip(dip), valid.rake(rake)
            )
        )
105
106

    @classmethod
107
    def from_focal_mechanism(cls, focal_mechanism: Dict) -> RuptureMechanism:
108
        """
109
110
        Constructs the focal mechanism from an evenly pair of nodal planes such as that
        of a focal mechanism
111
112

        Args:
113
114
115
            focal_mechanism: Dictionary containing two :class:
            `openquake.hazardlib.geo.NodalPlane` objects, each labeled as "nodal_plane_1"
            and "nodal_plane_2"
116
        """
Marius Kriegerowski's avatar
linting    
Marius Kriegerowski committed
117
118
119
        return cls(
            PMF(
                [
120
121
                    (0.5, focal_mechanism["nodal_plane_1"]),
                    (0.5, focal_mechanism["nodal_plane_2"]),
Marius Kriegerowski's avatar
linting    
Marius Kriegerowski committed
122
123
124
                ]
            )
        )
125

126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
    @classmethod
    def from_nodal_planes(cls, nodal_planes: List, probabilities: List) -> RuptureMechanism:
        """
        Constructs the rupture mechanism distribution from a list of nodal planes and their
        associated probabilities

        Args:
            nodal_planes: Set of nodal planes as a list of dictionaries, eac containing strike,
                          dip and rake
            probabilities: List of probabilities of the nodal planes (must sum to 1)
        """
        assert len(nodal_planes) == len(
            probabilities
        ), "Number of nodal planes not equal to number of probabilities"
        assert np.isclose(
            sum(probabilities), 1.0
        ), "Probabilities do not sum to 1.0 (sum = %.6f)" % sum(probabilities)
        mechanism_distribution = []
        for prob, npl in zip(probabilities, nodal_planes):
            mechanism = valid.mechanism(npl["strike"], npl["dip"], npl["rake"])
            mechanism_distribution.append((prob, mechanism))
        return cls(PMF(mechanism_distribution))

149
    @staticmethod
Marius Kriegerowski's avatar
linting    
Marius Kriegerowski committed
150
151
152
153
154
    def build_mechanism_distribution(
        strike: Optional[float] = None,
        dip: Optional[float] = None,
        rake: Optional[float] = None,
    ) -> PMF:
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
        """
        Builds a mechanism distribution from a partial (or complete) characterisation of the
        rupture mechanism

        Args:
            strike: Strike of fault in decimal degrees between 0 and 360, or None
            dip: Dip of fault in decimal degrees between 0 and 90, or None
            rake: Rake of fault in decimal degrees between -180 and 180, or None

        Returns:
            Probability mass function of the mechanism distribution as an instance of the
            :class:`openquake.hazardlib.pmf.PMF`
        """
        if rake is None:
            # If rake is not defined then describe a complete range of rakes every 15 degrees
            # from -165.0 to 180.0
            rakes = np.arange(-165.0, 181.0, 15.0)
            weight_r = (1.0 / len(rakes)) * np.ones(len(rakes))
        else:
            rakes = [rake]
            weight_r = [1.0]

        if strike is None:
            # If strike is not defined then describe a complete range of strikes every 15
            # degrees between 0 and 360 (not included)
            strikes = np.arange(0.0, 359.0, 15.0)
            weight_s = (1.0 / len(strikes)) * np.ones(len(strikes))
        else:
            strikes = [strike]
            weight_s = [1.0]

        mechanisms = []
        for wght1, rake_i in zip(weight_r, rakes):
            if dip is None:
                # If dip is undefined then sample uniformly in the range
                # appropriate to the style of faulting
                if rake_i >= 45.0 and rake_i < 135.0:
                    # Reverse
                    dip_range = DIP_RANGES.R.value
                elif rake_i < -45.0 and rake_i >= -135.0:
                    # Normal
                    dip_range = DIP_RANGES.N.value
                else:
                    # Strike-slip
                    dip_range = DIP_RANGES.SS.value
                dips = np.arange(dip_range[0], dip_range[1], 5.0)
                weight_d = (1.0 / len(dips)) * np.ones(len(dips))
            else:
                dips = [dip]
                weight_d = [1.0]
            for wght2, dip_i in zip(weight_d, dips):
                for wght3, strike_i in zip(weight_s, strikes):
                    prob = wght1 * wght2 * wght3
                    mechanisms.append((prob, NodalPlane(strike_i, dip_i, rake_i)))
        return PMF(mechanisms)