1# http://pyrocko.org - GPLv3
2#
3# The Pyrocko Developers, 21st Century
4# ---|P------/S----------~Lg----------
6import logging
7import numpy as num
8import hashlib
9import base64
11from pyrocko import util, moment_tensor
13from pyrocko.guts import Float, String, Timestamp, Unicode, \
14 StringPattern, List, Dict, Any
15from .location import Location
17logger = logging.getLogger('pyrocko.model.event')
19guts_prefix = 'pf'
21d2r = num.pi / 180.
24def cmp(a, b):
25 return (a > b) - (a < b)
28def ehash(s):
29 return str(base64.urlsafe_b64encode(
30 hashlib.sha1(s.encode('utf8')).digest()).decode('ascii'))
33def float_or_none_to_str(x, prec=9):
34 return 'None' if x is None else '{:.{prec}e}'.format(x, prec=prec)
37class FileParseError(Exception):
38 pass
41class EventExtrasDumpError(Exception):
42 pass
45class EOF(Exception):
46 pass
49class EmptyEvent(Exception):
50 pass
53class Tag(StringPattern):
54 pattern = r'^[A-Za-z][A-Za-z0-9._]{0,128}(:[A-Za-z0-9._-]*)?$'
57class Event(Location):
58 '''
59 Representation of a seismic event.
61 :param lat: latitude of hypocenter (default 0.0)
62 :param lon: longitude of hypocenter (default 0.0)
63 :param time: origin time system timestamp
64 :param name: event identifier as string (optional)
65 :param depth: source depth (optional)
66 :param magnitude: magnitude of event (optional)
67 :param region: source region (optional)
68 :param catalog: name of catalog that lists this event (optional)
69 :param moment_tensor: moment tensor as
70 :py:class:`moment_tensor.MomentTensor` instance (optional)
71 :param duration: source duration as float (optional)
72 :param tags: list of tags describing event (optional)
73 :param extras: dictionary for user defined event attributes (optional).
74 Keys must be strings, values must be YAML serializable.
75 '''
77 time = Timestamp.T(default=Timestamp.D('1970-01-01 00:00:00'))
78 depth = Float.T(optional=True)
79 name = String.T(default='', optional=True, yamlstyle="'")
80 magnitude = Float.T(optional=True)
81 magnitude_type = String.T(optional=True, yamlstyle="'")
82 region = Unicode.T(optional=True, yamlstyle="'")
83 catalog = String.T(optional=True, yamlstyle="'")
84 moment_tensor = moment_tensor.MomentTensor.T(optional=True)
85 duration = Float.T(optional=True)
86 tags = List.T(Tag.T(), default=[])
87 extras = Dict.T(String.T(), Any.T(), default={})
89 def __init__(
90 self, lat=0., lon=0., north_shift=0., east_shift=0., time=0.,
91 name='', depth=None, elevation=None,
92 magnitude=None, magnitude_type=None, region=None, load=None,
93 loadf=None, catalog=None, moment_tensor=None, duration=None,
94 tags=None, extras=None):
96 if tags is None:
97 tags = []
99 if extras is None:
100 extras = {}
102 vals = None
103 if load is not None:
104 vals = Event.oldload(load)
105 elif loadf is not None:
106 vals = Event.oldloadf(loadf)
108 if vals:
109 lat, lon, north_shift, east_shift, time, name, depth, magnitude, \
110 magnitude_type, region, catalog, moment_tensor, duration, \
111 tags = vals
113 Location.__init__(
114 self, lat=lat, lon=lon,
115 north_shift=north_shift, east_shift=east_shift,
116 time=time, name=name, depth=depth,
117 elevation=elevation,
118 magnitude=magnitude, magnitude_type=magnitude_type,
119 region=region, catalog=catalog,
120 moment_tensor=moment_tensor, duration=duration, tags=tags,
121 extras=extras)
123 def time_as_string(self):
124 return util.time_to_str(self.time)
126 def set_name(self, name):
127 self.name = name
129 def olddump(self, filename):
130 file = open(filename, 'w')
131 self.olddumpf(file)
132 file.close()
134 def olddumpf(self, file):
135 if self.extras:
136 raise EventExtrasDumpError(
137 'Event user-defined extras attributes cannot be dumped in the '
138 '"basic" event file format. Use '
139 'dump_events(..., format="yaml").')
141 file.write('name = %s\n' % self.name)
142 file.write('time = %s\n' % util.time_to_str(self.time))
144 if self.lat != 0.0:
145 file.write('latitude = %.12g\n' % self.lat)
146 if self.lon != 0.0:
147 file.write('longitude = %.12g\n' % self.lon)
149 if self.north_shift != 0.0:
150 file.write('north_shift = %.12g\n' % self.north_shift)
151 if self.east_shift != 0.0:
152 file.write('east_shift = %.12g\n' % self.east_shift)
154 if self.magnitude is not None:
155 file.write('magnitude = %g\n' % self.magnitude)
156 file.write('moment = %g\n' %
157 moment_tensor.magnitude_to_moment(self.magnitude))
158 if self.magnitude_type is not None:
159 file.write('magnitude_type = %s\n' % self.magnitude_type)
160 if self.depth is not None:
161 file.write('depth = %.10g\n' % self.depth)
162 if self.region is not None:
163 file.write('region = %s\n' % self.region)
164 if self.catalog is not None:
165 file.write('catalog = %s\n' % self.catalog)
166 if self.moment_tensor is not None:
167 m = self.moment_tensor.m()
168 sdr1, sdr2 = self.moment_tensor.both_strike_dip_rake()
169 file.write((
170 'mnn = %g\nmee = %g\nmdd = %g\nmne = %g\nmnd = %g\nmed = %g\n'
171 'strike1 = %g\ndip1 = %g\nrake1 = %g\n'
172 'strike2 = %g\ndip2 = %g\nrake2 = %g\n') % (
173 (m[0, 0], m[1, 1], m[2, 2], m[0, 1], m[0, 2], m[1, 2]) +
174 sdr1 + sdr2))
176 if self.duration is not None:
177 file.write('duration = %g\n' % self.duration)
179 if self.tags:
180 file.write('tags = %s\n' % ', '.join(self.tags))
182 @staticmethod
183 def unique(events, deltat=10., group_cmp=(lambda a, b:
184 cmp(a.catalog, b.catalog))):
185 groups = Event.grouped(events, deltat)
187 events = []
188 for group in groups:
189 if group:
190 group.sort(group_cmp)
191 events.append(group[-1])
193 return events
195 @staticmethod
196 def grouped(events, deltat=10.):
197 events = list(events)
198 groups = []
199 for ia, a in enumerate(events):
200 groups.append([])
201 haveit = False
202 for ib, b in enumerate(events[:ia]):
203 if abs(b.time - a.time) < deltat:
204 groups[ib].append(a)
205 haveit = True
206 break
208 if not haveit:
209 groups[ia].append(a)
211 groups = [g for g in groups if g]
212 groups.sort(key=lambda g: sum(e.time for e in g) // len(g))
213 return groups
215 @staticmethod
216 def dump_catalog(events, filename=None, stream=None):
217 if filename is not None:
218 file = open(filename, 'w')
219 else:
220 file = stream
221 try:
222 i = 0
223 for ev in events:
225 ev.olddumpf(file)
227 file.write('--------------------------------------------\n')
228 i += 1
230 finally:
231 if filename is not None:
232 file.close()
234 @staticmethod
235 def oldload(filename):
236 with open(filename, 'r') as file:
237 return Event.oldloadf(file)
239 @staticmethod
240 def oldloadf(file):
241 d = {}
242 try:
243 for line in file:
244 if line.lstrip().startswith('#'):
245 continue
247 toks = line.split(' = ', 1)
248 if len(toks) == 2:
249 k, v = toks[0].strip(), toks[1].strip()
250 if k in ('name', 'region', 'catalog', 'magnitude_type'):
251 d[k] = v
252 if k in (('latitude longitude magnitude depth duration '
253 'north_shift east_shift '
254 'mnn mee mdd mne mnd med strike1 dip1 rake1 '
255 'strike2 dip2 rake2 duration').split()):
256 d[k] = float(v)
257 if k == 'time':
258 d[k] = util.str_to_time(v)
259 if k == 'tags':
260 d[k] = [x.strip() for x in v.split(',')]
262 if line.startswith('---'):
263 d['have_separator'] = True
264 break
266 except Exception as e:
267 raise FileParseError(e)
269 if not d:
270 raise EOF()
272 if 'have_separator' in d and len(d) == 1:
273 raise EmptyEvent()
275 mt = None
276 m6 = [d[x] for x in 'mnn mee mdd mne mnd med'.split() if x in d]
277 if len(m6) == 6:
278 mt = moment_tensor.MomentTensor(m=moment_tensor.symmat6(*m6))
279 else:
280 sdr = [d[x] for x in 'strike1 dip1 rake1'.split() if x in d]
281 if len(sdr) == 3:
282 moment = 1.0
283 if 'moment' in d:
284 moment = d['moment']
285 elif 'magnitude' in d:
286 moment = moment_tensor.magnitude_to_moment(d['magnitude'])
288 mt = moment_tensor.MomentTensor(
289 strike=sdr[0], dip=sdr[1], rake=sdr[2],
290 scalar_moment=moment)
292 return (
293 d.get('latitude', 0.0),
294 d.get('longitude', 0.0),
295 d.get('north_shift', 0.0),
296 d.get('east_shift', 0.0),
297 d.get('time', 0.0),
298 d.get('name', ''),
299 d.get('depth', None),
300 d.get('magnitude', None),
301 d.get('magnitude_type', None),
302 d.get('region', None),
303 d.get('catalog', None),
304 mt,
305 d.get('duration', None),
306 d.get('tags', []))
308 @staticmethod
309 def load_catalog(filename):
311 with open(filename, 'r') as file:
312 try:
313 while True:
314 try:
315 ev = Event(loadf=file)
316 yield ev
317 except EmptyEvent:
318 pass
320 except EOF:
321 pass
323 def get_hash(self):
324 e = self
325 if isinstance(e.time, float):
326 stime = util.time_to_str(e.time, format='%Y-%m-%d %H:%M:%S.3FRAC')
327 else:
328 stime = util.time_to_str(e.time, format='%Y-%m-%d %H:%M:%S.6FRAC')
330 s = float_or_none_to_str
332 to_hash = ', '.join((
333 stime,
334 s(e.lat), s(e.lon), s(e.depth),
335 float_or_none_to_str(e.magnitude, 5),
336 str(e.catalog), str(e.name or ''),
337 str(e.region)))
339 return ehash(to_hash)
341 def human_str(self):
342 s = [
343 'Latitude [deg]: %g' % self.lat,
344 'Longitude [deg]: %g' % self.lon,
345 'Time [UTC]: %s' % util.time_to_str(self.time)]
347 if self.name:
348 s.append('Name: %s' % self.name)
350 if self.depth is not None:
351 s.append('Depth [km]: %g' % (self.depth / 1000.))
353 if self.magnitude is not None:
354 s.append('Magnitude [%s]: %3.1f' % (
355 self.magnitude_type or 'M?', self.magnitude))
357 if self.region:
358 s.append('Region: %s' % self.region)
360 if self.catalog:
361 s.append('Catalog: %s' % self.catalog)
363 if self.moment_tensor:
364 s.append(str(self.moment_tensor))
366 return '\n'.join(s)
368 @property
369 def summary(self):
370 return '%s: %s, %s, %s, %s' % (
371 self.__class__.__name__,
372 self.name,
373 util.time_to_str(self.time),
374 '%-3s %3.1f' % (self.magnitude_type or ' ', self.magnitude)
375 if self.magnitude else 'M ---',
376 self.region)
379def detect_format(filename):
380 with open(filename, 'r') as f:
381 for line in f:
382 line = line.strip()
383 if not line or line.startswith('#') or line.startswith('%'):
384 continue
385 if line.startswith('--- !pf.Event'):
386 return 'yaml'
387 else:
388 return 'basic'
390 return 'basic'
393def load_events(filename, format='detect'):
394 '''
395 Read events file.
397 :param filename: name of file as str
398 :param format: file format: ``'detect'``, ``'basic'``, or ``'yaml'``
399 :returns: list of :py:class:`Event` objects
400 '''
402 if format == 'detect':
403 format = detect_format(filename)
405 if format == 'yaml':
406 from pyrocko import guts
407 events = [
408 ev for ev in guts.load_all(filename=filename)
409 if isinstance(ev, Event)]
411 return events
412 elif format == 'basic':
413 return list(Event.load_catalog(filename))
414 else:
415 from pyrocko.io.io_common import FileLoadError
416 raise FileLoadError('unknown event file format: %s' % format)
419class OneEventRequired(Exception):
420 pass
423def load_one_event(filename, format='detect'):
424 events = load_events(filename)
425 if len(events) != 1:
426 raise OneEventRequired(
427 'exactly one event is required in "%s"' % filename)
429 return events[0]
432def dump_events(events, filename=None, stream=None, format='basic'):
433 '''
434 Write events file.
436 :param events: list of :py:class:`Event` objects
437 :param filename: name of file as str
438 :param format: file format: ``'basic'``, or ``'yaml'``
439 '''
441 if format == 'basic':
442 Event.dump_catalog(events, filename=filename, stream=stream)
444 elif format == 'yaml':
445 from pyrocko import guts
446 events = [ev for ev in events if isinstance(ev, Event)]
447 guts.dump_all(object=events, filename=filename, stream=None)
449 else:
450 from pyrocko.io.io_common import FileSaveError
451 raise FileSaveError('unknown event file format: %s' % format)
454def load_kps_event_list(filename):
455 elist = []
456 f = open(filename, 'r')
457 for line in f:
458 toks = line.split()
459 if len(toks) < 7:
460 continue
462 tim = util.to_time_float(util.ctimegm(toks[0]+' '+toks[1]))
463 lat, lon, depth, magnitude = [float(x) for x in toks[2:6]]
464 duration = float(toks[10])
465 region = toks[-1]
466 name = util.gmctime_fn(tim)
467 e = Event(
468 lat, lon, tim,
469 name=name,
470 depth=depth,
471 magnitude=magnitude,
472 duration=duration,
473 region=region)
475 elist.append(e)
477 f.close()
478 return elist
481def load_gfz_event_list(filename):
482 from pyrocko import catalog
483 cat = catalog.Geofon()
485 elist = []
486 f = open(filename, 'r')
487 for line in f:
488 e = cat.get_event(line.strip())
489 elist.append(e)
491 f.close()
492 return elist