Coverage for /usr/local/lib/python3.11/dist-packages/pyrocko/model/event.py: 66%

294 statements  

« prev     ^ index     » next       coverage.py v6.5.0, created at 2023-10-06 15:01 +0000

1# http://pyrocko.org - GPLv3 

2# 

3# The Pyrocko Developers, 21st Century 

4# ---|P------/S----------~Lg---------- 

5 

6''' 

7Simple representation of a seismic event. 

8''' 

9 

10import re 

11import logging 

12import numpy as num 

13import hashlib 

14import base64 

15 

16from pyrocko import util, moment_tensor 

17 

18from pyrocko.guts import Float, String, Timestamp, Unicode, \ 

19 StringPattern, List, Dict, Any 

20from .location import Location 

21 

22logger = logging.getLogger('pyrocko.model.event') 

23 

24guts_prefix = 'pf' 

25 

26d2r = num.pi / 180. 

27 

28 

29def cmp(a, b): 

30 return (a > b) - (a < b) 

31 

32 

33def ehash(s): 

34 return str(base64.urlsafe_b64encode( 

35 hashlib.sha1(s.encode('utf8')).digest()).decode('ascii')) 

36 

37 

38def float_or_none_to_str(x, prec=9): 

39 return 'None' if x is None else '{:.{prec}e}'.format(x, prec=prec) 

40 

41 

42class FileParseError(Exception): 

43 pass 

44 

45 

46class EventExtrasDumpError(Exception): 

47 pass 

48 

49 

50class EOF(Exception): 

51 pass 

52 

53 

54class EmptyEvent(Exception): 

55 pass 

56 

57 

58class Tag(StringPattern): 

59 pattern = r'^([A-Za-z][A-Za-z0-9._]{0,128})(:([A-Za-z0-9._-]*))?$' 

60 

61 

62def opportunistic_cast(v): 

63 try: 

64 return int(v) 

65 except ValueError: 

66 pass 

67 

68 try: 

69 return float(v) 

70 except ValueError: 

71 pass 

72 

73 return v 

74 

75 

76class Event(Location): 

77 ''' 

78 Representation of a seismic event. 

79 ''' 

80 

81 time = Timestamp.T( 

82 default=Timestamp.D('1970-01-01 00:00:00'), 

83 help='Origin time (UTC system timestamp) [s].') 

84 depth = Float.T( 

85 optional=True, 

86 help='Depth below surface [m].') 

87 name = String.T( 

88 default='', 

89 optional=True, 

90 yamlstyle="'", 

91 help='Event identifier.') 

92 magnitude = Float.T( 

93 optional=True, 

94 help='Magnitude of the event.') 

95 magnitude_type = String.T( 

96 optional=True, 

97 yamlstyle="'", 

98 help='Magnitude type :py:gattr:`magnitude` is given in.') 

99 region = Unicode.T( 

100 optional=True, 

101 yamlstyle="'", 

102 help='Source region.') 

103 catalog = String.T( 

104 optional=True, 

105 yamlstyle="'", 

106 help='Name of catalog that lists this event.') 

107 moment_tensor = moment_tensor.MomentTensor.T( 

108 optional=True, 

109 help='Moment tensor of the event.') 

110 duration = Float.T( 

111 optional=True, 

112 help='Source duration [s].') 

113 tags = List.T( 

114 Tag.T(), 

115 default=[], 

116 help='Auxiliary tags.') 

117 extras = Dict.T( 

118 String.T(), 

119 Any.T(), 

120 default={}, 

121 help='Additional user defined event attributes. The given values must ' 

122 'be YAML-serializable.') 

123 

124 def __init__( 

125 self, lat=0., lon=0., north_shift=0., east_shift=0., time=0., 

126 name='', depth=None, elevation=None, 

127 magnitude=None, magnitude_type=None, region=None, load=None, 

128 loadf=None, catalog=None, moment_tensor=None, duration=None, 

129 tags=None, extras=None): 

130 

131 if tags is None: 

132 tags = [] 

133 

134 if extras is None: 

135 extras = {} 

136 

137 vals = None 

138 if load is not None: 

139 vals = Event.oldload(load) 

140 elif loadf is not None: 

141 vals = Event.oldloadf(loadf) 

142 

143 if vals: 

144 lat, lon, north_shift, east_shift, time, name, depth, magnitude, \ 

145 magnitude_type, region, catalog, moment_tensor, duration, \ 

146 tags = vals 

147 

148 Location.__init__( 

149 self, lat=lat, lon=lon, 

150 north_shift=north_shift, east_shift=east_shift, 

151 time=time, name=name, depth=depth, 

152 elevation=elevation, 

153 magnitude=magnitude, magnitude_type=magnitude_type, 

154 region=region, catalog=catalog, 

155 moment_tensor=moment_tensor, duration=duration, tags=tags, 

156 extras=extras) 

157 

158 def tags_as_dict(self): 

159 d = {} 

160 for tag in self.tags: 

161 m = re.match(Tag.pattern, tag) 

162 if m: 

163 k, v = m.group(1), opportunistic_cast(m.group(3)) 

164 d[k] = None if m.group(2) == '' else v 

165 else: 

166 logger.warning('Invalid event tag: %s' % tag) 

167 

168 return d 

169 

170 def time_as_string(self): 

171 return util.time_to_str(self.time) 

172 

173 def set_name(self, name): 

174 self.name = name 

175 

176 def olddump(self, filename): 

177 file = open(filename, 'w') 

178 self.olddumpf(file) 

179 file.close() 

180 

181 def olddumpf(self, file): 

182 if self.extras: 

183 raise EventExtrasDumpError( 

184 'Event user-defined extras attributes cannot be dumped in the ' 

185 '"basic" event file format. Use ' 

186 'dump_events(..., format="yaml").') 

187 

188 file.write('name = %s\n' % self.name) 

189 file.write('time = %s\n' % util.time_to_str(self.time)) 

190 

191 if self.lat != 0.0: 

192 file.write('latitude = %.12g\n' % self.lat) 

193 if self.lon != 0.0: 

194 file.write('longitude = %.12g\n' % self.lon) 

195 

196 if self.north_shift != 0.0: 

197 file.write('north_shift = %.12g\n' % self.north_shift) 

198 if self.east_shift != 0.0: 

199 file.write('east_shift = %.12g\n' % self.east_shift) 

200 

201 if self.magnitude is not None: 

202 file.write('magnitude = %g\n' % self.magnitude) 

203 file.write('moment = %g\n' % 

204 moment_tensor.magnitude_to_moment(self.magnitude)) 

205 if self.magnitude_type is not None: 

206 file.write('magnitude_type = %s\n' % self.magnitude_type) 

207 if self.depth is not None: 

208 file.write('depth = %.10g\n' % self.depth) 

209 if self.region is not None: 

210 file.write('region = %s\n' % self.region) 

211 if self.catalog is not None: 

212 file.write('catalog = %s\n' % self.catalog) 

213 if self.moment_tensor is not None: 

214 m = self.moment_tensor.m() 

215 sdr1, sdr2 = self.moment_tensor.both_strike_dip_rake() 

216 file.write(( 

217 'mnn = %g\nmee = %g\nmdd = %g\nmne = %g\nmnd = %g\nmed = %g\n' 

218 'strike1 = %g\ndip1 = %g\nrake1 = %g\n' 

219 'strike2 = %g\ndip2 = %g\nrake2 = %g\n') % ( 

220 (m[0, 0], m[1, 1], m[2, 2], m[0, 1], m[0, 2], m[1, 2]) + 

221 sdr1 + sdr2)) 

222 

223 if self.duration is not None: 

224 file.write('duration = %g\n' % self.duration) 

225 

226 if self.tags: 

227 file.write('tags = %s\n' % ', '.join(self.tags)) 

228 

229 @staticmethod 

230 def unique(events, deltat=10., group_cmp=(lambda a, b: 

231 cmp(a.catalog, b.catalog))): 

232 groups = Event.grouped(events, deltat) 

233 

234 events = [] 

235 for group in groups: 

236 if group: 

237 group.sort(group_cmp) 

238 events.append(group[-1]) 

239 

240 return events 

241 

242 @staticmethod 

243 def grouped(events, deltat=10.): 

244 events = list(events) 

245 groups = [] 

246 for ia, a in enumerate(events): 

247 groups.append([]) 

248 haveit = False 

249 for ib, b in enumerate(events[:ia]): 

250 if abs(b.time - a.time) < deltat: 

251 groups[ib].append(a) 

252 haveit = True 

253 break 

254 

255 if not haveit: 

256 groups[ia].append(a) 

257 

258 groups = [g for g in groups if g] 

259 groups.sort(key=lambda g: sum(e.time for e in g) // len(g)) 

260 return groups 

261 

262 @staticmethod 

263 def dump_catalog(events, filename=None, stream=None): 

264 if filename is not None: 

265 file = open(filename, 'w') 

266 else: 

267 file = stream 

268 try: 

269 i = 0 

270 for ev in events: 

271 

272 ev.olddumpf(file) 

273 

274 file.write('--------------------------------------------\n') 

275 i += 1 

276 

277 finally: 

278 if filename is not None: 

279 file.close() 

280 

281 @staticmethod 

282 def oldload(filename): 

283 with open(filename, 'r') as file: 

284 return Event.oldloadf(file) 

285 

286 @staticmethod 

287 def oldloadf(file): 

288 d = {} 

289 try: 

290 for line in file: 

291 if line.lstrip().startswith('#'): 

292 continue 

293 

294 toks = line.split(' = ', 1) 

295 if len(toks) == 2: 

296 k, v = toks[0].strip(), toks[1].strip() 

297 if k in ('name', 'region', 'catalog', 'magnitude_type'): 

298 d[k] = v 

299 if k in (('latitude longitude magnitude depth duration ' 

300 'north_shift east_shift ' 

301 'mnn mee mdd mne mnd med strike1 dip1 rake1 ' 

302 'strike2 dip2 rake2 duration').split()): 

303 d[k] = float(v) 

304 if k == 'time': 

305 d[k] = util.str_to_time(v) 

306 if k == 'tags': 

307 d[k] = [x.strip() for x in v.split(',')] 

308 

309 if line.startswith('---'): 

310 d['have_separator'] = True 

311 break 

312 

313 except Exception as e: 

314 raise FileParseError(e) 

315 

316 if not d: 

317 raise EOF() 

318 

319 if 'have_separator' in d and len(d) == 1: 

320 raise EmptyEvent() 

321 

322 mt = None 

323 m6 = [d[x] for x in 'mnn mee mdd mne mnd med'.split() if x in d] 

324 if len(m6) == 6: 

325 mt = moment_tensor.MomentTensor(m=moment_tensor.symmat6(*m6)) 

326 else: 

327 sdr = [d[x] for x in 'strike1 dip1 rake1'.split() if x in d] 

328 if len(sdr) == 3: 

329 moment = 1.0 

330 if 'moment' in d: 

331 moment = d['moment'] 

332 elif 'magnitude' in d: 

333 moment = moment_tensor.magnitude_to_moment(d['magnitude']) 

334 

335 mt = moment_tensor.MomentTensor( 

336 strike=sdr[0], dip=sdr[1], rake=sdr[2], 

337 scalar_moment=moment) 

338 

339 return ( 

340 d.get('latitude', 0.0), 

341 d.get('longitude', 0.0), 

342 d.get('north_shift', 0.0), 

343 d.get('east_shift', 0.0), 

344 d.get('time', 0.0), 

345 d.get('name', ''), 

346 d.get('depth', None), 

347 d.get('magnitude', None), 

348 d.get('magnitude_type', None), 

349 d.get('region', None), 

350 d.get('catalog', None), 

351 mt, 

352 d.get('duration', None), 

353 d.get('tags', [])) 

354 

355 @staticmethod 

356 def load_catalog(filename): 

357 

358 with open(filename, 'r') as file: 

359 try: 

360 while True: 

361 try: 

362 ev = Event(loadf=file) 

363 yield ev 

364 except EmptyEvent: 

365 pass 

366 

367 except EOF: 

368 pass 

369 

370 def get_hash(self): 

371 ''' 

372 Get a pseudo-unique hash over the main attributes of the event. 

373 

374 The following attributes are hashed: :py:gattr:`time`, 

375 :py:gattr:`~pyrocko.model.location.Location.lat`, 

376 :py:gattr:`~pyrocko.model.location.Location.lon`, :py:gattr:`depth`, 

377 :py:gattr:`magnitude`, :py:gattr:`catalog`, :py:gattr:`name`, 

378 :py:gattr:`region`. 

379 

380 :returns: 

381 URL-safe base64 encoded SHA1 hash. 

382 ''' 

383 e = self 

384 if isinstance(e.time, float): 

385 stime = util.time_to_str(e.time, format='%Y-%m-%d %H:%M:%S.3FRAC') 

386 else: 

387 stime = util.time_to_str(e.time, format='%Y-%m-%d %H:%M:%S.6FRAC') 

388 

389 s = float_or_none_to_str 

390 

391 to_hash = ', '.join(( 

392 stime, 

393 s(e.lat), s(e.lon), s(e.depth), 

394 float_or_none_to_str(e.magnitude, 5), 

395 str(e.catalog), str(e.name or ''), 

396 str(e.region))) 

397 

398 return ehash(to_hash) 

399 

400 def human_str(self): 

401 s = [ 

402 'Latitude [deg]: %g' % self.lat, 

403 'Longitude [deg]: %g' % self.lon, 

404 'Time [UTC]: %s' % util.time_to_str(self.time)] 

405 

406 if self.name: 

407 s.append('Name: %s' % self.name) 

408 

409 if self.depth is not None: 

410 s.append('Depth [km]: %g' % (self.depth / 1000.)) 

411 

412 if self.magnitude is not None: 

413 s.append('Magnitude [%s]: %3.1f' % ( 

414 self.magnitude_type or 'M?', self.magnitude)) 

415 

416 if self.region: 

417 s.append('Region: %s' % self.region) 

418 

419 if self.catalog: 

420 s.append('Catalog: %s' % self.catalog) 

421 

422 if self.moment_tensor: 

423 s.append(str(self.moment_tensor)) 

424 

425 return '\n'.join(s) 

426 

427 @property 

428 def summary(self): 

429 return '%s: %s, %s, %s, %s' % ( 

430 self.__class__.__name__, 

431 self.name, 

432 util.time_to_str(self.time), 

433 '%-3s %3.1f' % (self.magnitude_type or ' ', self.magnitude) 

434 if self.magnitude else 'M ---', 

435 self.region) 

436 

437 

438def detect_format(filename): 

439 with open(filename, 'r') as f: 

440 for line in f: 

441 line = line.strip() 

442 if not line or line.startswith('#') or line.startswith('%'): 

443 continue 

444 if line.startswith('--- !pf.Event'): 

445 return 'yaml' 

446 else: 

447 return 'basic' 

448 

449 return 'basic' 

450 

451 

452def load_events(filename, format='detect'): 

453 ''' 

454 Read events file. 

455 

456 :param filename: name of file as str 

457 :param format: file format: ``'detect'``, ``'basic'``, or ``'yaml'`` 

458 :returns: list of :py:class:`Event` objects 

459 ''' 

460 

461 if filename.startswith('http://') or filename.startswith('https://'): 

462 import tempfile 

463 with tempfile.NamedTemporaryFile() as fp: 

464 util.download_file(filename, fp.name) 

465 return load_events(fp.name, format=format) 

466 

467 if format == 'detect': 

468 format = detect_format(filename) 

469 

470 if format == 'yaml': 

471 from pyrocko import guts 

472 events = [ 

473 ev for ev in guts.load_all(filename=filename) 

474 if isinstance(ev, Event)] 

475 

476 return events 

477 elif format == 'basic': 

478 return list(Event.load_catalog(filename)) 

479 else: 

480 from pyrocko.io.io_common import FileLoadError 

481 raise FileLoadError('unknown event file format: %s' % format) 

482 

483 

484class OneEventRequired(Exception): 

485 pass 

486 

487 

488def load_one_event(filename, format='detect'): 

489 events = load_events(filename) 

490 if len(events) != 1: 

491 raise OneEventRequired( 

492 'exactly one event is required in "%s"' % filename) 

493 

494 return events[0] 

495 

496 

497def dump_events(events, filename=None, stream=None, format='basic'): 

498 ''' 

499 Write events file. 

500 

501 :param events: list of :py:class:`Event` objects 

502 :param filename: name of file as str 

503 :param format: file format: ``'basic'``, or ``'yaml'`` 

504 ''' 

505 

506 if format == 'basic': 

507 Event.dump_catalog(events, filename=filename, stream=stream) 

508 

509 elif format == 'yaml': 

510 from pyrocko import guts 

511 events = [ev for ev in events if isinstance(ev, Event)] 

512 guts.dump_all(object=events, filename=filename, stream=None) 

513 

514 else: 

515 from pyrocko.io.io_common import FileSaveError 

516 raise FileSaveError('unknown event file format: %s' % format) 

517 

518 

519def load_kps_event_list(filename): 

520 elist = [] 

521 f = open(filename, 'r') 

522 for line in f: 

523 toks = line.split() 

524 if len(toks) < 7: 

525 continue 

526 

527 tim = util.to_time_float(util.ctimegm(toks[0]+' '+toks[1])) 

528 lat, lon, depth, magnitude = [float(x) for x in toks[2:6]] 

529 duration = float(toks[10]) 

530 region = toks[-1] 

531 name = util.gmctime_fn(tim) 

532 e = Event( 

533 lat, lon, tim, 

534 name=name, 

535 depth=depth, 

536 magnitude=magnitude, 

537 duration=duration, 

538 region=region) 

539 

540 elist.append(e) 

541 

542 f.close() 

543 return elist 

544 

545 

546def load_gfz_event_list(filename): 

547 from pyrocko import catalog 

548 cat = catalog.Geofon() 

549 

550 elist = [] 

551 f = open(filename, 'r') 

552 for line in f: 

553 e = cat.get_event(line.strip()) 

554 elist.append(e) 

555 

556 f.close() 

557 return elist