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