1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6import re
7import logging
8import numpy as num
9import hashlib
10import base64
12from pyrocko import util, moment_tensor
14from pyrocko.guts import Float, String, Timestamp, Unicode, \
15 StringPattern, List, Dict, Any
16from .location import Location
18logger = logging.getLogger('pyrocko.model.event')
20guts_prefix = 'pf'
22d2r = num.pi / 180.
25def cmp(a, b):
26 return (a > b) - (a < b)
29def ehash(s):
30 return str(base64.urlsafe_b64encode(
31 hashlib.sha1(s.encode('utf8')).digest()).decode('ascii'))
34def float_or_none_to_str(x, prec=9):
35 return 'None' if x is None else '{:.{prec}e}'.format(x, prec=prec)
38class FileParseError(Exception):
39 pass
42class EventExtrasDumpError(Exception):
43 pass
46class EOF(Exception):
47 pass
50class EmptyEvent(Exception):
51 pass
54class Tag(StringPattern):
55 pattern = r'^([A-Za-z][A-Za-z0-9._]{0,128})(:([A-Za-z0-9._-]*))?$'
58def opportunistic_cast(v):
59 try:
60 return int(v)
61 except ValueError:
62 pass
64 try:
65 return float(v)
66 except ValueError:
67 pass
69 return v
72class Event(Location):
73 '''
74 Representation of a seismic event.
76 :param lat: latitude of hypocenter (default 0.0)
77 :param lon: longitude of hypocenter (default 0.0)
78 :param time: origin time system timestamp
79 :param name: event identifier as string (optional)
80 :param depth: source depth (optional)
81 :param magnitude: magnitude of event (optional)
82 :param region: source region (optional)
83 :param catalog: name of catalog that lists this event (optional)
84 :param moment_tensor: moment tensor as
85 :py:class:`moment_tensor.MomentTensor` instance (optional)
86 :param duration: source duration as float (optional)
87 :param tags: list of tags describing event (optional)
88 :param extras: dictionary for user defined event attributes (optional).
89 Keys must be strings, values must be YAML serializable.
90 '''
92 time = Timestamp.T(default=Timestamp.D('1970-01-01 00:00:00'))
93 depth = Float.T(optional=True)
94 name = String.T(default='', optional=True, yamlstyle="'")
95 magnitude = Float.T(optional=True)
96 magnitude_type = String.T(optional=True, yamlstyle="'")
97 region = Unicode.T(optional=True, yamlstyle="'")
98 catalog = String.T(optional=True, yamlstyle="'")
99 moment_tensor = moment_tensor.MomentTensor.T(optional=True)
100 duration = Float.T(optional=True)
101 tags = List.T(Tag.T(), default=[])
102 extras = Dict.T(String.T(), Any.T(), default={})
104 def __init__(
105 self, lat=0., lon=0., north_shift=0., east_shift=0., time=0.,
106 name='', depth=None, elevation=None,
107 magnitude=None, magnitude_type=None, region=None, load=None,
108 loadf=None, catalog=None, moment_tensor=None, duration=None,
109 tags=None, extras=None):
111 if tags is None:
112 tags = []
114 if extras is None:
115 extras = {}
117 vals = None
118 if load is not None:
119 vals = Event.oldload(load)
120 elif loadf is not None:
121 vals = Event.oldloadf(loadf)
123 if vals:
124 lat, lon, north_shift, east_shift, time, name, depth, magnitude, \
125 magnitude_type, region, catalog, moment_tensor, duration, \
126 tags = vals
128 Location.__init__(
129 self, lat=lat, lon=lon,
130 north_shift=north_shift, east_shift=east_shift,
131 time=time, name=name, depth=depth,
132 elevation=elevation,
133 magnitude=magnitude, magnitude_type=magnitude_type,
134 region=region, catalog=catalog,
135 moment_tensor=moment_tensor, duration=duration, tags=tags,
136 extras=extras)
138 def tags_as_dict(self):
139 d = {}
140 for tag in self.tags:
141 m = re.match(Tag.pattern, tag)
142 if m:
143 k, v = m.group(1), opportunistic_cast(m.group(3))
144 d[k] = None if m.group(2) == '' else v
145 else:
146 logger.warning('Invalid event tag: %s' % tag)
148 return d
150 def time_as_string(self):
151 return util.time_to_str(self.time)
153 def set_name(self, name):
154 self.name = name
156 def olddump(self, filename):
157 file = open(filename, 'w')
158 self.olddumpf(file)
159 file.close()
161 def olddumpf(self, file):
162 if self.extras:
163 raise EventExtrasDumpError(
164 'Event user-defined extras attributes cannot be dumped in the '
165 '"basic" event file format. Use '
166 'dump_events(..., format="yaml").')
168 file.write('name = %s\n' % self.name)
169 file.write('time = %s\n' % util.time_to_str(self.time))
171 if self.lat != 0.0:
172 file.write('latitude = %.12g\n' % self.lat)
173 if self.lon != 0.0:
174 file.write('longitude = %.12g\n' % self.lon)
176 if self.north_shift != 0.0:
177 file.write('north_shift = %.12g\n' % self.north_shift)
178 if self.east_shift != 0.0:
179 file.write('east_shift = %.12g\n' % self.east_shift)
181 if self.magnitude is not None:
182 file.write('magnitude = %g\n' % self.magnitude)
183 file.write('moment = %g\n' %
184 moment_tensor.magnitude_to_moment(self.magnitude))
185 if self.magnitude_type is not None:
186 file.write('magnitude_type = %s\n' % self.magnitude_type)
187 if self.depth is not None:
188 file.write('depth = %.10g\n' % self.depth)
189 if self.region is not None:
190 file.write('region = %s\n' % self.region)
191 if self.catalog is not None:
192 file.write('catalog = %s\n' % self.catalog)
193 if self.moment_tensor is not None:
194 m = self.moment_tensor.m()
195 sdr1, sdr2 = self.moment_tensor.both_strike_dip_rake()
196 file.write((
197 'mnn = %g\nmee = %g\nmdd = %g\nmne = %g\nmnd = %g\nmed = %g\n'
198 'strike1 = %g\ndip1 = %g\nrake1 = %g\n'
199 'strike2 = %g\ndip2 = %g\nrake2 = %g\n') % (
200 (m[0, 0], m[1, 1], m[2, 2], m[0, 1], m[0, 2], m[1, 2]) +
201 sdr1 + sdr2))
203 if self.duration is not None:
204 file.write('duration = %g\n' % self.duration)
206 if self.tags:
207 file.write('tags = %s\n' % ', '.join(self.tags))
209 @staticmethod
210 def unique(events, deltat=10., group_cmp=(lambda a, b:
211 cmp(a.catalog, b.catalog))):
212 groups = Event.grouped(events, deltat)
214 events = []
215 for group in groups:
216 if group:
217 group.sort(group_cmp)
218 events.append(group[-1])
220 return events
222 @staticmethod
223 def grouped(events, deltat=10.):
224 events = list(events)
225 groups = []
226 for ia, a in enumerate(events):
227 groups.append([])
228 haveit = False
229 for ib, b in enumerate(events[:ia]):
230 if abs(b.time - a.time) < deltat:
231 groups[ib].append(a)
232 haveit = True
233 break
235 if not haveit:
236 groups[ia].append(a)
238 groups = [g for g in groups if g]
239 groups.sort(key=lambda g: sum(e.time for e in g) // len(g))
240 return groups
242 @staticmethod
243 def dump_catalog(events, filename=None, stream=None):
244 if filename is not None:
245 file = open(filename, 'w')
246 else:
247 file = stream
248 try:
249 i = 0
250 for ev in events:
252 ev.olddumpf(file)
254 file.write('--------------------------------------------\n')
255 i += 1
257 finally:
258 if filename is not None:
259 file.close()
261 @staticmethod
262 def oldload(filename):
263 with open(filename, 'r') as file:
264 return Event.oldloadf(file)
266 @staticmethod
267 def oldloadf(file):
268 d = {}
269 try:
270 for line in file:
271 if line.lstrip().startswith('#'):
272 continue
274 toks = line.split(' = ', 1)
275 if len(toks) == 2:
276 k, v = toks[0].strip(), toks[1].strip()
277 if k in ('name', 'region', 'catalog', 'magnitude_type'):
278 d[k] = v
279 if k in (('latitude longitude magnitude depth duration '
280 'north_shift east_shift '
281 'mnn mee mdd mne mnd med strike1 dip1 rake1 '
282 'strike2 dip2 rake2 duration').split()):
283 d[k] = float(v)
284 if k == 'time':
285 d[k] = util.str_to_time(v)
286 if k == 'tags':
287 d[k] = [x.strip() for x in v.split(',')]
289 if line.startswith('---'):
290 d['have_separator'] = True
291 break
293 except Exception as e:
294 raise FileParseError(e)
296 if not d:
297 raise EOF()
299 if 'have_separator' in d and len(d) == 1:
300 raise EmptyEvent()
302 mt = None
303 m6 = [d[x] for x in 'mnn mee mdd mne mnd med'.split() if x in d]
304 if len(m6) == 6:
305 mt = moment_tensor.MomentTensor(m=moment_tensor.symmat6(*m6))
306 else:
307 sdr = [d[x] for x in 'strike1 dip1 rake1'.split() if x in d]
308 if len(sdr) == 3:
309 moment = 1.0
310 if 'moment' in d:
311 moment = d['moment']
312 elif 'magnitude' in d:
313 moment = moment_tensor.magnitude_to_moment(d['magnitude'])
315 mt = moment_tensor.MomentTensor(
316 strike=sdr[0], dip=sdr[1], rake=sdr[2],
317 scalar_moment=moment)
319 return (
320 d.get('latitude', 0.0),
321 d.get('longitude', 0.0),
322 d.get('north_shift', 0.0),
323 d.get('east_shift', 0.0),
324 d.get('time', 0.0),
325 d.get('name', ''),
326 d.get('depth', None),
327 d.get('magnitude', None),
328 d.get('magnitude_type', None),
329 d.get('region', None),
330 d.get('catalog', None),
331 mt,
332 d.get('duration', None),
333 d.get('tags', []))
335 @staticmethod
336 def load_catalog(filename):
338 with open(filename, 'r') as file:
339 try:
340 while True:
341 try:
342 ev = Event(loadf=file)
343 yield ev
344 except EmptyEvent:
345 pass
347 except EOF:
348 pass
350 def get_hash(self):
351 e = self
352 if isinstance(e.time, float):
353 stime = util.time_to_str(e.time, format='%Y-%m-%d %H:%M:%S.3FRAC')
354 else:
355 stime = util.time_to_str(e.time, format='%Y-%m-%d %H:%M:%S.6FRAC')
357 s = float_or_none_to_str
359 to_hash = ', '.join((
360 stime,
361 s(e.lat), s(e.lon), s(e.depth),
362 float_or_none_to_str(e.magnitude, 5),
363 str(e.catalog), str(e.name or ''),
364 str(e.region)))
366 return ehash(to_hash)
368 def human_str(self):
369 s = [
370 'Latitude [deg]: %g' % self.lat,
371 'Longitude [deg]: %g' % self.lon,
372 'Time [UTC]: %s' % util.time_to_str(self.time)]
374 if self.name:
375 s.append('Name: %s' % self.name)
377 if self.depth is not None:
378 s.append('Depth [km]: %g' % (self.depth / 1000.))
380 if self.magnitude is not None:
381 s.append('Magnitude [%s]: %3.1f' % (
382 self.magnitude_type or 'M?', self.magnitude))
384 if self.region:
385 s.append('Region: %s' % self.region)
387 if self.catalog:
388 s.append('Catalog: %s' % self.catalog)
390 if self.moment_tensor:
391 s.append(str(self.moment_tensor))
393 return '\n'.join(s)
395 @property
396 def summary(self):
397 return '%s: %s, %s, %s, %s' % (
398 self.__class__.__name__,
399 self.name,
400 util.time_to_str(self.time),
401 '%-3s %3.1f' % (self.magnitude_type or ' ', self.magnitude)
402 if self.magnitude else 'M ---',
403 self.region)
406def detect_format(filename):
407 with open(filename, 'r') as f:
408 for line in f:
409 line = line.strip()
410 if not line or line.startswith('#') or line.startswith('%'):
411 continue
412 if line.startswith('--- !pf.Event'):
413 return 'yaml'
414 else:
415 return 'basic'
417 return 'basic'
420def load_events(filename, format='detect'):
421 '''
422 Read events file.
424 :param filename: name of file as str
425 :param format: file format: ``'detect'``, ``'basic'``, or ``'yaml'``
426 :returns: list of :py:class:`Event` objects
427 '''
429 if filename.startswith('http://') or filename.startswith('https://'):
430 import tempfile
431 with tempfile.NamedTemporaryFile() as fp:
432 util.download_file(filename, fp.name)
433 return load_events(fp.name, format=format)
435 if format == 'detect':
436 format = detect_format(filename)
438 if format == 'yaml':
439 from pyrocko import guts
440 events = [
441 ev for ev in guts.load_all(filename=filename)
442 if isinstance(ev, Event)]
444 return events
445 elif format == 'basic':
446 return list(Event.load_catalog(filename))
447 else:
448 from pyrocko.io.io_common import FileLoadError
449 raise FileLoadError('unknown event file format: %s' % format)
452class OneEventRequired(Exception):
453 pass
456def load_one_event(filename, format='detect'):
457 events = load_events(filename)
458 if len(events) != 1:
459 raise OneEventRequired(
460 'exactly one event is required in "%s"' % filename)
462 return events[0]
465def dump_events(events, filename=None, stream=None, format='basic'):
466 '''
467 Write events file.
469 :param events: list of :py:class:`Event` objects
470 :param filename: name of file as str
471 :param format: file format: ``'basic'``, or ``'yaml'``
472 '''
474 if format == 'basic':
475 Event.dump_catalog(events, filename=filename, stream=stream)
477 elif format == 'yaml':
478 from pyrocko import guts
479 events = [ev for ev in events if isinstance(ev, Event)]
480 guts.dump_all(object=events, filename=filename, stream=None)
482 else:
483 from pyrocko.io.io_common import FileSaveError
484 raise FileSaveError('unknown event file format: %s' % format)
487def load_kps_event_list(filename):
488 elist = []
489 f = open(filename, 'r')
490 for line in f:
491 toks = line.split()
492 if len(toks) < 7:
493 continue
495 tim = util.to_time_float(util.ctimegm(toks[0]+' '+toks[1]))
496 lat, lon, depth, magnitude = [float(x) for x in toks[2:6]]
497 duration = float(toks[10])
498 region = toks[-1]
499 name = util.gmctime_fn(tim)
500 e = Event(
501 lat, lon, tim,
502 name=name,
503 depth=depth,
504 magnitude=magnitude,
505 duration=duration,
506 region=region)
508 elist.append(e)
510 f.close()
511 return elist
514def load_gfz_event_list(filename):
515 from pyrocko import catalog
516 cat = catalog.Geofon()
518 elist = []
519 f = open(filename, 'r')
520 for line in f:
521 e = cat.get_event(line.strip())
522 elist.append(e)
524 f.close()
525 return elist