Waveform plot

This will download continuous seismic waveforms & plot them and requires ObsPy

Import ObsPy module

In [1]:
from obspy import read
from obspy import UTCDateTime
from obspy.clients.fdsn import Client

import obspy as ob
print("# obspy version = ",ob.__version__)
# obspy version =  1.2.2

Set client (Data Center)

This example uses NCEDC. We can use other dataceneter (e.g., SCEDC, IRIS)

In [2]:
client = Client("NCEDC") # data from NCEDC 
#client = Client("SCEDC") # data from SCEDC
#client = Client("IRIS") # data from IRIS

Set SNCL

Which SNCL (Station, Network, Component, Location)? This example uses BRK.BK.LHZ.00 data

In [3]:
# BRK BHZ data (https://seismo.berkeley.edu/station_book/brk.html; Haviland Hall)
sta = "BRK" # station
#com = "HHZ" # component (HHZ 100sps; BHZ 40sps, LHZ 1sps)
com = "LHZ"
net = "BK" # network
loc = "00" # location "--" for blank location code

Set time window

This example uses 2-hour data for the 2020 M6.8 Philippines earthquake

In [4]:
# M 6.4 - 10 km WSW of Polloc, Philippines
# https://earthquake.usgs.gov/earthquakes/eventpage/us6000b80p/executive
#M 6.4 - 10 km WSW of Polloc, Philippines
#2020-08-01 17:09:01 (UTC)7.304°N 124.142°E479.6 km depth
start_day = "2020-08-01T17:09:01"
end_day = "2020-08-01T19:08:01"
starttime = UTCDateTime(start_day)
endtime = UTCDateTime(end_day)

Download seismic data

use get_waveforms to download data and do st.plot() for plotting

In [5]:
st = client.get_waveforms(network=net, station=sta, location=loc, channel=com,
                     starttime=starttime, endtime=endtime, 
                     attach_response=True)
_plot = st.plot()

Correct instrument response

use remove_response to correct the instrument response. We can select output unit (displacement, velocity or accerelation)

In [6]:
st.detrend() # remove liner trend
st.taper(max_percentage=0.01) # apply taper
Out[6]:
1 Trace(s) in Stream:
BK.BRK.00.LHZ | 2020-08-01T17:09:01.069538Z - 2020-08-01T19:08:00.069538Z | 1.0 Hz, 7140 samples
In [7]:
st = st.remove_response( output="VEL" ) # get velocity data (m/s)
#st = st.remove_response( output="DISP" ) # get displacement data (m)
#st = st.remove_response( output="ACC" ) # get acceleration data (m/s^2)
In [8]:
_plot = st.plot()

Filtering seismic data

first remove liner trend, apply a cosin taper, and then do filtering

In [9]:
st.detrend() # remove liner trend
st.taper(max_percentage=0.01) # apply taper

fl = 0.05 # in Hz 
fh = 0.10 # in Hz
st.filter(type='bandpass', freqmin=fl, freqmax=fh, corners=6, zerophase=False)
Out[9]:
1 Trace(s) in Stream:
BK.BRK.00.LHZ | 2020-08-01T17:09:01.069538Z - 2020-08-01T19:08:00.069538Z | 1.0 Hz, 7140 samples
In [10]:
_plot = st.plot()
In [ ]:
 
In [ ]:
 
In [ ]:
 
In [ ]:
 
In [ ]: