Download seismic waveforms using Obspy¶

Reference: https://docs.obspy.org/packages/obspy.clients.fdsn.html

In [1]:

import obspy
import cartopy
from datetime import datetime
import numpy as np
from obspy import UTCDateTime
from obspy.clients.fdsn import Client
client = Client("IRIS")

import obspy import cartopy from datetime import datetime import numpy as np from obspy import UTCDateTime from obspy.clients.fdsn import Client client = Client("IRIS")

Explore Seismic Networks¶

As an example, we can find all CI stations of Southern California Seismic Network within a region of 33.48 N - 33.84 N, 117.19 - 118.94 W from IRIS

We then can find station information using Obspy:

In [2]:

minlatitude=33.48
maxlatitude=34.84
minlongitude=-118.94
maxlongitude=-117.19
stations = client.get_stations(network="CI", station="*", location="*", channel="HH*",
                               minlatitude=minlatitude, maxlatitude=maxlatitude,
                               minlongitude=minlongitude, maxlongitude=maxlongitude)
stations.plot(projection="local", resolution="f");

minlatitude=33.48 maxlatitude=34.84 minlongitude=-118.94 maxlongitude=-117.19 stations = client.get_stations(network="CI", station="*", location="*", channel="HH*", minlatitude=minlatitude, maxlatitude=maxlatitude, minlongitude=minlongitude, maxlongitude=maxlongitude) stations.plot(projection="local", resolution="f");

Explore Earthquakes¶

We can also use Obspy to find out the events occured in the past month:

In [3]:

now = UTCDateTime(datetime.now())
# starttime = now - 3600.0*24*30#s
# endtime = now
starttime = UTCDateTime("2022-05-26T00:00:00.000")
endtime = UTCDateTime("2022-06-26T00:00:00.000")
events = client.get_events(starttime=starttime, endtime=endtime,
                           minmagnitude=2.5,
                           minlatitude=minlatitude, maxlatitude=maxlatitude,
                           minlongitude=minlongitude, maxlongitude=maxlongitude)
events.plot(projection="local", resolution="f");

now = UTCDateTime(datetime.now()) # starttime = now - 3600.0*24*30#s # endtime = now starttime = UTCDateTime("2022-05-26T00:00:00.000") endtime = UTCDateTime("2022-06-26T00:00:00.000") events = client.get_events(starttime=starttime, endtime=endtime, minmagnitude=2.5, minlatitude=minlatitude, maxlatitude=maxlatitude, minlongitude=minlongitude, maxlongitude=maxlongitude) events.plot(projection="local", resolution="f");

In [4]:

print(events)

print(events)

1 Event(s) in Catalog:
2022-06-03T12:05:21.050000Z | +33.901, -118.414 | 2.51 Ml

Download Waveforms¶

Based on the stations and events information, we now can download the waveforms for analysis.

We can first select one event and one station to look at the waveform:

In [5]:

# We select the station FMP
station = stations.select(station="FMP")
print(station[0][0])

# We select the station FMP station = stations.select(station="FMP") print(station[0][0])

Station FMP (Fort Macarthur Park)
	Station Code: FMP
	Channel Count: 0/175 (Selected/Total)
	2000-08-30T00:00:00.000000Z - 
	Access: open 
	Latitude: 33.71, Longitude: -118.29, Elevation: 89.0 m
	Available Channels:

In [6]:

# The event information:
print(events[0].origins[0])

# The event information: print(events[0].origins[0])

Origin
	   resource_id: ResourceIdentifier(id="smi:service.iris.edu/fdsnws/event/1/query?originid=46955690")
	          time: UTCDateTime(2022, 6, 3, 12, 5, 21, 50000)
	     longitude: -118.413667
	      latitude: 33.9005
	         depth: 11270.0
	 creation_info: CreationInfo(author='ci,us')

In [7]:

# download and plot the waveforms
waveforms = client.get_waveforms(network=station[0].code,
                                 station=station[0][0].code,
                                 location="*",
                                 channel="HH*", 
                                 starttime=events[0].origins[0].time, 
                                 endtime=events[0].origins[0].time+60)  
waveforms.plot();

# download and plot the waveforms waveforms = client.get_waveforms(network=station[0].code, station=station[0][0].code, location="*", channel="HH*", starttime=events[0].origins[0].time, endtime=events[0].origins[0].time+60) waveforms.plot();

Save waveforms for subsequent processing:

Compare waveforms from seismic network and raspberry shake network¶

In [8]:

# download and plot the waveforms from the seismic network
import obspy
from obspy.clients.fdsn import Client
client = Client("NCEDC") ## northern california
# client = Client("SCEDC") ## sourthern california
# client = Client("IRIS") ## global
starttime = UTCDateTime("2022-06-29T14:08:04")
waveforms = client.get_waveforms(network="NC",
                                 station="KPR",
                                 location="*",
                                 channel="*", 
                                 starttime=starttime, 
                                 endtime=starttime+120) #s
waveforms[0].plot();

# download and plot the waveforms from the seismic network import obspy from obspy.clients.fdsn import Client client = Client("NCEDC") ## northern california # client = Client("SCEDC") ## sourthern california # client = Client("IRIS") ## global starttime = UTCDateTime("2022-06-29T14:08:04") waveforms = client.get_waveforms(network="NC", station="KPR", location="*", channel="*", starttime=starttime, endtime=starttime+120) #s waveforms[0].plot();

In [9]:

# download and plot the waveforms from the raspberry shake network
client = Client("RASPISHAKE")
waveforms = client.get_waveforms(network="AM",
                                 station="R5E62", ## change this!
                                 location="*",
                                 channel="*", 
                                 starttime=UTCDateTime("2022-06-29T14:08:04"), 
                                 endtime=UTCDateTime("2022-06-29T14:08:04")+120)  
waveforms[0].plot();

# download and plot the waveforms from the raspberry shake network client = Client("RASPISHAKE") waveforms = client.get_waveforms(network="AM", station="R5E62", ## change this! location="*", channel="*", starttime=UTCDateTime("2022-06-29T14:08:04"), endtime=UTCDateTime("2022-06-29T14:08:04")+120) waveforms[0].plot();

In [10]:

waveforms.write("waveforms.mseed")

waveforms.write("waveforms.mseed")

Download waveforms of multiple stations¶

In [11]:

from obspy.clients.fdsn import RoutingClient
client = RoutingClient("iris-federator")

from obspy.clients.fdsn import RoutingClient client = RoutingClient("iris-federator")

In [12]:

stations = client.get_stations(network="CI",
                               channel="HHZ", 
                               starttime=events[0].origins[0].time, 
                               endtime=events[0].origins[0].time+60,
                               latitude=events[0].origins[0].latitude,
                               longitude=events[0].origins[0].longitude,
                               maxradius=0.5, level="channel", )  
# stations.plot(projection="local", resolution="f");

stations = client.get_stations(network="CI", channel="HHZ", starttime=events[0].origins[0].time, endtime=events[0].origins[0].time+60, latitude=events[0].origins[0].latitude, longitude=events[0].origins[0].longitude, maxradius=0.5, level="channel", ) # stations.plot(projection="local", resolution="f");

In [13]:

waveforms = client.get_waveforms(network="CI",
                                 channel="HHZ", 
                                 starttime=events[0].origins[0].time, 
                                 endtime=events[0].origins[0].time+60,
                                 latitude=events[0].origins[0].latitude,
                                 longitude=events[0].origins[0].longitude,
                                 maxradius=0.5, level="channel", )  
# waveforms.plot();

waveforms = client.get_waveforms(network="CI", channel="HHZ", starttime=events[0].origins[0].time, endtime=events[0].origins[0].time+60, latitude=events[0].origins[0].latitude, longitude=events[0].origins[0].longitude, maxradius=0.5, level="channel", ) # waveforms.plot();

In [14]:

degree2meter = 111.32 * 1e3 #m
for waveform in waveforms:
    station = stations.select(station=waveform.stats.station)[0][0]
    waveform.stats.distance = np.sqrt((station.longitude-events[0].origins[0].longitude)**2 + (station.latitude-events[0].origins[0].latitude)**2) * degree2meter
waveforms.plot(type='section');

degree2meter = 111.32 * 1e3 #m for waveform in waveforms: station = stations.select(station=waveform.stats.station)[0][0] waveform.stats.distance = np.sqrt((station.longitude-events[0].origins[0].longitude)**2 + (station.latitude-events[0].origins[0].latitude)**2) * degree2meter waveforms.plot(type='section');

Download Waveforms from the Raspberry Shake network¶

We can also using the Raspberry Shake network to conduct a similar experiment:

Now we need to use the client of "RASPISHAKE" instead of "IRIS"

And the network code of Raspberry Shake is AM

You can find more details from their manual.

• Explore station of the Raspberry Shake network using stationview

• Or download station information using Obspy

In [15]:

client = Client("RASPISHAKE")
stations = client.get_stations(network="AM", station="*", location="*", channel="*",
                               minlatitude=minlatitude, maxlatitude=maxlatitude,
                               minlongitude=minlongitude, maxlongitude=maxlongitude)
stations.plot(projection="local", resolution="f");

client = Client("RASPISHAKE") stations = client.get_stations(network="AM", station="*", location="*", channel="*", minlatitude=minlatitude, maxlatitude=maxlatitude, minlongitude=minlongitude, maxlongitude=maxlongitude) stations.plot(projection="local", resolution="f");

• Downlaod waveform of the same event using the station R2285 that is close to the event:

In [16]:

client = Client("RASPISHAKE")
waveforms = client.get_waveforms(network="AM",
                                 station="R2285",
                                 location="*",
                                 channel="*", 
                                 starttime=events[0].origins[0].time, 
                                 endtime=events[0].origins[0].time+60)  
waveforms.plot();

client = Client("RASPISHAKE") waveforms = client.get_waveforms(network="AM", station="R2285", location="*", channel="*", starttime=events[0].origins[0].time, endtime=events[0].origins[0].time+60) waveforms.plot();

Compare waveforms from seismic network and raspberry shake network¶

In [17]:

# download and plot the waveforms from the seismic network
import obspy
from obspy.clients.fdsn import Client
client = Client("NCEDC") ## northern california
# client = Client("SCEDC") ## sourthern california
# client = Client("IRIS") ## global
starttime = UTCDateTime("2022-06-29T14:08:04")
waveforms = client.get_waveforms(network="NC",
                                 station="KPR",
                                 location="*",
                                 channel="*", 
                                 starttime=starttime, 
                                 endtime=starttime+120) #s
waveforms[0].plot();

# download and plot the waveforms from the seismic network import obspy from obspy.clients.fdsn import Client client = Client("NCEDC") ## northern california # client = Client("SCEDC") ## sourthern california # client = Client("IRIS") ## global starttime = UTCDateTime("2022-06-29T14:08:04") waveforms = client.get_waveforms(network="NC", station="KPR", location="*", channel="*", starttime=starttime, endtime=starttime+120) #s waveforms[0].plot();

In [18]:

# download and plot the waveforms from the raspberry shake network
client = Client("RASPISHAKE")
waveforms = client.get_waveforms(network="AM",
                                 station="R5E62", ## change this!
                                 location="*",
                                 channel="*", 
                                 starttime=UTCDateTime("2022-06-29T14:08:04"), 
                                 endtime=UTCDateTime("2022-06-29T14:08:04")+120)  
waveforms[0].plot();

# download and plot the waveforms from the raspberry shake network client = Client("RASPISHAKE") waveforms = client.get_waveforms(network="AM", station="R5E62", ## change this! location="*", channel="*", starttime=UTCDateTime("2022-06-29T14:08:04"), endtime=UTCDateTime("2022-06-29T14:08:04")+120) waveforms[0].plot();

Use multiple stations to detect earthquakes¶

For example, we can look at the waveforms of an earthquake occurred on 2022-07-08:

In [19]:

client = Client("RASPISHAKE")
waveforms = client.get_waveforms(network="AM",
                                 station="R8E9F,RF5B0,R42DE,R6D1A,R891A,RE924,R98B8,RAA32,R06D3,RB07E,RE569",
                                 location="*",
                                 channel="*", 
                                 starttime=UTCDateTime("2022-07-08 13:10:43"), 
                                 endtime=UTCDateTime("2022-07-08 13:10:43")+20)  
waveforms.plot();

client = Client("RASPISHAKE") waveforms = client.get_waveforms(network="AM", station="R8E9F,RF5B0,R42DE,R6D1A,R891A,RE924,R98B8,RAA32,R06D3,RB07E,RE569", location="*", channel="*", starttime=UTCDateTime("2022-07-08 13:10:43"), endtime=UTCDateTime("2022-07-08 13:10:43")+20) waveforms.plot();

In [22]:

waveforms = client.get_waveforms(network="AM",
                                 station="R8E9F,RF5B0,R42DE,R6D1A,R891A,RE924,R98B8,RAA32,R06D3,RB07E,RE569",
                                 location="*",
                                 channel="*", 
                                 starttime=UTCDateTime("2022-07-15 01:19:07"), 
                                 endtime=UTCDateTime("2022-07-15 01:19:07")+120)  
waveforms.normalize()
waveforms.plot();

waveforms = client.get_waveforms(network="AM", station="R8E9F,RF5B0,R42DE,R6D1A,R891A,RE924,R98B8,RAA32,R06D3,RB07E,RE569", location="*", channel="*", starttime=UTCDateTime("2022-07-15 01:19:07"), endtime=UTCDateTime("2022-07-15 01:19:07")+120) waveforms.normalize() waveforms.plot();

In [ ]: