obspy.signal - Signal processing routines for ObsPy

Capabilities include filtering, triggering, rotation, instrument correction and coordinate transformations.

copyright:

The ObsPy Development Team (devs@obspy.org)

license:

GNU Lesser General Public License, Version 3 (https://www.gnu.org/copyleft/lesser.html)

Filter

The filter module provides various filters, including different bandpass, lowpass, highpass, bandstop and FIR filter.

Warning

Before filtering you should make sure that data is demeaned/detrended, e.g. using detrend(). Otherwise there can be massive artifacts from filtering.

The following example shows how to highpass a seismogram at 1.0Hz. In the example only the first trace is processed to see the changes in comparison with the other traces in the plot.

Note

The filter takes the data explicitly as argument (i.e. an numpy.ndarray) and therefore the sampling_rate needs to be also specified. It returns the filtered data. For Stream and Trace objects simply use their respective filtering methods Stream.filter() and Trace.filter().

>>> from obspy import read
>>> import obspy.signal
>>> st = read()
>>> tr = st[0]
>>> tr.data = obspy.signal.filter.highpass(
...     tr.data, 1.0, corners=1, zerophase=True, df=tr.stats.sampling_rate)
>>> st.plot()  

Working with the convenience methods implemented on Stream/Trace works similar:

>>> tr.filter('highpass', freq=1.0, corners=1, zerophase=True)
... 
<...Trace object at 0x...>

(Source code, png)

../_images/obspy-signal-1.png

Instrument Correction

The response of the instrument can be removed by the invsim module. The following example shows how to remove the the instrument response of a STS2 and simulate an instrument with 2Hz corner frequency.

>>> from obspy import read
>>> st = read()
>>> st.plot() 

(Source code, png)

../_images/obspy-signal-2.png

Now we apply the instrument correction and simulation:

>>> from obspy.signal.invsim import simulate_seismometer, corn_freq_2_paz
>>> inst2hz = corn_freq_2_paz(2.0)
>>> sts2 = {'gain': 60077000.0,
...         'poles': [(-0.037004+0.037016j),
...                   (-0.037004-0.037016j),
...                   (-251.33+0j),
...                   (-131.04-467.29j),
...                   (-131.04+467.29j)],
...         'sensitivity': 2516778400.0,
...         'zeros': [0j, 0j]}
>>> for tr in st:
...     df = tr.stats.sampling_rate
...     tr.data = simulate_seismometer(tr.data, df, paz_remove=sts2,
...                                    paz_simulate=inst2hz,
...                                    water_level=60.0)
>>> st.plot()  

Again, there are convenience methods implemented on Stream/Trace:

>>> tr.simulate(paz_remove=sts2, paz_simulate=inst2hz, water_level=60.0)
... 
<...Trace object at 0x...>

(Source code, png)

../_images/obspy-signal-3.png

Trigger

The trigger module provides various triggering algorithms, including different STA/LTA routines, Z-Detector, AR picker and the P-picker by M. Bear. The implementation is based on [Withers1998] and [Baer1987].

The following example demonstrates a recursive STA/LTA triggering:

>>> from obspy import read
>>> from obspy.signal.trigger import recursive_sta_lta, plot_trigger
>>> st = read()
>>> tr = st.select(component="Z")[0]
>>> tr.filter("bandpass", freqmin=1, freqmax=20)  
<...Trace object at 0x...>
>>> sta = 0.5
>>> lta = 4
>>> cft = recursive_sta_lta(tr.data, int(sta * tr.stats.sampling_rate),
...                        int(lta * tr.stats.sampling_rate))
>>> thrOn = 4
>>> thrOff = 0.7
>>> plot_trigger(tr, cft, thrOn, thrOff) 

(Source code, png)

../_images/obspy-signal-4.png

There is also a convenience method implemented on Stream/Trace. It works on and overwrites the traces waveform data and is intended for batch processing rather than for interactive determination of triggering parameters. But it also means that the trace’s built-in methods can be used.

>>> tr.trigger("recstalta", sta=0.5, lta=4)  
<...Trace object at 0x...>
>>> tr.plot()  

(Source code, png)

../_images/obspy-signal-5.png

For more examples check out the trigger in the Tutorial. For network coincidence refer to obspy.signal.trigger.coincidence_trigger() and the same page in the Tutorial. For automated use see the following stalta example scripts.

Classes & Functions

array_processing

Method for Seismic-Array-Beamforming/FK-Analysis/Capon

array_rotation_strain

This routine calculates the best-fitting rigid body rotation and uniform strain as functions of time, and their formal errors, given three-component ground motion time series recorded on a seismic array.

ar_pick

Pick P and S arrivals with an AR-AIC + STA/LTA algorithm.

bandpass

Butterworth-Bandpass Filter.

bandstop

Butterworth-Bandstop Filter.

carl_sta_trig

Computes the carlSTAtrig characteristic function.

classic_sta_lta

Computes the standard STA/LTA from a given input array a.

coincidence_trigger

Perform a network coincidence trigger.

corn_freq_2_paz

Convert corner frequency and damping to poles and zeros.

cosine_taper

Cosine Taper.

delayed_sta_lta

Delayed STA/LTA.

envelope

Envelope of a function.

interpolate_1d

Wrapper around some scipy interpolation functions.

weighted_average_slopes

Implements the weighted average slopes interpolation scheme proposed in [Wiggins1976] for evenly sampled data.

estimate_magnitude

Estimate local magnitude.

evalresp

Get the instrument response from a SEED RESP-file.

highpass

Butterworth-Highpass Filter.

lowpass

Butterworth-Lowpass Filter.

paz_to_freq_resp

Convert Poles and Zeros (PAZ) to frequency response.

pk_baer

Wrapper for P-picker routine by M.

polarization_analysis

Method carrying out polarization analysis with the [Flinn1965b], [Jurkevics1988], ParticleMotion, or [Vidale1986] algorithm.

linear_regression

Use linear least squares to fit a function, f, to data.

PPSD

Class to compile probabilistic power spectral densities for one combination of network/station/location/channel/sampling_rate.

MSEEDMetadata

A container for MiniSEED specific metadata, including quality control parameters.

recursive_sta_lta

Recursive STA/LTA.

rotate_ne_rt

Rotates horizontal components of a seismogram.

simulate_seismometer

Simulate/Correct seismometer.

util_geo_km

Transform lon, lat to km with reference to orig_lon and orig_lat on the elliptic Earth.

util_lon_lat

Transform x, y [km] to decimal degree in reference to orig_lon and orig_lat

correlate

Cross-correlation of two signals up to a specified maximal shift.

z_detect

Z-detector.

Modules

array_analysis

Functions for Array Analysis

calibration

Functions for relative calibration.

cpxtrace

Complex Trace Analysis

cross_correlation

Signal processing routines based on cross correlation techniques.

detrend

Python module containing detrend methods.

differentiate_and_integrate

Integration and differentiation routines.

filter

Various Seismogram Filtering Functions

freqattributes

Frequency Attributes

hoctavbands

Half Octave Bands

interpolation

Some Seismogram Interpolating Functions.

invsim

Python Module for Instrument Correction (Seismology).

konnoohmachismoothing

Functions to smooth spectra with the so called Konno & Ohmachi method.

polarization

Functions for polarization analysis.

quality_control

Quality control module for ObsPy.

regression

Python Module for (Weighted) Linear Regression.

rotate

Various Seismogram Rotation Functions

spectral_estimation

Various Routines Related to Spectral Estimation

tf_misfit

Various Time Frequency Misfit Functions based on [Kristekova2006] and [Kristekova2009].

trigger

Various routines related to triggering/picking

util

Various additional utilities for obspy.signal.