1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
5from __future__ import absolute_import, division, print_function
7import numpy as num
9from pyrocko.guts import Float, Int
10from pyrocko import moment_tensor, gf, util
12from .base import SourceGenerator
14km = 1e3
15guts_prefix = 'pf.scenario'
18class RectangularSourceGenerator(SourceGenerator):
19 depth_min = Float.T(
20 default=0.0)
21 depth_max = Float.T(
22 default=5*km)
23 decimation_factor = Int.T(
24 default=4)
26 strike = Float.T(
27 optional=True)
28 dip = Float.T(
29 optional=True)
30 rake = Float.T(
31 optional=True)
32 depth = Float.T(
33 optional=True)
35 width = Float.T(
36 optional=True)
37 length = Float.T(
38 optional=True)
40 def get_source(self, ievent):
41 rstate = self.get_rstate(ievent)
42 time = self.time_min + rstate.uniform(
43 0., float(self.time_max - self.time_min)) # hptime aware
44 lat, lon = self.get_latlon(ievent)
45 depth = rstate.uniform(self.depth_min, self.depth_max)
47 magnitude = self.draw_magnitude(rstate)
48 moment = moment_tensor.magnitude_to_moment(magnitude)
50 # After Mai and Beroza (2000)
51 length = num.exp(-6.27 + 0.4*num.log(moment))
52 width = num.exp(-4.24 + 0.32*num.log(moment))
54 length = length if not self.length else self.length
55 width = width if not self.width else self.width
56 depth = depth if not self.depth else self.depth
58 if self.strike is None and self.dip is None and self.rake is None:
59 strike, rake = rstate.uniform(-180., 180., 2)
60 dip = rstate.uniform(0., 90.)
61 else:
62 if None in (self.strike, self.dip, self.rake):
63 raise ValueError(
64 'RectangularFaultGenerator: '
65 'strike, dip, rake must be used in combination.')
67 strike = self.strike
68 dip = self.dip
69 rake = self.rake
71 source = gf.RectangularSource(
72 time=util.to_time_float(time),
73 lat=float(lat),
74 lon=float(lon),
75 anchor='top',
76 depth=float(depth),
77 length=float(length),
78 width=float(width),
79 strike=float(strike),
80 dip=float(dip),
81 rake=float(rake),
82 magnitude=magnitude,
83 decimation_factor=self.decimation_factor)
85 return source
87 def add_map_artists(self, automap):
88 for source in self.get_sources():
89 automap.gmt.psxy(
90 in_rows=source.outline(cs='lonlat'),
91 L='+p2p,black',
92 W='1p,black',
93 G='black',
94 t=50,
95 *automap.jxyr)