1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
5from __future__ import absolute_import, division
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._-]*)?$'
58class Event(Location):
59 '''
60 Representation of a seismic event.
62 :param lat: latitude of hypocenter (default 0.0)
63 :param lon: longitude of hypocenter (default 0.0)
64 :param time: origin time system timestamp
65 :param name: event identifier as string (optional)
66 :param depth: source depth (optional)
67 :param magnitude: magnitude of event (optional)
68 :param region: source region (optional)
69 :param catalog: name of catalog that lists this event (optional)
70 :param moment_tensor: moment tensor as
71 :py:class:`moment_tensor.MomentTensor` instance (optional)
72 :param duration: source duration as float (optional)
73 :param tags: list of tags describing event (optional)
74 :param extras: dictionary for user defined event attributes (optional).
75 Keys must be strings, values must be YAML serializable.
76 '''
78 time = Timestamp.T(default=Timestamp.D('1970-01-01 00:00:00'))
79 depth = Float.T(optional=True)
80 name = String.T(default='', optional=True, yamlstyle="'")
81 magnitude = Float.T(optional=True)
82 magnitude_type = String.T(optional=True, yamlstyle="'")
83 region = Unicode.T(optional=True, yamlstyle="'")
84 catalog = String.T(optional=True, yamlstyle="'")
85 moment_tensor = moment_tensor.MomentTensor.T(optional=True)
86 duration = Float.T(optional=True)
87 tags = List.T(Tag.T(), default=[])
88 extras = Dict.T(String.T(), Any.T(), default={})
90 def __init__(
91 self, lat=0., lon=0., north_shift=0., east_shift=0., time=0.,
92 name='', depth=None, elevation=None,
93 magnitude=None, magnitude_type=None, region=None, load=None,
94 loadf=None, catalog=None, moment_tensor=None, duration=None,
95 tags=None, extras=None):
97 if tags is None:
98 tags = []
100 if extras is None:
101 extras = {}
103 vals = None
104 if load is not None:
105 vals = Event.oldload(load)
106 elif loadf is not None:
107 vals = Event.oldloadf(loadf)
109 if vals:
110 lat, lon, north_shift, east_shift, time, name, depth, magnitude, \
111 magnitude_type, region, catalog, moment_tensor, duration, \
112 tags = vals
114 Location.__init__(
115 self, lat=lat, lon=lon,
116 north_shift=north_shift, east_shift=east_shift,
117 time=time, name=name, depth=depth,
118 elevation=elevation,
119 magnitude=magnitude, magnitude_type=magnitude_type,
120 region=region, catalog=catalog,
121 moment_tensor=moment_tensor, duration=duration, tags=tags,
122 extras=extras)
124 def time_as_string(self):
125 return util.time_to_str(self.time)
127 def set_name(self, name):
128 self.name = name
130 def olddump(self, filename):
131 file = open(filename, 'w')
132 self.olddumpf(file)
133 file.close()
135 def olddumpf(self, file):
136 if self.extras:
137 raise EventExtrasDumpError(
138 'Event user-defined extras attributes cannot be dumped in the '
139 '"basic" event file format. Use '
140 'dump_events(..., format="yaml").')
142 file.write('name = %s\n' % self.name)
143 file.write('time = %s\n' % util.time_to_str(self.time))
145 if self.lat != 0.0:
146 file.write('latitude = %.12g\n' % self.lat)
147 if self.lon != 0.0:
148 file.write('longitude = %.12g\n' % self.lon)
150 if self.north_shift != 0.0:
151 file.write('north_shift = %.12g\n' % self.north_shift)
152 if self.east_shift != 0.0:
153 file.write('east_shift = %.12g\n' % self.east_shift)
155 if self.magnitude is not None:
156 file.write('magnitude = %g\n' % self.magnitude)
157 file.write('moment = %g\n' %
158 moment_tensor.magnitude_to_moment(self.magnitude))
159 if self.magnitude_type is not None:
160 file.write('magnitude_type = %s\n' % self.magnitude_type)
161 if self.depth is not None:
162 file.write('depth = %.10g\n' % self.depth)
163 if self.region is not None:
164 file.write('region = %s\n' % self.region)
165 if self.catalog is not None:
166 file.write('catalog = %s\n' % self.catalog)
167 if self.moment_tensor is not None:
168 m = self.moment_tensor.m()
169 sdr1, sdr2 = self.moment_tensor.both_strike_dip_rake()
170 file.write((
171 'mnn = %g\nmee = %g\nmdd = %g\nmne = %g\nmnd = %g\nmed = %g\n'
172 'strike1 = %g\ndip1 = %g\nrake1 = %g\n'
173 'strike2 = %g\ndip2 = %g\nrake2 = %g\n') % (
174 (m[0, 0], m[1, 1], m[2, 2], m[0, 1], m[0, 2], m[1, 2]) +
175 sdr1 + sdr2))
177 if self.duration is not None:
178 file.write('duration = %g\n' % self.duration)
180 if self.tags:
181 file.write('tags = %s\n' % ', '.join(self.tags))
183 @staticmethod
184 def unique(events, deltat=10., group_cmp=(lambda a, b:
185 cmp(a.catalog, b.catalog))):
186 groups = Event.grouped(events, deltat)
188 events = []
189 for group in groups:
190 if group:
191 group.sort(group_cmp)
192 events.append(group[-1])
194 return events
196 @staticmethod
197 def grouped(events, deltat=10.):
198 events = list(events)
199 groups = []
200 for ia, a in enumerate(events):
201 groups.append([])
202 haveit = False
203 for ib, b in enumerate(events[:ia]):
204 if abs(b.time - a.time) < deltat:
205 groups[ib].append(a)
206 haveit = True
207 break
209 if not haveit:
210 groups[ia].append(a)
212 groups = [g for g in groups if g]
213 groups.sort(key=lambda g: sum(e.time for e in g) // len(g))
214 return groups
216 @staticmethod
217 def dump_catalog(events, filename=None, stream=None):
218 if filename is not None:
219 file = open(filename, 'w')
220 else:
221 file = stream
222 try:
223 i = 0
224 for ev in events:
226 ev.olddumpf(file)
228 file.write('--------------------------------------------\n')
229 i += 1
231 finally:
232 if filename is not None:
233 file.close()
235 @staticmethod
236 def oldload(filename):
237 with open(filename, 'r') as file:
238 return Event.oldloadf(file)
240 @staticmethod
241 def oldloadf(file):
242 d = {}
243 try:
244 for line in file:
245 if line.lstrip().startswith('#'):
246 continue
248 toks = line.split(' = ', 1)
249 if len(toks) == 2:
250 k, v = toks[0].strip(), toks[1].strip()
251 if k in ('name', 'region', 'catalog', 'magnitude_type'):
252 d[k] = v
253 if k in (('latitude longitude magnitude depth duration '
254 'north_shift east_shift '
255 'mnn mee mdd mne mnd med strike1 dip1 rake1 '
256 'strike2 dip2 rake2 duration').split()):
257 d[k] = float(v)
258 if k == 'time':
259 d[k] = util.str_to_time(v)
260 if k == 'tags':
261 d[k] = [x.strip() for x in v.split(',')]
263 if line.startswith('---'):
264 d['have_separator'] = True
265 break
267 except Exception as e:
268 raise FileParseError(e)
270 if not d:
271 raise EOF()
273 if 'have_separator' in d and len(d) == 1:
274 raise EmptyEvent()
276 mt = None
277 m6 = [d[x] for x in 'mnn mee mdd mne mnd med'.split() if x in d]
278 if len(m6) == 6:
279 mt = moment_tensor.MomentTensor(m=moment_tensor.symmat6(*m6))
280 else:
281 sdr = [d[x] for x in 'strike1 dip1 rake1'.split() if x in d]
282 if len(sdr) == 3:
283 moment = 1.0
284 if 'moment' in d:
285 moment = d['moment']
286 elif 'magnitude' in d:
287 moment = moment_tensor.magnitude_to_moment(d['magnitude'])
289 mt = moment_tensor.MomentTensor(
290 strike=sdr[0], dip=sdr[1], rake=sdr[2],
291 scalar_moment=moment)
293 return (
294 d.get('latitude', 0.0),
295 d.get('longitude', 0.0),
296 d.get('north_shift', 0.0),
297 d.get('east_shift', 0.0),
298 d.get('time', 0.0),
299 d.get('name', ''),
300 d.get('depth', None),
301 d.get('magnitude', None),
302 d.get('magnitude_type', None),
303 d.get('region', None),
304 d.get('catalog', None),
305 mt,
306 d.get('duration', None),
307 d.get('tags', []))
309 @staticmethod
310 def load_catalog(filename):
312 with open(filename, 'r') as file:
313 try:
314 while True:
315 try:
316 ev = Event(loadf=file)
317 yield ev
318 except EmptyEvent:
319 pass
321 except EOF:
322 pass
324 def get_hash(self):
325 e = self
326 if isinstance(e.time, float):
327 stime = util.time_to_str(e.time, format='%Y-%m-%d %H:%M:%S.3FRAC')
328 else:
329 stime = util.time_to_str(e.time, format='%Y-%m-%d %H:%M:%S.6FRAC')
331 s = float_or_none_to_str
333 to_hash = ', '.join((
334 stime,
335 s(e.lat), s(e.lon), s(e.depth),
336 float_or_none_to_str(e.magnitude, 5),
337 str(e.catalog), str(e.name or ''),
338 str(e.region)))
340 return ehash(to_hash)
342 def human_str(self):
343 s = [
344 'Latitude [deg]: %g' % self.lat,
345 'Longitude [deg]: %g' % self.lon,
346 'Time [UTC]: %s' % util.time_to_str(self.time)]
348 if self.name:
349 s.append('Name: %s' % self.name)
351 if self.depth is not None:
352 s.append('Depth [km]: %g' % (self.depth / 1000.))
354 if self.magnitude is not None:
355 s.append('Magnitude [%s]: %3.1f' % (
356 self.magnitude_type or 'M?', self.magnitude))
358 if self.region:
359 s.append('Region: %s' % self.region)
361 if self.catalog:
362 s.append('Catalog: %s' % self.catalog)
364 if self.moment_tensor:
365 s.append(str(self.moment_tensor))
367 return '\n'.join(s)
369 @property
370 def summary(self):
371 return '%s: %s, %s, %s, %s' % (
372 self.__class__.__name__,
373 self.name,
374 util.time_to_str(self.time),
375 '%-3s %3.1f' % (self.magnitude_type or ' ', self.magnitude)
376 if self.magnitude else 'M ---',
377 self.region)
380def detect_format(filename):
381 with open(filename, 'r') as f:
382 for line in f:
383 line = line.strip()
384 if not line or line.startswith('#') or line.startswith('%'):
385 continue
386 if line.startswith('--- !pf.Event'):
387 return 'yaml'
388 else:
389 return 'basic'
391 return 'basic'
394def load_events(filename, format='detect'):
395 '''
396 Read events file.
398 :param filename: name of file as str
399 :param format: file format: ``'detect'``, ``'basic'``, or ``'yaml'``
400 :returns: list of :py:class:`Event` objects
401 '''
403 if format == 'detect':
404 format = detect_format(filename)
406 if format == 'yaml':
407 from pyrocko import guts
408 events = [
409 ev for ev in guts.load_all(filename=filename)
410 if isinstance(ev, Event)]
412 return events
413 elif format == 'basic':
414 return list(Event.load_catalog(filename))
415 else:
416 from pyrocko.io.io_common import FileLoadError
417 raise FileLoadError('unknown event file format: %s' % format)
420class OneEventRequired(Exception):
421 pass
424def load_one_event(filename, format='detect'):
425 events = load_events(filename)
426 if len(events) != 1:
427 raise OneEventRequired(
428 'exactly one event is required in "%s"' % filename)
430 return events[0]
433def dump_events(events, filename=None, stream=None, format='basic'):
434 '''
435 Write events file.
437 :param events: list of :py:class:`Event` objects
438 :param filename: name of file as str
439 :param format: file format: ``'basic'``, or ``'yaml'``
440 '''
442 if format == 'basic':
443 Event.dump_catalog(events, filename=filename, stream=stream)
445 elif format == 'yaml':
446 from pyrocko import guts
447 events = [ev for ev in events if isinstance(ev, Event)]
448 guts.dump_all(object=events, filename=filename, stream=None)
450 else:
451 from pyrocko.io.io_common import FileSaveError
452 raise FileSaveError('unknown event file format: %s' % format)
455def load_kps_event_list(filename):
456 elist = []
457 f = open(filename, 'r')
458 for line in f:
459 toks = line.split()
460 if len(toks) < 7:
461 continue
463 tim = util.to_time_float(util.ctimegm(toks[0]+' '+toks[1]))
464 lat, lon, depth, magnitude = [float(x) for x in toks[2:6]]
465 duration = float(toks[10])
466 region = toks[-1]
467 name = util.gmctime_fn(tim)
468 e = Event(
469 lat, lon, tim,
470 name=name,
471 depth=depth,
472 magnitude=magnitude,
473 duration=duration,
474 region=region)
476 elist.append(e)
478 f.close()
479 return elist
482def load_gfz_event_list(filename):
483 from pyrocko import catalog
484 cat = catalog.Geofon()
486 elist = []
487 f = open(filename, 'r')
488 for line in f:
489 e = cat.get_event(line.strip())
490 elist.append(e)
492 f.close()
493 return elist