Coverage for /usr/local/lib/python3.11/dist-packages/pyrocko/gui/snuffler/snufflings/seismosizer.py: 48%

1# https://pyrocko.org - GPLv3

3# The Pyrocko Developers, 21st Century

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

6import math

7import numpy as num

8import os

10from pyrocko import moment_tensor, model

11from ..snuffling import Snuffling, Param, Choice, Switch, EventMarker

12from pyrocko import gf

14km = 1000.

17class Seismosizer(Snuffling):

19 def help(self):

20 return '''

21Generate synthetic traces on the fly

22====================================

24Activate an event (press `e`) to generate synthetic waveforms for it.

25If no stations have been loaded pripor to execution, two stations will be

26generated at lat/lon = (5., 0) and (-5., 0.).

28All geometrical units are kilometers (if not stated otherwise).

29`GF Stores` will be loaded on start from `gf_store_superdirs` defined

30in your pyrocko config file located at `$HOME/.pyrocko/config.pf`.

31'''

33 def __init__(self):

34 Snuffling.__init__(self)

35 self.stf_types = ['half sin', 'triangular', 'boxcar', 'None']

36 self.stf_instances = [gf.HalfSinusoidSTF(), gf.TriangularSTF(),

37 gf.BoxcarSTF(), None]

39 def get_source(self, event):

40 raise NotImplementedError()

42 def panel_visibility_changed(self, bool):

43 if bool:

44 if self._engine is None:

45 self.set_engine()

47 def set_engine(self):

48 self._engine = None

49 self.store_ids = self.get_store_ids()

50 if self.store_ids == []:

51 return

53 self.set_parameter_choices('store_id', self.store_ids)

54 self.store_id = self.store_ids[0]

56 def get_engine(self):

57 if not self._engine:

58 self._engine = gf.LocalEngine(use_config=True)

60 return self._engine

62 def get_store_ids(self):

63 return self.get_engine().get_store_ids()

65 def get_stf(self):

66 stf = dict(zip(self.stf_types, self.stf_instances))[self.stf_type]

67 if stf is not None:

68 stf.duration = self.stf_duration

70 return stf

72 def setup(self):

73 self.add_parameter(

74 Choice('GF Store', 'store_id',

75 '<not loaded yet>', ['<not loaded yet>']))

76 self.add_parameter(

77 Choice('Waveform type', 'waveform_type', 'Displacement [m]',

78 ['Displacement [m]',

79 'Displacement [nm]',

80 'Velocity [m/s]',

81 'Velocity [nm/s]',

82 'Acceleration [m/s^2]',

83 'Acceleration [nm/s^2]']))

85 self.add_trigger('Set Engine', self.set_engine)

86 self.add_trigger('Set Params from Event', self.mechanism_from_event)

87 self.add_trigger('Add Stores', self.add_store)

89 self.store_ids = None

90 self.offline_config = None

91 self._engine = None

93 def call(self):

94 '''

95 Main work routine of the snuffling.

96 '''

97 self.cleanup()

99 # get time range visible in viewer

100 viewer = self.get_viewer()

101

102 event = viewer.get_active_event()

103 if event:

104 event, stations = self.get_active_event_and_stations(

105 missing='warn')

106 else:

107 # event = model.Event(lat=self.lat, lon=self.lon)

108 event = model.Event(lat=0., lon=0.)

109 stations = []

110

111 stations = self.get_stations()

112

113 s2c = {}

114 for traces in self.chopper_selected_traces(fallback=True,

115 mode='visible'):

116 for tr in traces:

117 net, sta, loc, cha = tr.nslc_id

118 ns = net, sta

119 if ns not in s2c:

120 s2c[ns] = set()

121

122 s2c[ns].add((loc, cha))

123

124 if not stations:

125 stations = []

126 for (lat, lon) in [(5., 0.), (-5., 0.)]:

127 s = model.Station(station='(%g, %g)' % (lat, lon),

128 lat=lat, lon=lon)

129 stations.append(s)

130 viewer.add_stations(stations)

131

132 for s in stations:

133 ns = s.nsl()[:2]

134 if ns not in s2c:

135 s2c[ns] = set()

136

137 for cha in 'NEZ':

138 s2c[ns].add(('', cha))

139

140 source = self.get_source(event)

141 source.regularize()

142

143 m = EventMarker(source.pyrocko_event())

144 self.add_marker(m)

145

146 targets = []

147

148 if self.store_id == '<not loaded yet>':

149 self.fail('Select a GF Store first')

150

151 for station in stations:

152

153 nsl = station.nsl()

154 if nsl[:2] not in s2c:

155 continue

156

157 for loc, cha in s2c[nsl[:2]]:

158

159 target = gf.Target(

160 codes=(

161 station.network,

162 station.station,

163 loc + '-syn',

164 cha),

165 quantity='displacement',

166 lat=station.lat,

167 lon=station.lon,

168 depth=station.depth,

169 store_id=self.store_id,

170 optimization='enable',

171 interpolation='nearest_neighbor')

172

173 _, bazi = source.azibazi_to(target)

174

175 if cha.endswith('T'):

176 dip = 0.

177 azi = bazi + 270.

178 elif cha.endswith('R'):

179 dip = 0.

180 azi = bazi + 180.

181 elif cha.endswith('1'):

182 dip = 0.

183 azi = 0.

184 elif cha.endswith('2'):

185 dip = 0.

186 azi = 90.

187 else:

188 dip = None

189 azi = None

190

191 target.azimuth = azi

192 target.dip = dip

193

194 targets.append(target)

195

196 req = gf.Request(

197 sources=[source],

198 targets=targets)

199

200 req.regularize()

201

202 try:

203 resp = self.get_engine().process(req, nthreads=0)

204 except (gf.meta.OutOfBounds, gf.store_ext.StoreExtError)as e:

205 self.fail(e)

206

207 traces = resp.pyrocko_traces()

208

209 if self.waveform_type.startswith('Velocity'):

210 for tr in traces:

211 tr.set_ydata(num.diff(tr.ydata) / tr.deltat)

212

213 elif self.waveform_type.startswith('Acceleration'):

214 for tr in traces:

215 tr.set_ydata(num.diff(num.diff(tr.ydata)) / tr.deltat**2)

216

217 if self.waveform_type.endswith('[nm]') or \

218 self.waveform_type.endswith('[nm/s]') or \

219 self.waveform_type.endswith('[nm/s^2]'):

220

221 for tr in traces:

222 tr.set_ydata(tr.ydata * 1e9)

223

224 self.add_traces(traces)

225

226 def mechanism_from_event(self):

227

228 event = self.get_viewer().get_active_event()

229

230 if event is None:

231 self.fail('No active event set.')

232

233 if event.moment_tensor is not None:

234 strike, dip, slip_rake = event.moment_tensor\

235 .both_strike_dip_rake()[0]

236 moment = event.moment_tensor.scalar_moment()

237 self.set_parameter('magnitude',

238 moment_tensor.moment_to_magnitude(moment))

239 self.set_parameter('strike', strike)

240 self.set_parameter('dip', dip)

241 self.set_parameter('rake', slip_rake)

242 m6 = event.moment_tensor.m6()

243 m9 = moment_tensor.symmat6(*m6)

244 m0_unscaled = math.sqrt(num.sum(num.array(m9)**2)) / math.sqrt(2.)

245 m9 /= m0_unscaled

246 rel_m6 = moment_tensor.to6(m9)

247 for m, val in zip('rmnn rmee rmdd rmne rmnd rmed'.split(), rel_m6):

248 self.set_parameter(m, val)

249

250 else:

251 self.warn(

252 'No source mechanism available for event %s. '

253 'Only setting location' % event.name)

254

255 if event.duration is not None:

256 self.set_parameter('stf_duration', event.duration)

257 else:

258 self.warn(

259 'No source duration available for event %s. ' % event.name)

260

261 self.set_parameter('lat', event.lat)

262 self.set_parameter('lon', event.lon)

263 self.set_parameter('depth_km', event.depth/km)

264

265 def add_store(self):

266 self._engine = self.get_engine()

267 superdir = self.input_directory()

268 if self.has_config(superdir):

269 self._engine.store_dirs.append(superdir)

270 else:

271 self._engine.store_superdirs.append(superdir)

272 self.store_ids = self._engine.get_store_ids()

273

274 self.set_parameter_choices('store_id', self.store_ids)

275

276 def has_config(self, directory):

277 return 'config' in os.listdir(directory)

278

279

280class DCSource(Seismosizer):

281

282 def setup(self):

283 '''Customization of the snuffling.'''

284

285 self.set_name('Seismosizer: DCSource')

286 self.add_parameter(

287 Param('Time', 'time', 0.0, -50., 50.))

288 # self.add_parameter(

289 # Param('Latitude', 'lat', 0.0, -90., 90.))

290 # self.add_parameter(

291 # Param('Longitude', 'lon', 0.0, -180., 180.))

292 self.add_parameter(

293 Param('North shift', 'north_km', 0.0, -50., 50.))

294 self.add_parameter(

295 Param('East shift', 'east_km', 0.0, -50., 50.))

296 self.add_parameter(

297 Param('Depth', 'depth_km', 10.0, -100.0, 600.0))

298 self.add_parameter(

299 Param('Magnitude', 'magnitude', 6.0, 0.0, 10.0))

300 self.add_parameter(

301 Param('Strike', 'strike', 0., -180., 180.))

302 self.add_parameter(

303 Param('Dip', 'dip', 90., 0., 90.))

304 self.add_parameter(

305 Param('Rake', 'rake', 0., -180., 180.))

306 self.add_parameter(

307 Param('STF duration', 'stf_duration', 0., 0., 50.))

308 self.add_parameter(

309 Choice('STF type', 'stf_type', self.stf_types[0], self.stf_types))

310

311 Seismosizer.setup(self)

312

313 def get_source(self, event):

314 return gf.DCSource(

315 time=event.time+self.time,

316 lat=event.lat,

317 lon=event.lon,

318 north_shift=self.north_km*km,

319 east_shift=self.east_km*km,

320 depth=self.depth_km*km,

321 magnitude=self.magnitude,

322 strike=self.strike,

323 dip=self.dip,

324 rake=self.rake,

325 stf=self.get_stf())

326

327

328class MTSource(Seismosizer):

329

330 def setup(self):

331 '''Customization of the snuffling.'''

332

333 self.set_name('Seismosizer: MTSource')

334 self.add_parameter(

335 Param('Time', 'time', 0.0, -50., 50.))

336 # self.add_parameter(

337 # Param('Latitude', 'lat', 0.0, -90., 90.))

338 # self.add_parameter(

339 # Param('Longitude', 'lon', 0.0, -180., 180.))

340 self.add_parameter(

341 Param('North shift', 'north_km', 0.0, -50., 50.))

342 self.add_parameter(

343 Param('East shift', 'east_km', 0.0, -50., 50.))

344 self.add_parameter(

345 Param('Depth', 'depth_km', 10.0, -100.0, 600.0))

346 self.add_parameter(

347 Param('Magnitude', 'magnitude', 5.0, -1.0, 9.0))

348 self.add_parameter(

349 Param('rel Mnn', 'rmnn', .0, -1.41421, 1.41421))

350 self.add_parameter(

351 Param('rel Mee', 'rmee', .0, -1.41421, 1.41421))

352 self.add_parameter(

353 Param('rel Mdd', 'rmdd', .0, -1.41421, 1.41421))

354 self.add_parameter(

355 Param('rel Mne', 'rmne', 1., -1.0, 1.0))

356 self.add_parameter(

357 Param('rel Mnd', 'rmnd', .0, -1.0, 1.0))

358 self.add_parameter(

359 Param('rel Med', 'rmed', .0, -1.0, 1.0))

360 self.add_parameter(

361 Param('STF duration', 'stf_duration', 0., 0., 50.))

362 self.add_parameter(

363 Choice('STF type', 'stf_type', self.stf_types[0], self.stf_types))

364

365 self.add_trigger('Normalize M', self.normalize)

366

367 Seismosizer.setup(self)

368

369 def normalize(self):

370 for m, val in zip(

371 'rmnn rmee rmdd rmne rmnd rmed'.split(),

372 self.get_normalized_m6()):

373 self.set_parameter(m, val)

374

375 def get_normalized_m6(self):

376 rm6 = num.array(

377 [self.rmnn, self.rmee, self.rmdd, self.rmne, self.rmnd, self.rmed],

378 dtype=float)

379

380 m9 = moment_tensor.symmat6(*rm6)

381 m0_unscaled = math.sqrt(num.sum(num.array(m9)**2)) / math.sqrt(2.)

382 m9 /= m0_unscaled

383 return moment_tensor.to6(m9)

384

385 def get_m6(self):

386 m0 = moment_tensor.magnitude_to_moment(self.magnitude)

387 return self.get_normalized_m6() * m0

388

389 def get_source(self, event):

390 return gf.MTSource(

391 time=event.time+self.time,

392 lat=event.lat,

393 lon=event.lon,

394 north_shift=self.north_km*km,

395 east_shift=self.east_km*km,

396 depth=self.depth_km*km,

397 m6=self.get_m6(),

398 stf=self.get_stf())

399

400

401class SFSource(Seismosizer):

402 def setup(self):

403 '''Customization of the snuffling.'''

404

405 self.set_name('Seismosizer: SFSource')

406 self.add_parameter(

407 Param('Time', 'time', 0.0, -50., 50.))

408 # self.add_parameter(

409 # Param('Latitude', 'lat', 0.0, -90., 90.))

410 # self.add_parameter(

411 # Param('Longitude', 'lon', 0.0, -180., 180.))

412 self.add_parameter(

413 Param('North shift', 'north_km', 0.0, -50., 50.))

414 self.add_parameter(

415 Param('East shift', 'east_km', 0.0, -50., 50.))

416 self.add_parameter(

417 Param('Depth', 'depth_km', 10.0, -100.0, 600.0))

418 self.add_parameter(

419 Param('North force', 'fn', 1e3, -1e9, 1e9))

420 self.add_parameter(

421 Param('East force', 'fe', 1e3, -1e9, 1e9))

422 self.add_parameter(

423 Param('Down force', 'fd', 1e3, -1e9, 1e9))

424 self.add_parameter(

425 Param('STF duration', 'stf_duration', 0., 0., 100.))

426 self.add_parameter(

427 Choice('STF type', 'stf_type', self.stf_types[0], self.stf_types))

428

429 Seismosizer.setup(self)

430

431 def get_source(self, event):

432 return gf.SFSource(

433 time=event.time+self.time,

434 lat=event.lat,

435 lon=event.lon,

436 north_shift=self.north_km*km,

437 east_shift=self.east_km*km,

438 depth=self.depth_km*km,

439 fn=self.fn,

440 fe=self.fe,

441 fd=self.fd,

442 stf=self.get_stf())

443

444

445class RectangularSource(Seismosizer):

446

447 def setup(self):

448 '''Customization of the snuffling.'''

449

450 self.set_name('Seismosizer: RectangularSource')

451 self.add_parameter(

452 Param('Time', 'time', 0.0, -50., 50.))

453 # self.add_parameter(

454 # Param('Latitude', 'lat', 0.0, -90., 90.))

455 # self.add_parameter(

456 # Param('Longitude', 'lon', 0.0, -180., 180.))

457 self.add_parameter(

458 Param('North shift', 'north_km', 0.0, -50., 50.))

459 self.add_parameter(

460 Param('East shift', 'east_km', 0.0, -50., 50.))

461 self.add_parameter(

462 Param('Depth', 'depth_km', 10.0, 0.0, 600.0))

463 self.add_parameter(

464 Param('Magnitude', 'magnitude', 6.0, 0.0, 10.0))

465 self.add_parameter(

466 Param('Strike', 'strike', 0., -180., 180.))

467 self.add_parameter(

468 Param('Dip', 'dip', 90., 0., 90.))

469 self.add_parameter(

470 Param('Rake', 'rake', 0., -180., 180.))

471 self.add_parameter(

472 Param('Length', 'length', 10.*km, .1*km, 100*km))

473 self.add_parameter(

474 Param('Width', 'width', 5.*km, .1*km, 50*km))

475 self.add_parameter(

476 Param('Nucleation X', 'nucleation_x', 0., -1., 1.))

477 self.add_parameter(

478 Param('Nucleation Y', 'nucleation_y', 0., -1., 1.))

479 self.add_parameter(

480 Param('Rupture velocity', 'velocity', 3500.0, 0.0, 5000.0))

481 self.add_parameter(

482 Param('STF duration', 'stf_duration', 0., 0., 20.))

483 self.add_parameter(

484 Choice('STF type', 'stf_type', self.stf_types[0], self.stf_types))

485

486 Seismosizer.setup(self)

487

488 def get_source(self, event):

489 return gf.RectangularSource(

490 time=event.time+self.time,

491 lat=event.lat,

492 lon=event.lon,

493 north_shift=self.north_km*km,

494 east_shift=self.east_km*km,

495 depth=self.depth_km*km,

496 magnitude=self.magnitude,

497 strike=self.strike,

498 dip=self.dip,

499 rake=self.rake,

500 length=self.length,

501 width=self.width,

502 nucleation_x=self.nucleation_x,

503 nucleation_y=self.nucleation_y,

504 velocity=self.velocity,

505 stf=self.get_stf())

506

507

508class PseudoDynamicRuptureSource(Seismosizer):

509

510 def setup(self):

511 '''Customization of the snuffling.'''

512

513 self.set_name('Seismosizer: PseudoDynamicRupture')

514 self.add_parameter(

515 Param('Time', 'time', 0.0, -50., 50.))

516 # self.add_parameter(

517 # Param('Latitude', 'lat', 0.0, -90., 90.))

518 # self.add_parameter(

519 # Param('Longitude', 'lon', 0.0, -180., 180.))

520 self.add_parameter(

521 Param('North shift', 'north_km', 0.0, -50., 50.))

522 self.add_parameter(

523 Param('East shift', 'east_km', 0.0, -50., 50.))

524 self.add_parameter(

525 Param('Depth', 'depth_km', 10.0, 0.0, 600.0))

526 self.add_parameter(

527 Param('Slip', 'slip', 1.0, 0.0, 20.0))

528 self.add_parameter(

529 Param('Strike', 'strike', 0., -180., 180.))

530 self.add_parameter(

531 Param('Dip', 'dip', 90., 0., 90.))

532 self.add_parameter(

533 Param('Rake', 'rake', 0., -180., 180.))

534 self.add_parameter(

535 Param('Length', 'length', 10.*km, .1*km, 300*km))

536 self.add_parameter(

537 Param('Width', 'width', 5.*km, .1*km, 50*km))

538 self.add_parameter(

539 Param('Nucleation X', 'nucleation_x', 0., -1., 1.))

540 self.add_parameter(

541 Param('Nucleation Y', 'nucleation_y', 0., -1., 1.))

542 self.add_parameter(

543 Param('Gamma', 'gamma', 0.8, 0.1, 1.5))

544 self.add_parameter(

545 Param('nx', 'nx', 5, 1, 20))

546 self.add_parameter(

547 Param('ny', 'ny', 5, 1, 20))

548 self.add_parameter(

549 Param('STF duration', 'stf_duration', 0., 0., 20.))

550 self.add_parameter(

551 Choice('STF type', 'stf_type', self.stf_types[0], self.stf_types))

552

553 self.add_parameter(Switch(

554 'Tapered tractions', 'tapered', True))

555

556 Seismosizer.setup(self)

557

558 def get_source(self, event):

559 source = gf.PseudoDynamicRupture(

560 time=event.time + self.time,

561 lat=event.lat,

562 lon=event.lon,

563 nx=int(self.nx),

564 ny=int(self.ny),

565 north_shift=self.north_km*km,

566 east_shift=self.east_km*km,

567 depth=self.depth_km*km,

568 slip=self.slip,

569 strike=self.strike,

570 dip=self.dip,

571 rake=self.rake,

572 length=self.length,

573 width=self.width,

574 nucleation_x=self.nucleation_x,

575 nucleation_y=self.nucleation_y,

576 gamma=self.gamma,

577 nthreads=5,

578 pure_shear=True,

579 smooth_rupture=True)

580

581 return source

582

583

584def __snufflings__():

585 '''Returns a list of snufflings to be exported by this module.'''

586 return [

587 DCSource(),

588 MTSource(),

589 SFSource(),

590 RectangularSource(),

591 PseudoDynamicRuptureSource()]