1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6import re
7import logging
8import numpy as num
9from os import path as op
10from collections import OrderedDict
11import json
13from pyrocko import config, util
16def parse_3tup(s):
17 m = re.match(r'^\(?([^,]+),([^,]*),([^,]*)\)$', s)
18 if m:
19 return [float(m.group(1)) if m.group(1) else None for i in range(3)]
20 else:
21 return [None, None, None]
24logger = logging.getLogger('ActiveFaults')
27class Fault(object):
28 __fields__ = OrderedDict()
29 __slots__ = list(__fields__.keys())
31 def get_property(self, fault_obj, attr):
32 try:
33 values = float(fault_obj['properties'][attr][1:3])
34 except KeyError:
35 if attr == 'lower_seis_depth' or attr == 'upper_seis_depth':
36 values = 0
37 else:
38 values = -999
39 return values
41 def __init__(self, f):
42 nodes = f['geometry']['coordinates']
43 props = f['properties']
44 lons = [p[0] for p in nodes]
45 lats = [p[1] for p in nodes]
46 self.lon = lons
47 self.lat = lats
48 self.slip_type = props.get('slip_type', 'Unknown')
50 for attr, attr_type in self.__fields__.items():
51 if attr in props:
52 if props[attr] is None or props[attr] == '':
53 continue
55 if isinstance(props[attr], attr_type):
56 v = props[attr]
58 elif attr_type is float:
59 try:
60 v = parse_3tup(props[attr])[0]
61 except TypeError:
62 v = float(props[attr])
64 elif attr_type is int:
65 v = int(props[attr])
67 else:
68 v = None
70 setattr(self, attr, v)
72 def get_surface_line(self):
73 arr = num.empty((len(self.lat), 2))
74 for i in range(len(self.lat)):
75 arr[i, 0] = self.lat[i]
76 arr[i, 1] = self.lon[i]
78 return arr
80 def __str__(self):
81 d = {attr: getattr(self, attr) for attr in self.__fields__.keys()}
82 return '\n'.join(['%s: %s' % (attr, val) for attr, val in d.items()])
85class ActiveFault(Fault):
86 __fields__ = OrderedDict([
88 ('lat', list),
89 ('lon', list),
90 ('average_dip', float),
91 ('average_rake', float),
92 ('lower_seis_depth', float),
93 ('upper_seis_depth', float),
94 ('slip_type', str),
95 ])
98class ActiveFaults(object):
99 URL_GEM_ACTIVE_FAULTS = 'https://raw.githubusercontent.com/cossatot/gem-global-active-faults/master/geojson/gem_active_faults.geojson' # noqa
101 def __init__(self):
102 self.fname_active_faults = op.join(
103 config.config().fault_lines_dir, 'gem_active_faults.geojson')
105 if not op.exists(self.fname_active_faults):
106 self.download()
108 self.active_faults = []
109 self._load_faults(self.fname_active_faults, ActiveFault)
111 def _load_faults(self, fname, cls):
112 with open(fname, 'r') as f:
113 gj = json.load(f)
114 faults = gj['features']
115 for f in faults:
116 fault = cls(f)
117 self.active_faults.append(fault)
118 logger.debug('loaded %d fault', self.nactive_faults)
120 def download(self):
121 logger.info('Downloading GEM active faults database...')
122 util.download_file(self.URL_GEM_ACTIVE_FAULTS,
123 self.fname_active_faults)
125 @property
126 def nactive_faults(self):
127 return len(self.active_faults)
129 def nactive_faults_nodes(self):
130 return int(sum(len(f.lat) for f in self.active_faults))
132 def get_coords(self):
133 return num.array([f['coordinates'] for f in self.active_faults])
136if __name__ == '__main__':
137 logging.basicConfig(level=logging.DEBUG)
138 activefaults = ActiveFaults()