Seismo-Live: http://seismo-live.org
%matplotlib inline
import matplotlib.pyplot as plt
plt.style.use("bmh")
plt.rcParams['figure.figsize'] = 10, 6
data/GR.FUR..BHN.D.2015.361
(station FUR
, LMU geophysical observatory in Fürstenfeldbruck)data/station_FUR.stationxml
from obspy import read, read_inventory
st = read("data/GR.FUR..BHN.D.2015.361")
inv = read_inventory("data/station_FUR.stationxml")
print(st)
print(inv)
inv.plot(projection="ortho");
1 Trace(s) in Stream: GR.FUR..BHN | 2015-12-27T00:00:09.769999Z - 2015-12-28T00:00:11.669999Z | 20.0 Hz, 1728039 samples Inventory created at 2016-05-03T10:14:19.000000Z Created by: JANE WEB SERVICE: fdsnws-station | Jane version: 0.0.0-gd6e1 http://jane/fdsnws/station/1/query?network=GR&station=FUR&starttime... Sending institution: Jane (Jane) Contains: Networks (1): GR Stations (1): GR.FUR (Fuerstenfeldbruck, Bavaria, GR-Net) Channels (2): GR.FUR..BHN, GR.FUR..HHN
/Users/lion/miniconda3/envs/seismo_live/lib/python3.7/site-packages/obspy/imaging/maps.py:45: UserWarning: basemap/pyproj with proj4 version >= 5 has a bug that results in inverted map axes. Your maps may be wrong. Please use another version of proj4, or use cartopy. warnings.warn(msg)
PPSD
class from obspy.signal, see http://docs.obspy.org/tutorial/code_snippets/probabilistic_power_spectral_density.html (but use the inventory you read from StationXML as metadata)PPSD
(plot()
method attached to PPSD
object)from obspy.signal import PPSD
tr = st[0]
ppsd = PPSD(stats=tr.stats, metadata=inv)
ppsd.add(tr)
ppsd.plot()
/Users/lion/miniconda3/envs/seismo_live/lib/python3.7/site-packages/obspy/signal/spectral_estimation.py:1830: MatplotlibDeprecationWarning: The set_clim function was deprecated in Matplotlib 3.1 and will be removed in 3.3. Use ScalarMappable.set_clim instead. cb.set_clim(*fig.ppsd.color_limits)
Since longer term stacks would need too much waveform data and take way too long to compute, we prepared one year continuous data preprocessed for a single channel of station FUR
to play with..
PPSD_FUR_HHN.npz
using PPSD
's load_npz()
staticmethod (i.e. it is called directly from the class, not an instance object of the class)max_percentage
option of plot()
option) (might take a couple of minutes..!)from obspy.signal import PPSD
ppsd = PPSD.load_npz("data/PPSD_FUR_HHN.npz")
ppsd.plot(max_percentage=10)
ppsd.plot(cumulative=True)
calculate_histogram()
(see docs!) method of PPSD
and visualize themcallback
option and use some crazy custom callback function in calculate_histogram()
, e.g. stack together all data from birthdays in your family.. or all German holidays + Sundays in the time span.. or from dates of some bands' concerts on a tour.. etc.ppsd.calculate_histogram(time_of_weekday=[(-1, 0, 2), (-1, 22, 24)])
ppsd.plot(max_percentage=10)
ppsd.calculate_histogram(time_of_weekday=[(-1, 8, 16)])
ppsd.plot(max_percentage=10)
calculate_histogram()
(see docs!) method of PPSD
and visualize themppsd.calculate_histogram(time_of_weekday=[(1, 0, 24), (2, 0, 24), (3, 0, 24), (4, 0, 24), (5, 0, 24)])
ppsd.plot(max_percentage=10)
ppsd.calculate_histogram(time_of_weekday=[(6, 0, 24), (7, 0, 24)])
ppsd.plot(max_percentage=10)
calculate_histogram()
(see docs!) method of PPSD
and visualize themppsd.calculate_histogram(month=[10, 11, 12, 1])
ppsd.plot(max_percentage=10)
ppsd.calculate_histogram(month=[4, 5, 6, 7])
ppsd.plot(max_percentage=10)
calculate_histogram()
(see docs!) method of PPSD
and visualize them