Seismo-Live: http://seismo-live.org
Authors:¶
- Martin van Driel (@martinvandriel)
- Lion Krischer (@krischer)
First Basic Exercise¶
Task: Calculate three component synthetic seismograms for the stations and events in the data/events and data/stations subdirectories and save them on disc.
Notes¶
- Receiver objects can also be created from StationXML, SEED, or STATIONS files as well as obpy inventories using
instaseis.Receiver.parse(); see the documentation for details. - Source objects can also be created from QuakeML, CMTSOLUTIONS, and in other ways using
instaseis.Source.parse(); see the documentation for details. The
get_seismograms()method has a couple of extra arguments:kind:displacement,velocity,accelerationremove_source_shift,reconvolve_stf,dt,
... see the documentation for details.
- You can use the properties of the Receiver and Source objects to create usefull filenames.
Basic lines to set up the notebook and some paths.
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import os
import obspy
plt.style.use('ggplot')
plt.rcParams['figure.figsize'] = (10, 8)
Import Instaseis and open the database:
import instaseis
db = instaseis.open_db("data/database")
1. Load Receivers¶
reminder: you can use ObsPy to load stations and plot a map:
from obspy import read_inventory
inventory = read_inventory('data/stations/all_stations.xml')
inventory.plot(projection="local", resolution="i");
This inventory can directly be used as input to generate a list of instaseis.Receiver objects:
receivers = instaseis.Receiver.parse(inventory)
for rec in receivers[:2]:
print(rec)
Alternatively, instaseis can directly open the station xml or STATIONS file (but then you don't have the nice plot):
receivers = instaseis.Receiver.parse('data/stations/all_stations.xml')
print(receivers[0])
receivers = instaseis.Receiver.parse('data/stations/STATIONS')
print(receivers[0])
2. Load Events¶
reminder: use ObsPy to load events from a QuakeML file containing all events and plot a map:
import glob # provides iterator to loop over files
cat = obspy.core.event.Catalog()
for filename in glob.iglob('data/events/quakeml/*.xml'):
cat += obspy.read_events(filename)
print(cat)
print(instaseis.Source.parse(cat.events[0]))
cat.plot();
Alternatively load QuakeML or CMTSOLUTION files directly using instaseis.Source.parse() and store the sources in a list:
sources = []
for filename in glob.iglob('data/events/quakeml/*.xml'):
sources.append(instaseis.Source.parse(filename))
print(sources[0])
for filename in glob.iglob('data/events/cmtsolutions/*'):
sources.append(instaseis.Source.parse(filename))
print(sources[0])
3. Extract Seismograms and Save to File¶
For the first solution using a ObsPy event catalog:
Solutions:¶
dt = 1.0
for event in cat:
src = instaseis.Source.parse(event)
srcname = '%s_Mw_%3.1f' % (src.origin_time.date, src.moment_magnitude)
for rec in receivers:
# create a usefull filename
recname = '%s_%s' % (rec.network, rec.station)
filename = '%s_%s' % (recname, srcname)
filename = filename.replace('.', '_')
# extract seismograms using instaseis
st = db.get_seismograms(source=src, receiver=rec, dt=dt)
# write to miniseed files in the data_out folder. Write as MiniSEED due to multi
# component support.
st.write(os.path.join('data_out', filename + '.mseed'), format='mseed')
For the second solution use a list of sources:
Solutions:¶
dt = 1.0
for src in sources:
srcname = '%s_Mw_%3.1f' % (src.origin_time.date, src.moment_magnitude)
for rec in receivers:
# create a usefull filename
recname = '%s_%s' % (rec.network, rec.station)
filename = '%s_%s' % (recname, srcname)
filename = filename.replace('.', '_')
# extract seismograms using instaseis
st = db.get_seismograms(source=src, receiver=rec, dt=dt)
# write to miniseed files in the data_out folder. Write as MiniSEED due to multi
# component support.
st.write(os.path.join('data_out', filename + '.mseed'), format='mseed')
Check the data:
ls data_out/