Skip to content

Instantly share code, notes, and snippets.

@jamm1985
Last active Jul 25, 2020
Embed
What would you like to do?
Simple waveforms (from various earthquakes) plot with trace differentiation, normalization.
"""
File: wave_plots.py
Author: Andrey Stepnov
Email: a.stepnov@geophystech.ru, myjamm@gmail.com
Github: https://github.com/jamm1985
Description: puts waveforms from different earthquakes in single plot using obspy and matplotlib
"""
from obspy import read
from obspy import UTCDateTime
import numpy as np
import matplotlib
import matplotlib.pyplot as plt
tracelength = 15
textsize = 12
# https://eqalert.ru/#/events/Wgvx4NMj
# LNSK SHZ - velocity
# read waveform file
ch_1 = read("../wave_data/2013-10-14-0836-05S.IMGG__021")
ch_1 = ch_1.select(station='LNSK', channel='SHZ')
# slice time
dt = UTCDateTime("2013-10-14T08:36:25.60")
ch_1 = ch_1.slice(dt, dt + tracelength)
ch_1.differentiate()
ch_1.normalize()
#ch_1.plot()
# https://eqalert.ru/#/events/GMmGPNgm
# ARGI SHZ - velocity
ch_2 = read("../wave_data/2016-06-21-1759-25M.IMGG__018")
ch_2 = ch_2.select(station='ARGI',channel='SHZ')
dt = UTCDateTime("2016-06-21T17:59:55.87")
ch_2 = ch_2.slice(dt, dt + tracelength)
ch_2.differentiate()
ch_2.normalize()
#ch_2.plot()
# https://eqalert.ru/#/events/WgvnwNgj
# STO11 Z - acceleration
ch_3 = read("../wave_data/2017-06-08-1907-39S.ST011_003")
ch_3 = ch_3.select(station='ST011',channel='C03')
dt = UTCDateTime("2017-06-08T19:08:07.80")
ch_3 = ch_3.slice(dt, dt + tracelength)
ch_3.normalize()
# custom fix zero shift :)
ch_3[0].data = ch_3[0].data - 0.08
#ch_3.plot()
# plot
matplotlib.rcParams.update({'font.size': 15})
ax3 = plt.subplot(313)
plt.plot(ch_3[0].times(), ch_3[0].data, "g-")
plt.yticks(np.arange(-1, 2, step=1))
ax3.set_title(str(ch_3[0].stats.starttime),
fontsize = textsize,
loc = 'right',
pad = 1)
plt.xlabel('Seconds')
ax1 = plt.subplot(311, sharex=ax3)
plt.plot(ch_1[0].times(), ch_1[0].data, "b-")
plt.setp(ax1.get_xticklabels(), visible=False)
plt.yticks(np.arange(-1, 2, step=1))
ax1.set_title(str(ch_1[0].stats.starttime),
fontsize = textsize,
loc = 'right',
pad = 1)
ax2 = plt.subplot(312, sharex=ax3)
plt.plot(ch_2[0].times(), ch_2[0].data, "r-")
plt.setp(ax2.get_xticklabels(), visible=False)
plt.yticks(np.arange(-1, 2, step=1))
ax2.set_title(str(ch_2[0].stats.starttime),
fontsize = textsize,
loc = 'right',
pad = 1)
plt.show()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment