1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
5from __future__ import absolute_import, division
7import time
8import re
9import logging
10import json
12from pyrocko import model, util
13from pyrocko.moment_tensor import MomentTensor, symmat6
14from .base_catalog import EarthquakeCatalog, NotFound
16from pyrocko.util import urlopen
18logger = logging.getLogger('pyrocko.client.geofon')
20km = 1000.
23class Geofon(EarthquakeCatalog):
24 '''
25 Access the GEOFON earthquake catalog.
26 '''
27 def __init__(self, get_moment_tensors=True):
28 self.events = {}
29 self._get_moment_tensors = get_moment_tensors
31 def flush(self):
32 self.events = {}
34 def iter_event_names(
35 self,
36 time_range=None,
37 nmax=1000,
38 magmin=None,
39 latmin=-90.,
40 latmax=90.,
41 lonmin=-180.,
42 lonmax=180.):
44 logger.debug('In Geofon.iter_event_names(...)')
46 dmin = time.strftime('%Y-%m-%d', time.gmtime(time_range[0]))
47 dmax = time.strftime('%Y-%m-%d', time.gmtime(time_range[1]+24*60*60))
49 if magmin is None:
50 magmin = ''
51 else:
52 magmin = '%g' % magmin
54 ipage = 1
55 while True:
56 url = ('http://geofon.gfz-potsdam.de/eqinfo/list.php?' + '&'.join([
57 'page=%i' % ipage,
58 'datemin=%s' % dmin,
59 'datemax=%s' % dmax,
60 'latmin=%g' % latmin,
61 'latmax=%g' % latmax,
62 'lonmin=%g' % lonmin,
63 'lonmax=%g' % lonmax,
64 'magmin=%s' % magmin,
65 'fmt=geojson',
66 'nmax=%i' % nmax]))
68 logger.debug('Opening URL: %s' % url)
69 page = urlopen(url).read()
70 logger.debug('Received page (%i bytes)' % len(page))
71 events = self._parse_events_page(page)
72 for ev in events:
73 if ev.moment_tensor is True:
74 ev.moment_tensor = self.get_mt(ev)
76 if not events:
77 break
79 for ev in events:
80 if time_range[0] <= ev.time and ev.time <= time_range[1]:
81 self.events[ev.name] = ev
82 yield ev.name
84 ipage += 1
86 def get_event(self, name):
87 logger.debug('In Geofon.get_event("%s")' % name)
89 if name not in self.events:
90 url = 'http://geofon.gfz-potsdam.de/eqinfo/event.php' \
91 '?id=%s&fmt=geojson' % name
92 logger.debug('Opening URL: %s' % url)
93 page = urlopen(url).read()
94 logger.debug('Received page (%i bytes)' % len(page))
96 try:
97 ev = self._parse_event_page(page)
98 self.events[name] = ev
100 except NotFound:
101 raise NotFound(url) # reraise with url
103 ev = self.events[name]
105 if ev.moment_tensor is True:
106 ev.moment_tensor = self.get_mt(ev)
108 return ev
110 def get_mt(self, ev):
111 syear = time.strftime('%Y', time.gmtime(ev.time))
112 url = 'http://geofon.gfz-potsdam.de/data/alerts/%s/%s/mt.txt' % (
113 syear, ev.name)
115 try:
116 logger.debug('Opening URL: %s' % url)
117 page = urlopen(url).read()
118 logger.debug('Received page (%i bytes)' % len(page))
119 except util.HTTPError:
120 logger.warning('No MT found for event "%s".' % ev.name)
121 return None
123 return self._parse_mt_page(page)
125 def _parse_mt_page(self, page):
126 d = {}
127 for k in 'Scale', 'Mrr', 'Mtt', 'Mpp', 'Mrt', 'Mrp', 'Mtp':
128 r = k.encode('ascii')+br'\s*=?\s*(\S+)'
129 m = re.search(r, page)
130 if m:
131 s = m.group(1).replace(b'10**', b'1e')
132 d[k.lower()] = float(s)
134 m = symmat6(*(d[x] for x in 'mrr mtt mpp mrt mrp mtp'.split()))
135 m *= d['scale']
136 mt = MomentTensor(m_up_south_east=m)
138 return mt
140 def _parse_events_page(self, page, limit=None):
141 j = json.loads(page.decode('utf-8'))
142 events = []
143 for ifeature, feature in enumerate(j['features']):
144 ev = self._json_feature_to_event(feature)
145 events.append(ev)
146 if limit and ifeature + 1 == limit:
147 break
149 return events
151 def _parse_event_page(self, page):
152 return self._parse_events_page(page, limit=1)[0]
154 def _json_feature_to_event(self, feature):
155 name = feature['id']
156 lon, lat, depth = feature['geometry']['coordinates']
157 depth *= 1000.
158 properties = feature['properties']
159 magnitude = properties['mag']
160 magnitude_type = properties['magType']
161 region = properties['place']
162 tevent = util.str_to_time(properties['time'].replace('T', ' '))
164 if ((properties.get('hasMT', 'no') == 'yes')
165 or properties['magType'] == 'Mw') and self._get_moment_tensors:
167 moment_tensor = True # flag for caller to query MT
168 else:
169 moment_tensor = None
171 status = properties['status'][:1]
172 tags = []
173 if status in 'AMC':
174 tags.append('geofon_status:%s' % status)
176 category = properties.get('evtype', '')
177 if re.match(r'^[a-zA-Z0-9]+$', category):
178 tags.append('geofon_category:%s' % category)
180 ev = model.Event(
181 lat=float(lat),
182 lon=float(lon),
183 time=tevent,
184 name=name,
185 depth=float(depth),
186 magnitude=float(magnitude),
187 magnitude_type=str(magnitude_type),
188 region=str(region),
189 moment_tensor=moment_tensor,
190 catalog='GEOFON',
191 tags=tags)
193 return ev