Classes

rfpy defines the following base classes:

RFData
HkStack
Harmonics
CCPimage

RFData

class rfpy.rfdata.RFData(sta=None)

A RFData object contains Class attributes that associate station information with a single event (i.e., earthquake) metadata, corresponding raw and rotated seismograms, spectra and receiver functions.

Note

The object is initialized with the sta field only, and other attributes are added to the object as the analysis proceeds.

Parameters:

sta (object) – Object containing station information - from stdb database. When RFData is initialized without sta, functions requiring station info coordinate will not be working.
meta (Meta) – Object of metadata information for single event (initially set to None)
data (Stream) – Stream object containing the three-component seismograms (either un-rotated or rotated by the method rotate())
rf (Stream) – Stream object containing the receiver functions (created when executing deconvlove())
specs (Stream) – Stream object containing the spectra of the data (created when executing calc_spectra())

add_event(event, gacmin=30.0, gacmax=90.0, phase='P', returned=False)

Adds event metadata to RFData object, including travel time info of P wave.

Parameters:

event (event) – Event XML object
returned (bool) – Whether or not to return the accept attribute

Returns:

accept – Whether or not the object is accepted for further analysis

Return type:

bool

add_data(stream, returned=False, new_sr=5.0)

Adds stream as object attribute

Parameters:

stream (Stream) – Stream container for NEZ seismograms
returned (bool) – Whether or not to return the accept attribute

zne_data

Stream container for NEZ seismograms

Type:: Stream

Returns:: accept – Whether or not the object is accepted for further analysis
Return type:: bool

download_data(client, stdata=[], dtype='SAC', ndval=nan, new_sr=5.0, dts=120.0, remove_response=False, local_response_dir='', returned=False, verbose=False)

Downloads seismograms based on event origin time and P phase arrival.

Parameters:

client (Client) – Client object
ndval (float) – Fill in value for missing data
new_sr (float) – New sampling rate (Hz)
dts (float) – Time duration (sec)
stdata (List) – Station list
remove_response (bool) – Remove instrument response from seismogram and resitute to true ground velocity (m/s) using obspy.core.trace.Trace.remove_response()
returned (bool) – Whether or not to return the accept attribute
verbose (bool) – Output diagnostics to screen

Returns:

accept – Whether or not the object is accepted for further analysis

Return type:

bool

data

Stream containing Trace objects

Type:: Stream

rotate(vp=None, vs=None, align=None)

Rotates 3-component seismograms from vertical (Z), east (E) and north (N) to longitudinal (L), radial (Q) and tangential (T) components of motion. Note that the method ‘rotate’ from obspy.core.stream.Stream is used for the rotation 'ZNE->ZRT' and 'ZNE->LQT'. Rotation 'ZNE->PVH' is implemented separately here due to different conventions.

Parameters:

vp (float) – P-wave velocity at surface (km/s)
vs (float) – S-wave velocity at surface (km/s)
align (str) – Alignment of coordinate system for rotation (‘ZRT’, ‘LQT’, or ‘PVH’)

Returns:

rotated – Whether or not the object has been rotated

Return type:

bool

Note

‘PVH’ following Equation 18 in: Kennett (1991). The removal of free surface interactions from three-component seismograms. Geophysical Journal International

calc_snr(dt=30.0, fmin=0.05, fmax=1.0)

Calculates signal-to-noise ratio on either Z, L or P component

Parameters:

dt (float) – Duration (sec)
fmin (float) – Minimum frequency corner for SNR filter (Hz)
fmax (float) – Maximum frequency corner for SNR filter (Hz)

snr

Signal-to-noise ratio on vertical component (dB)

Type:: float

snrh

Signal-to-noise ratio on radial component (dB)

Type:: float

calc_spectra(phase='P', vp=None, vs=None, align=None, method='wiener', wavelet='complete', envelope_threshold=0.05, time=5.0, pre_filt=None, writeto=False, norm=False, length=None)

Calculate the numerators and denominator of the spectral division for receiver function calculation using one component as the source wavelet. The source component is always taken as the dominant compressional component, which can be either ‘Z’, ‘L’, or ‘P’.

Parameters:

vp (float) – P-wave velocity at surface (km/s)
vs (float) – S-wave velocity at surface (km/s)
align (str) – Alignment of coordinate system for rotation (‘ZRT’, ‘LQT’, or ‘PVH’)
method (str) – Method for deconvolution. Options are ‘wiener’, ‘water’, ‘water2’ or ‘multitaper’
wavelet (str) – Type of wavelet for deconvolution. Options are ‘complete’, ‘time’ or ‘envelope’
envelope_threshold (float) – Threshold [0-1] used in wavelet='envelope'.
time (float) – Window length used in wavelet='time'. Minimum window length for wavelet='envelope'.
pre_filt (list of 2 floats) – Low and High frequency corners of bandpass filter applied before deconvolution
gfilt (float) – Center frequency of Gaussian filter (Hz).
wlevel (float) – Water level used in method='water' and method='water2'
writeto (str or None) – Write wavelets for deconvolution to file.
length (float or None) – Time length (s) of the numerator. If None, use entire traces.

specs

Stream containing the cross spectral quantities

Type:: Stream

deconvolve(align=None, method='wiener', gfilt=None, wlevel=0.01)

Deconvolves three-component data using one component as the source wavelet. The source component is always taken as the dominant compressional component, which can be either ‘Z’, ‘L’, or ‘P’.

Parameters:

align (str) – Alignment of coordinate system for rotation (‘ZRT’, ‘LQT’, or ‘PVH’)
method (str) – Method for deconvolution. Options are ‘wiener’, ‘water’ (aka ‘wlevel’), ‘water2’ (aka ‘wlevel2’), or ‘multitaper’
gfilt (float) – Center frequency of Gaussian filter (Hz).
wlevel (float) – Water level used in method='water'.

rf

Stream containing the receiver function traces

Type:: Stream

to_stream(store=None): Method to switch from RFData object to Stream object. This allows easier manipulation of the receiver functions for post-processing.

save(file)

Saves RFData object to file

Parameters:: file (str) – File name for RFData object

Class Meta

class rfpy.rfdata.Meta(sta, event, gacmin=30.0, gacmax=90.0, phase='P')

A Meta object contains attributes associated with the metadata for a single receiver function analysis.

Parameters:

time (UTCDateTime) – Origin time of earthquake
dep (float) – Depth of hypocenter (km)
lon (float) – Longitude coordinate of epicenter
lat (float) – Latitude coordinate of epicenter
mag (float) – Magnitude of earthquake
gac (float) – Great arc circle between station and epicenter (degrees)
epi_dist (float) – Epicentral distance between station and epicenter (km)
baz (float) – Back-azimuth - pointing to earthquake from station (degrees)
az (float) – Azimuth - pointing to station from earthquake (degrees)
ttime (float) – Predicted arrival time (sec)
ph (str) – Phase name
slow (float) – Horizontal slowness of phase
inc (float) – Incidence angle of phase at surface
vp (float) – P-wave velocity at surface (km/s)
vs (float) – S-wave velocity at surface (km/s)
align (str) – Alignment of coordinate system for rotation (‘ZRT’, ‘LQT’, or ‘PVH’)
rotated (bool) – Whether or not data have been rotated to align coordinate system

HkStack

class rfpy.hk.HkStack(rfV1, rfV2=None, strike=None, dip=None, vp=6.0)

A HkStack object contains attributes and methods to stack radial receiver functions along moveout curves for point measurements of the depth to Moho (H) and P-to-S velocity ratio (k) beneath a seismic stations. The object is initialized with at least one Stream containing observed (or synthetised) radial receiver functions. The methods available can produce linear weighted stacks or product stacks, and can be used in the presence of flat or dipping Moho (with known strike and dip).

Note

The object is initialized with the rfV1 field only, and other attributes are added to the object as the analysis proceeds. A second rfV2 can be included, which is typically a copy of rfV1 filtered at different corner frequencies and is used to stack along the Pps and Pss moveout curves.

Parameters:

rfV1 (Stream) – Stream object containing the radial-component receiver function seismograms
rfV2 (Stream) – Stream object containing the radial-component receiver function seismograms (typically filtered at lower frequencies)
strike (float) – Strike angle of dipping Moho (has to be known or estimated a priori)
dip (float) – Dip angle of Moho (has to be known or estimated a priori)
vp (float) – Mean P-wave velocity of the crust (km/s)
kbound (list) – List of 2 floats that determine the range of Vp/Vs values to search
dk (float) – Spacing between adjacent Vp/Vs search values
hbound (list) – List of 2 floats that determine the range of Moho depth values to search
dh (float) – Spacing between adjacent Moho depth search values
weights (list) – List of 3 floats that determine the relative weights of the individual phase stacks to be used in the final stack. The third weight is negative since Pss amplitudes are opposite to those of the Ps and Pps phases.
phases (list) – List of 3 strings (‘ps’, ‘pps’, ‘pss’) corresponding to the thre phases of interest (do not modify this attribute)

stack(vp=None)

Method to calculate Hk stacks from radial receiver functions. The stacks are calculated using phase-weighted stacking for individual phases and take the median of the weighted stack to avoid bias by outliers.

Note

If two streams are available as attributes, the method will assume that the second stream should be used for stacking along the Pps and Pss move out curves (e.g., if the second stream contains lower frequency signals). Furthermore, If the vp argument is not specified, the method will use the value set during initialization (vp=6.0 km/s)

Parameters:: vp (float) – Mean crust P-wave velocity (km/s).

pws

Array of phase stacks, where the outer dimension corresponds to the phase index (shape nH, nk, nph)

Type:: ndarray

sig

Variance of phase stacks, where the outer dimension corresponds to the phase index (shape nH, nk, nph)

Type:: ndarray

stack_dip(vp=None, strike=None, dip=None)

Method to calculate Hk stacks from radial receiver functions using known stike and dip angles of the Moho. The stacks are calculated using phase-weighted stacking for individual phases and take the median of the weighted stack to avoid bias by outliers.

Note

If two streams are available as attributes, the method will assume that the second stream should be used for stacking along the Pps and Pss move out curves (e.g., if the second stream contains lower frequency signals). Furthermore, If the arguments are not specified, the method will use the values set during initialization (vp=6.0 km/s, strike=0., dip=0.)

Parameters:

vp (float) – Mean crust P-wave velocity (km/s).
strike (float) – Strike angle of dipping Moho (has to be known or estimated a priori)
dip (float) – Dip angle of Moho (has to be known or estimated a priori)

pws

Array of phase stacks, where the outer dimension corresponds to the phase index (shape nH, nk, nph)

Type:: ndarray

sig

Variance of phase stacks, where the outer dimension corresponds to the phase index (shape nH, nk, nph)

Type:: ndarray

average(typ='sum', q=0.05, err_method='amp')

Method to combine the phase-weighted stacks to produce a final stack, from which to estimate the H and k parameters and their associated errors.

Parameters:

typ (str) – How the phase-weigthed stacks should be combined to produce a final stack. Available options are: weighted sum (typ=sum) or product (typ=product).
q (float) – Confidence level for the error estimate
err_method (str) – How errors should be estimated. Options are err_method='amp' to estimate errors from amplitude, or err_method='stats' to use a statistical F test from the residuals.

error(q=0.05, err_method='amp')

Method to determine the error on H and k estimates.

From Walsh, JGR, 2013

Parameters:

q (float) – Confidence level for the error estimate
err_method (str) – How errors should be estimated. Options are err_method='amp' to estimate errors from amplitude, or err_method='stats' to use a statistical F test from the residuals.

plot(save=False, title=None, form='png')

Method to plot H-K stacks. By default all 4 panels are plotted: The ps, pps and pss stacks, and the final (averaged) stack. Error contours are also plotted, as well as the position of the maximum stack values.

Parameters:

save (bool) – Whether or not to save the Figure
title (str) – Title of plot

save(file)

Saves HkStack object to file

Parameters:: file (str) – File name for HkStack object

Harmonics

class rfpy.harmonics.Harmonics(radialRF, transvRF=None, azim=0, xmin=0.0, xmax=10.0)

A Harmonics object contains attributes and methods to decompose radial and transverse component receiver functions into back-azimuth harmonics. The object is initialized with two Stream objects containing observed (or synthetised) radial and transverse receiver functions. The methods available can decompose the receiver functions along a fixed azimuth, or search for the optimal azimuth within a time range by minimizing one component.

Note

The object is initialized with the rfV1 field only, and other attributes are added to the object as the analysis proceeds. A second rfV2 can be included, which is typically a copy of rfV1 filtered at different corner frequencies and is used to stack along the Pps and Pss moveout curves.

Parameters:

radialRF (Stream) – Stream object containing the radial-component receiver function seismograms
transvRF (Stream) – Stream object containing the transverse-component receiver function seismograms
azim (float) – Direction (azimuth) along which the B1 component of the stream is minimized (between xmin and xmax)
xmin (float) – Minimum x axis value over which to calculate azim
xmax (float) – Maximum x axis value over which to calculate azim
hstream (Stream) – Stream containing the 5 harmonics, oriented in direction azim
radial_forward (Stream, optional) – Stream containing the radial receiver functions
transv_forward (Stream, optional) – Stream containing the transverse receiver functions

dcomp_find_azim(xmin=None, xmax=None)

Method to decompose radial and transverse receiver function streams into back-azimuth harmonics and determine the main orientation azim, obtained by minimizing the B1 component between xmin and xmax (i.e., time or depth).

Parameters:

xmin (float) – Minimum x axis value over which to calculate azim
xmax (float) – Maximum x axis value over which to calculate azim

hstream

Stream containing the 5 harmonics, oriented in direction azim

Type:: Stream

azim

Direction (azimuth) along which the B1 component of the stream is minimized (between xmin and xmax)

Type:: float

var

Variance of the 5 harmonics between xmin and xmax

Type:: ndarray

dcomp_fix_azim(azim=None)

Method to decompose radial and transverse receiver function streams into back-azimuth harmonics along direction azim.

Parameters:: azim (float) – Direction (azimuth) along which the B1 component of the stream is minimized (between xmin and xmax)

hstream

Stream containing the 5 harmonics, oriented in direction azim

Type:: Stream

forward(baz_list=None)

Method to forward calculate radial and transverse component receiver functions given the 5 pre-determined harmonics and a list of back-azimuth values. The receiver function signal parameters (length, sampling rate, etc.) will be identical to those in the stream of harmonic components.

Parameters:: baz_list (list) – List of back-azimuth directions over which to calculate the receiver functions. If no list is specified, the method will use the same back-azimuths as those in the original receiver function streams

radial_forward

Stream containing the radial receiver functions

Type:: Stream

transv_forward

Stream containing the transverse receiver functions

Type:: Stream

plot(ymax=30.0, scale=10.0, save=False, title=None, form='png')

Method to plot the 5 harmonic components.

Parameters:

ymax (float) – Maximum y axis value (time or depth) over which to plot the harmonic components
scale (float) – Scaling factor for the amplitudes (typically > 1)
save (bool) – Whether or not to save the plot
title (str) – Title of plot, to be used in the Figure and the file name (if save==True)

save(file)

Saves harmonics object to file

Parameters:: file (str) – File name for Harmonics object

CCPimage

class rfpy.ccp.CCPimage(coord_start=[None, None], coord_end=[None, None], weights=[1.0, 3.0, -3.0], dep=array([0., 4., 8., 14., 30., 35., 45., 110.]), vp=array([4., 5.9, 6.2, 6.3, 6.8, 7.2, 8., 8.1]), vs=None, vpvs=1.73, dx=2.5, dz=1.0)

A CCPimage object contains attributes and methods to produce Common Conversion Point (CCP) stacks for each of the main three main phases (Ps, Pps and Pss) using radial-component receiver functions. The object is used to project the stacks along a linear profile, specified by start and end geographical coordinate locations, which are subsequently averaged to produce a final CCP image. The averaging can be done using a linear weighted sum, or a phase-weighted sum. Methods should be used in the appropriate sequence (see rfpy_ccp.py for details).

Note

By default, the object initializes with un-defined coordinate locations for the profile. If not specified during initialization, make sure they are specified later, before the other methods are used, e.g. ccpimage.xs_lat1 = 10.; ccpimage.xs_lon1 = 110., etc. Note also that the default 1D velocity model may not be applicable to your region of interest and a different model can be implemented during initialization or later during processing.

Parameters:

coord_start (list) – List of two floats corresponding to the (latitude, longitude) pair for the start point of the profile
coord_end (list) – List of two floats corresponding to the (latitude, longitude) pair for the end point of the profile
weights (list) – List of three floats with corresponding weights for the Ps, Pps and Pss phases used during linear, weighted averaging
dep (ndarray) – Array of depth values defining the 1D background seismic velocity model. Note that the maximum depth defined here sets the maximum depth in each of the CCP stacks and the final CCP image.
vp (ndarray) – Array of Vp values defining the 1D background seismic velocity model
vpvs (float) – Constant Vp/Vs ratio for the 1D model.
radialRF (list) – List of Stream objects containing the radial receiver functions along the line. Each item in the list contains the streams for one single station.
vs (ndarray) – Array of Vp values defining the 1D background seismic velocity model
xs_lat1 (float) – Latitude of start point defining the linear profile.
xs_lon1 (float) – Longitude of start point defining the linear profile.
xs_lat2 (float) – Latitude of end point defining the linear profile.
xs_lon2 (float) – Longitude of end point defining the linear profile.
is_ready_for_prep (boolean) – Whether or not the object is ready for the method prep_data
is_ready_for_prestack (boolean) – Whether or not the object is ready for the method prestack
is_ready_for_ccp (boolean) – Whether or not the object is ready for the method ccp
is_ready_for_gccp (boolean) – Whether or not the object is ready for the method gccp

add_rfstream(rfstream)

Method to add a Stream object to the list radialRF. With at least one stream in radialRF, the object is ready for the prep_data method, and the corresponding flag is updated.

Parameters:: rfstream (Stream) – Stream object containing radial receiver functions for one station

prep_data(f1=0.05, f2ps=0.5, f2pps=0.25, f2pss=0.2, nbaz=37, nslow=41)

Method to pre-process the data and calculate the CCP points for each of the receiver functions. Pre-processing includes the binning to back-azimuth and slowness bins to reduce processing time and generate cleaner stacks, as well as filtering to emphasize the energy of the various phases as a function of depth. As a general rule, the high frequency corner of the 2 reverberated phases (Pps and Pss) should be approximately half of the high frequency corner of the direct (Ps) phase. At the end of this step, the object is updated with the amplitude at the lat, lon and depth values corresponding to the raypath for each receiver function. The object is now ready for the method prestack, with the corresponding flag updated.

Parameters:

f1 (float) – Low-frequency corner of the bandpass filter used for all phases (Hz)
f2ps (float) – High-frequency corner of the bandpass filter used for the Ps phase (Hz)
f2pps (float) – High-frequency corner of the bandpass filter used for the Pps phase (Hz)
f2pss (float) – High-frequency corner of the bandpass filter used for the Pss phase (Hz)
nbaz (int) – Number of increments in the back-azimuth bins
nslow (int) – Number of increments in the slowness bins
object (The following attributes are added to the) –
amp_ps_depth (numpy.ndarray) – 2D array of amplitudes as a function of depth for the Ps phase
amp_pps_depth (numpy.ndarray) – 2D array of amplitudes as a function of depth for the Pps phase
amp_pss_depth (numpy.ndarray) – 2D array of amplitudes as a function of depth for the Pss phase
lon_depth (numpy.ndarray) – 2D array of longitude as a function of depth (i.e., piercing points)
lat_depth (numpy.ndarray) – 2D array of latitude as a function of depth (i.e., piercing points)
is_ready_for_presstack (boolean) – Flag specifying that the object is ready for the presstack() method
n_traces (float) – Total number of traces used in creating the CCP image.

prestack()

Method to project the raypaths onto the 2D profile for each of the three phases. The final grid is defined here, using the parameter dx in km. The horizontal extent is pre-determined from the start and end points of the profile. At the end of this step, the object contains the set of amplitudes at each of the 2D grid points, for each of the three phases. The object is now ready for the methods ccp and/or gccp, with the corresponding flag updated.

The following attributes are added to the object:

Parameters:

xs_amps_ps (numpy.ndarray) – 3D array of amplitudes (1D array of amplitudes at each grid cell) for the Ps phase
xs_amps_pps (numpy.ndarray) – 3D array of amplitudes (1D array of amplitudes at each grid cell) for the Pps phase
xs_amps_pss (numpy.ndarray) – 3D array of amplitudes (1D array of amplitudes at each grid cell) for the Pss phase
is_ready_for_ccp (boolean) – Flag specifying that the object is ready for the ccp() method
is_ready_for_gccp (boolean) – Flag specifying that the object is ready for the gccp() method

ccp()

Method to average the amplitudes at each grid point to produce 2D images for each of the three phases. At the end of this step, the object contains the three 2D arrays that can be further averaged into a single final image.

The following attributes are added to the object:

Parameters:

xs_ps_avg (numpy.ndarray) – 2D array of stacked amplitudes for the Ps phase
xs_pps_avg (numpy.ndarray) – 2D array of stacked amplitudes for the Pps phase
xs_pss_avg (numpy.ndarray) – 2D array of stacked amplitudes for the Pss phase

gccp(wlen=15.0)

Method to average the amplitudes at each grid point to produce 2D images for each of the three phases. In this method, the grid points are further smoothed in the horizontal direction using a Gaussian function to simulate P-wave sensitivity kernels. At the end of this step, the object contains the three 2D smoothed arrays that can be further averaged into a single final image.

Parameters:

wlen (float) – Wavelength of the P-wave for smoothing (km).
object (The following attributes are added to the) –
xs_gauss_ps (numpy.ndarray) – 2D array of stacked and Gaussian-filtered amplitudes for the Ps phase
xs_gauss_pps (numpy.ndarray) – 2D array of stacked and Gaussian-filtered amplitudes for the Pps phase
xs_gauss_pss (numpy.ndarray) – 2D array of stacked and Gaussian-filtered amplitudes for the Pss phase

linear_stack(typ='ccp')

Method to average the three 2D images into a final, weighted CCP image using the weights defined in the attribute.

Parameters:

typ (str) – Type of phase stacks to use (either ccp or gccp)
object (The following attributes are added to the) –
tot_trace (numpy.ndarray) – 2D array of amplitudes for the linearly combined Ps, Pps and Pss phases

phase_weighted_stack(typ='gccp')

Method to average the three 2D smoothed images into a final, phase-weighted CCP image.

Parameters:

typ (str) – Type of phase stacks to use (either ccp or gccp)
object (The following attributes are added to the) –
tot_trace (numpy.ndarray) – 2D array of amplitudes for the phase-weighted, combined Ps, Pps and Pss phases

save(title)

Method to save the ccpimage object to file.

Parameters:: title (str) – File name for the saved object

plot_ccp(vmin=-0.05, vmax=0.05, save=False, title='', fmt='png')

Method to plot the final CCP stacks along the line.

Parameters:

vmin (float) – Maximum negative value for the color scale
vmax (float) – Maximum positive value for the color scale
save (boolean) – Whether or not to save the figure.
fmt (str) – Format of saved figure
xlen (float) – Total length of profile

plot_gccp(vmin=-0.015, vmax=0.015, save=False, title='', fmt='png')

Method to plot the final GCCP stacks along the line.

Parameters:

vmin (float) – Maximum negative value for the color scale
vmax (float) – Maximum positive value for the color scale
save (boolean) – Whether or not to save the figure.
fmt (str) – Format of saved figure
xlen (float) – Total length of profile

Modules

binning

Functions to bin receiver functions along back-azimuth or slowness bins, or both, or for all data, regardless of the x axis (time or depth).

rfpy.binning.bin(stream1, stream2=None, typ='baz', nbin=37, pws=False, include_empty=False)

Function to stack receiver functions into (baz or slow) bins This can be done using a linear stack (i.e., simple mean), or using phase-weighted stacking.

Parameters:

stream1 (Stream) – Stream of equal-length seismograms to be stacked into a single trace.
stream2 (Stream) – Optionally stack a second stream in the same operation.
typ (str) – Attribute to bin ‘baz’: backazimuth (degree) ‘slow’: Horizontal slowness (s/km) ‘dist’: epicentral distance (degree)
nbin (int) – Number of equally sized bins within data range
pws (bool) – Whether or not to perform phase-weighted stacking
include_empty (bool) – Return empty bins as null traces (default omits them)

Returns:

stack – Stream containing one or two stacked traces, depending on the number of input streams

Return type:

Stream

Note

Sets the following attributes of the stack:: nbin: Number of bins index_list: Indices of constituent traces in source stream

rfpy.binning.bin_baz_slow(stream1, stream2=None, nbaz=37, nslow=21, baz_range=(0, 360), slow_range=(0.04, 0.08), pws=False, include_empty=False)

Function to stack receiver functions into back-azimuth and slowness bins. This can be done using a linear stack (i.e., simple mean), or using phase-weighted stacking.

Parameters:

stream1 (Stream) – Stream of equal-length seismograms to be stacked into a single trace.
stream2 (Stream) – Optionally stack a second stream in the same operation.
nbaz (int) – Number of back-azimuth samples in bins
nslow (int) – Number of slowness samples in bins
baz_range (2-tuple) – Minimum and maximum of the back-azimuth bins
slow_range (2-tuple) – Minimum and maximum of the slowness bins
pws (bool) – Whether or not to perform phase-weighted stacking
include_empty (bool) – Return empty bins as null traces (default omits them)

Returns:

stack – Stream containing one or two stacked traces, depending on the number of input streams

Return type:

Stream

Note

Sets the following attributes of the stack:: nbin: Number of bins index_list: Indices of constituent traces in source stream

rfpy.binning.stack(rfdatas, which='rf', pws=False, attributes={'bin': None, 'index_list': None, 'slow': None})

Function to stack receiver functions stored in a RFData list. This can be done using a linear stack (i.e., simple mean), or using phase-weighted stacking.

Parameters:

rfdatas (list of rfpy.rfdata.RFdata) – List of RFData objects to be stacked
which (str) – ‘rf’ stack receiver functions ‘specs’ stack spectra
pws (bool) – Whether or not to perform phase-weighted stacking (which=rf only)
attributes (dict) – Attributes passed to the traces of stack

Returns:

stack – Stream containing stacked traces

Return type:

Stream

rfpy.binning.baz_slow_bin_indices(rfs, nbaz=37, nslow=21, baz_range=(0, 360), slow_range=(0.04, 0.08), include_empty=False)

Sort traces of streams into backazimuth (nbaz) and slowness (nslow) bins.

Parameters:

rfs (list of RFData) – List of receiver functions to be sorted in bins
bazs (int) – Number of back-azimuth samples in bins
slows (int) – Number of slowness samples in bins
baz_range (2-tuple) – Minimum and maximum of the back-azimuth bins
slow_range (2-tuple) – Minimum and maximum of the slowness bins
include_empty (bool) – Include empty bins (default omits them)

Returns:

index_lists (list of lists) – Indices that sorts traces of stream into bins defined by nbaz and nslow
bazslow_list (list of 2*tuple) – Backazimuth and slowness values of each element of index list

rfpy.binning.bin_all(stream1, stream2=None, pws=False)

Function to bin all streams into a single trace. This can be done using a linear stack (i.e., simple mean), or using phase-weighted stacking.

Parameters:

stream1 (Stream) – Stream of equal-length seismograms to be stacked into a single trace.
stream2 (Stream) – Optionally stack a second stream in the same operation.
pws (bool) – Whether or not to perform phase-weighted stacking

Returns:

stack – Stream containing one or two stacked traces, depending on the number of input streams

Return type:

Stream

plotting

Functions to plot single station P-receiver functions as wiggle plots.

Options to plot by receiver functions # vs time/depth, by back-azimuth vs time/depth or slowness vs time.

rfpy.plotting.wiggle(stream1, stream2=None, sort=None, tmin=0.0, tmax=30, normalize=True, save=False, title=None, form='png')

Function to plot receiver function traces by index in stream. By default, only one stream is required, which produces a Figure with a single panel. Optionally, if a second stream is present, the Figure will contain two panels. The stream(s) can be sorted by stats attribute sort, normalized, and the Figure can be saved as .eps.

Parameters:

stream1 (Stream) – Stream of receiver function traces
stream2 (Stream) – Stream of receiver function traces, for a two-panel Figure
sort (str) – Name of attribute in the stats attribute of indivudial traces, used to sort the traces in the Stream(s).
xmax (float) – Maximum x-axis value displayed in the Figure.
normalize (bool) – Whether or not to normalize the traces according to the max amplitude in stream1
save (bool) – Whether or not to save the Figure
title (str) – Title of plot

rfpy.plotting.wiggle_bins(stream1, stream2=None, tr1=None, tr2=None, btyp='baz', tmin=0.0, tmax=30.0, xtyp='time', scale=None, norm=None, save=False, title=None, form='png')

Function to plot receiver function according to either baz or slowness bins. By default, only one stream is required, which produces a Figure with a single panel. Optionally, if a second stream is present, the Figure will contain two panels. If the single trace arguments are present, they will be plotted at the top.

Parameters:

stream1 (Stream) – Stream of receiver function traces
stream2 (Stream) – Stream of receiver function traces, for a two-panel Figure
tr1 (Trace) – Trace to plot at the top of stream1
tr2 (Trace) – Trace to plot at the top of stream2
btyp (str) – Type of plot to produce (either ‘baz’, ‘slow’, or ‘dist’)
tmax (float) – Maximum x-axis value displayed in the Figure.
xtyp (str) – Type of x-axis label (either ‘time’ or ‘depth’)
scale (float) – Scale factor applied to trace amplitudes for plotting
save (bool) – Whether or not to save the Figure
title (str) – Title of plot