API Docs

rfstream Module

Classes and functions for receiver function calculation.

class rf.rfstream.RFStream(traces=None)

Bases: Stream

Class providing the necessary functions for receiver function calculation.

Parameters:

traces – list of traces, single trace or stream object

To initialize a RFStream from a Stream object use

>>> rfstream = RFStream(stream)

To initialize a RFStream from a file use

>>> rfstream = read_rf('test.SAC')

Format specific headers are loaded into the stats object of all traces.

deconvolve(*args, **kwargs)

Deconvolve source component of stream.

All args and kwargs are passed to the function deconvolve().

harmonics(*args, **kwargs)

Perform harmonic decomposition on stream.

All args and kwargs are passed to the function harmonics().

property method

Property for used rf method, ‘P’ or ‘S’

moveout(phase=None, ref=6.4, model='iasp91')

In-place moveout correction to a reference slowness.

Needs stats attributes slowness and onset.

Parameters:

• phase – ‘Ps’, ‘Sp’, ‘Ppss’ or other multiples, if None is set to ‘Ps’ for P receiver functions or ‘Sp’ for S receiver functions

• ref – reference ray parameter in s/deg

• model – Path to model file (see SimpleModel, default: iasp91)

plot_profile(*args, **kwargs)

Create receiver function profile plot.

See imaging.plot_profile() for help on arguments.

plot_rf(*args, **kwargs)

Create receiver function plot.

See imaging.plot_rf() for help on arguments.

ppoints(pp_depth, pp_phase=None, model='iasp91')

Return coordinates of piercing point calculated by 1D ray tracing.

Piercing point coordinates are stored in the stats attributes plat and plon. Needs stats attributes station_latitude, station_longitude, slowness and back_azimuth.

Parameters:

• pp_depth – depth of interface in km

• pp_phase – ‘P’ for piercing points of P wave, ‘S’ for piercing points of S wave or multiples, if None will be set to ‘S’ for P receiver functions or ‘P’ for S receiver functions

• model – path to model file (see SimpleModel, default: iasp91)

Returns:

NumPy array with coordinates of piercing points

profile(*args, **kwargs)

Return profile of receiver functions in the stream.

See profile.profile() for help on arguments.

rf(method=None, filter=None, trim=None, downsample=None, rotate='ZNE->LQT', deconvolve='time', source_components=None, **kwargs)

Calculate receiver functions in-place.

Parameters:

• method – ‘P’ for P receiver functions, ‘S’ for S receiver functions, if None method will be determined from the phase

• filter (dict) – filter stream with its filter method and given kwargs

• trim (tuple (start, end)) – trim stream relative to P- or S-onset with trim2() (seconds)

• downsample (float) – downsample stream with its decimate() method to the given frequency

• rotate ('ZNE->LQT' or 'NE->RT') – rotate stream with its rotate() method with the angles given by the back_azimuth and inclination attributes of the traces stats objects. You can set these to your needs or let them be computed by rfstats().

• deconvolve – ‘time’, ‘waterlevel’, ‘iterative’ or ‘multitaper’ for time domain damped, frequency domain water level, time domain iterative, or frequency domain multitaper deconvolution using the stream’s deconvolve() method. See deconvolve(), deconv_time(), deconv_waterlevel(), deconv_iterative(), and deconv_multitaper() for further documentation.

• source_components – parameter is passed to deconvolve. If None, source components will be chosen depending on method.

• **kwargs – all other kwargs not mentioned here are passed to deconvolve

After performing the deconvolution the Q/R and T components are multiplied by -1 to get a positive phase for a Moho-like positive velocity contrast. Furthermore for method=’S’ all components are mirrored at t=0 for a better comparison with P receiver functions. See source code of this function for the default deconvolution windows.

slice2(starttime=None, endtime=None, reftime=None, keep_empty_traces=False, **kwargs)

Alternative slice method accepting relative times.

See slice() and trim2().

stack()

Return stack of traces with same seed ids.

Traces with same id need to have the same number of datapoints. Each trace in the returned stream will correspond to one unique seed id.

trim2(starttime=None, endtime=None, reftime=None, **kwargs)

Alternative trim method accepting relative times.

See trim().

Parameters:

• starttime,endtime – accept UTCDateTime or seconds relative to reftime

• reftime – reference time, can be an UTCDateTime object or a string. The string will be looked up in the stats dictionary (e.g. ‘starttime’, ‘endtime’, ‘onset’).

property type

Property for the type of stream, ‘rf’, ‘profile’ or None

write(filename, format, sh_compat=False, **kwargs)

Save stream to file including format specific headers.

See Stream.write() in ObsPy.

Parameters:

sh_compat – Ensure files in Q format can be read with SeismicHandler (default: False). If set to True the COMMENT field will be deleted which might result in loss of some metadata. (see issue #37).

class rf.rfstream.RFTrace(data=array([], dtype=float64), header=None, trace=None)

Bases: Trace

Class providing the Trace object for receiver function calculation.

write(filename, format, **kwargs)

Save current trace into a file including format specific headers.

See Trace.write() in ObsPy.

rf.rfstream.obj2stats(event=None, station=None)

Map event and station object to stats with attributes.

Parameters:

• event – ObsPy Event object

• station – station object with attributes latitude, longitude and elevation

Returns:

stats object with station and event attributes

rf.rfstream.read_rf(pathname_or_url=None, format=None, **kwargs)

Read waveform files into RFStream object.

See read() in ObsPy.

rf.rfstream.rfstats(obj=None, event=None, station=None, phase='P', dist_range='default', tt_model='iasp91', pp_depth=None, pp_phase=None, model='iasp91')

Calculate ray specific values like slowness for given event and station.

Parameters:

• obj – Stats object with event and/or station attributes. Can be None if both event and station are given. It is possible to specify a stream object, too. Then, rfstats will be called for each Trace.stats object and traces outside dist_range will be discarded.

• event – ObsPy Event object

• station – dictionary like object with items latitude, longitude and elevation

• phase – string with phase. Usually this will be ‘P’ or ‘S’ for P and S receiver functions, respectively.

• dist_range (tuple of length 2) –

  if epicentral of event is not in this intervall, None is returned by this function,

  if phase == ‘P’ defaults to (30, 90),

  if phase == ‘S’ defaults to (50, 85)

• tt_model – model for travel time calculation. (see the obspy.taup module, default: iasp91)

• pp_depth – Depth for piercing point calculation (in km, default: None -> No calculation)

• pp_phase – Phase for pp calculation (default: ‘S’ for P-receiver function and ‘P’ for S-receiver function)

• model – Path to model file for pp calculation (see SimpleModel, default: iasp91)

Returns:

Stats object with event and station attributes, distance, back_azimuth, inclination, onset and slowness or None if epicentral distance is not in the given interval. Stream instance if stream was specified instead of stats.

deconvolve Module

Frequency and time domain deconvolution.

rf.deconvolve.deconv_iterative(rsp, src, sampling_rate, tshift=10, gauss=0.5, itmax=400, minderr=0.001, mute_shift=False, normalize=0)

Iterative deconvolution.

Deconvolve src from arrays in rsp. Iteratively construct a spike train based on least-squares minimzation of the difference between one component of an observed seismogram and a predicted signal generated by convolving the spike train with an orthogonal component of the seismogram.

Reference: Ligorria, J. P., & Ammon, C. J. (1999). Iterative Deconvolution and Receiver-Function Estimation. Bulletin of the Seismological Society of America, 89, 5.

Parameters:

• rsp – either a list of arrays containing the response functions or a single array

• src – array with source function

• sampling_rate – sampling rate of the data

• tshift – delay time 0s will be at time tshift afterwards

• gauss – Gauss parameter (standard deviation) of the Gaussian Low-pass filter, corresponds to cut-off frequency in Hz for a response value of exp(-0.5)=0.607.

• itmax – limit on number of iterations/spikes to add

• minderr – stop iteration when the change in error from adding another spike drops below this threshold

• mute_shift – Mutes all samples at beginning of trace (lenght given by time shift). For len(src)==len(rsp) this mutes all samples before the onset.

• normalize – normalize all results so that the maximum of the trace with the supplied index is 1. Set normalize to None for no normalization.

Returns:

(list of) array(s) with deconvolution(s)

rf.deconvolve.deconv_multitaper(rsp, src, nse, sampling_rate, tshift, gauss=0.5, K=3, tband=4, T=10, olap=0.75, normalize=0)

Multitaper frequency domain deconvolution

Deconvolve src from arrays in rsp.

Based on mtdecon.f by G. Helffrich with some slight modifications

References: Helffrich, G (2006). Extended-time multitaper frequency domain cross-correlation receiver function estimation. Bulletin of the Seismological Society of America, 96 (1). Shibutani, T., Ueno, T., & Hirahara, K. (2008). Improvement in the Extended-Time Multitaper Receiver Function Estimation Technique. Bulletin of the Seismological Society of America, 98 (2).

Parameters:

• rsp – a list of arrays containing the response functions

• src – array of source function

• nse – a list of arrays containing samples of pre-event noise

• sampling_rate – sampling rate of the data

• tshift – shift the source by that amount of samples to the left side to get onset in RF at the desired time (src time window length pre-onset)

• gauss – Gauss parameter (standard deviation) of the Gaussian Low-pass filter, corresponds to cut-off frequency in Hz for a response value of exp(-0.5)=0.607.

• K – number of Slepian tapers to use (default: 3)

• tband – time-bandwidth product (default: 4)

• T – time length of taper window in seconds (default: 10)

• olap – window overlap between 0 and 1 (default: 0.75)

• normalize – normalize all results so that the maximum of the trace with supplied index is 1. Set normalize to None for no normalization.

Returns:

(list of) array(s) with deconvolution(s)

rf.deconvolve.deconv_time(rsp_list, src, shift, spiking=1.0, length=None, normalize=0, solve_toeplitz='toeplitz')

Time domain deconvolution.

Deconvolve src from arrays in rsp_list. Calculate Toeplitz auto-correlation matrix of source, invert it, add noise and multiply it with cross-correlation vector of response and source.

In one formula:

RF = (STS + spiking*I)^-1 * STR

N... length
    ( S0   S-1  S-2 ... S-N+1 )
    ( S1   S0   S-1 ... S-N+2 )
S = ( S2   ...                )
    ( ...                     )
    ( SN-1 ...          S0    )
R = (R0 R1 ... RN-1)^T
RF = (RF0 RF1 ... RFN-1)^T
S... source matrix (shape N*N)
R... response vector (length N)
RF... receiver function (deconvolution) vector (length N)
STS = S^T*S = symmetric Toeplitz autocorrelation matrix
STR = S^T*R = cross-correlation vector
I... Identity

Parameters:

• rsp_list – either a list of arrays containing the response functions or a single array

• src – array of source function

• shift –

  shift the source by that amount of samples to the left side to get onset in RF at the desired time (negative -> shift source to the right side)

  shift = (middle of rsp window - middle of src window) + (0 - middle rf window)

• spiking – random noise added to autocorrelation (eg. 1.0, 0.1)

• length – number of data points in results

• normalize – normalize all results so that the maximum of the trace with supplied index is 1. Set normalize to None for no normalization.

• solve_toeplitz – Python function to use as a solver for Toeplitz system. Possible values are 'toeplitz' for using the toeplitz package (default), 'scipy' for using solve_toeplitz from the SciPy package, or a custom function.

Returns:

(list of) array(s) with deconvolution(s)

rf.deconvolve.deconv_waterlevel(rsp_list, src, sampling_rate, waterlevel=0.05, gauss=0.5, tshift=10.0, length=None, normalize=0, nfft=None, return_info=False)

Frequency-domain deconvolution using waterlevel method.

Deconvolve src from arrays in rsp_list.

Parameters:

• rsp_list – either a list of arrays containing the response functions or a single array

• src – array with source function

• sampling_rate – sampling rate of the data

• waterlevel – waterlevel to stabilize the deconvolution

• gauss – Gauss parameter (standard deviation) of the Gaussian Low-pass filter, corresponds to cut-off frequency in Hz for a response value of exp(-0.5)=0.607.

• tshift – delay time 0s will be at time tshift afterwards

• nfft – explicitely set number of samples for the fft

• length – number of data points in results, optional

• normalize – normalize all results so that the maximum of the trace with supplied index is 1. Set normalize to ‘src’ to normalize for the maximum of the prepared source. Set normalize to None for no normalization.

• return_info – return additionally a lot of different parameters in a dict for debugging purposes

Returns:

(list of) array(s) with deconvolution(s)

rf.deconvolve.deconvolve(stream, method='time', func=None, source_components='LZ', response_components=None, winsrc='P', **kwargs)

Deconvolve one component of a stream from other components.

The deconvolutions are written to the data arrays of the stream. To keep the original data use the copy method of Stream. The stats dictionaries of the traces inside stream must have an ‘onset’ entry with a UTCDateTime object. This will be used for determining the data windows.

Parameters:

• stream – stream including responses and source

• method –

  ‘time’ -> use time domain deconvolution, see deconv_time(),

  ’waterlevel’ -> use frequency domain deconvolution with water level, see deconv_waterlevel()

  ’iterative’ -> use iterative time domain deconvolution, see deconv_iterative()

  ’multitaper’ -> use frequency domain multitaper deconvolution, see deconv_multitaper()

  ’func’ -> user defined function (func keyword)

• func –

  Custom deconvolution function with the following signature

  def deconv_custom(rsp: RFStream, src: RFTrace, tshift=10,

  **other_kwargs_possible) -> RFStream:

• source_components – names of components identifying the source traces, e.g. ‘LZ’ for P receiver functions and ‘QR’ for S receiver functions

• response_components – names of components identifying the response traces (default None: all traces are used as response)

• winsrc (tuple (start, end, taper)) –

  data window for source function, in seconds relative to onset. The function will be cosine tapered at both ends (seconds).

  winsrc can also be a string (‘P’ or ‘S’). In this case the function defines a source time window appropriate for this type of receiver function and deconvolution method (see source code for details).

• **kwargs – other kwargs are passed to the underlying deconvolution functions

Note

If parameter normalize is not present in kwargs and source component is not excluded from the results by response_components, results will be normalized such that the maximum of the deconvolution of the trimmed and tapered source from the untouched source is 1. If the source is excluded from the results, the normalization will performed against the first trace in results.

Note

If multitaper deconvolution is used and a stream of (pre-event) noise is not present in kwargs, noise will be sampled from the data in a pre-event window whose length depends on the trace length prior to the onset time.

imaging Module

Functions for receiver function plotting.

rf.imaging.plot_harmonics(hd, hd2=None, fillcolors=('b', 'r'), trim=None)

Plot components from harmonic decomposition.

The plot will have two panels. In all cases the left panel will show the modeled components of hd. If hd2 is supplied, the right panel will show the modeled components of hd2. If hd2 is None and hd has unmodeled components, those will be on the right. If hd2 is None and hd does not have unmodeled components, the right panel will be empty.

Parameters:

• hd – RFStream of harmonics, returned from rf.harmonics.harmonics() or the corresponding harmonics() method. The stream should contain modeled components, and may also include unmodeled components

• hd2 – Optional second RFStream of harmonics. If used, it should contain modeled components

• fillcolors – fill colors for positive and negative wiggles

• trim – trim stream relative to onset before plotting using slice2()

rf.imaging.plot_ppoints(ppoints, inventory=None, label_stations=True, ax=None, crs=None, **kwargs)

Plot piercing points with stations.

Parameters:

• ppoints – list of (lat, lon) tuples of piercing points

• label_stations (inventory,) – plot stations, see plot_stations

• ax – geoaxes (default None: new ax will be created)

• crs – coordinate reference system for new geoaxis, (default: None, then AzimuthalEquidistant projection with appropriate center is used.)

• **kwargs – other kwargs are passed to ax.scatter() call

rf.imaging.plot_profile(profile, fname=None, figsize=None, dpi=None, scale=1, fillcolors=('C3', 'C0'), trim=None, top=None, moveout_model='iasp91')

Plot receiver function profile.

Parameters:

• profile – stream holding the profile

• fname – filename to save plot to. Can be None. In this case the figure is left open.

• figsize – figsize of the created figure

• dpi – dots per inch for the created figure

• scale – scale for individual traces

• fillcolors – fill colors for positive and negative wiggles

• trim – trim stream relative to onset before plotting using slice2()

• top – show second axes on top of profile with additional information. Valid values: ‘hist’ - Plot histogram showing the number of receiver functions stacked in the corresponding bin

• moveout_model – string with model filename. Will be loaded into a SimpleModel object to calculate depths for tick labels.

rf.imaging.plot_profile_map(boxes, inventory=None, label_stations=True, ppoints=None, ax=None, crs=None, **kwargs)

Plot profile map with stations and piercing points.

Parameters:

• boxes – boxes created with get_profile_boxes()

• label_stations (inventory,) – plot stations, see plot_stations

• ppoints – list of (lat, lon) tuples of piercing points, see plot_ppoints

• ax – geoaxes (default None: new ax will be created)

• crs – coordinate reference system for new geoaxis, (default: None, then AzimuthalEquidistant projection with appropriate center is used.)

• **kwargs – other kwargs are passed to ax.add_geometries() call

rf.imaging.plot_rf(stream, fname=None, fig_width=7.0, trace_height=0.5, stack_height=0.5, dpi=None, trace_scale=False, scale_factor=1, fillcolors=(None, None), trim=None, info=(('back_azimuth', 'baz (°)', 'C0'), ('distance', 'dist (°)', 'C3')), show_traces=True, show_vlines=False)

Plot receiver functions.

Parameters:

• stream – stream to plot

• fname – filename to save plot to. Can be None. In this case the figure is left open.

• fig_width – width of figure in inches

• trace_height – height of one trace in inches

• stack_height – height of stack axes in inches

• dpi – dots per inch for the created figure

• trace_scale – if True, scale each trace individually; if false, scale traces by global max amplitude

• scale_factor – factor for scaling traces, either individually or globally. A value of 1 generally works well.

• fillcolors – fill colors for positive and negative wiggles

• trim – trim stream relative to onset before plotting using slice2()

• info – Plot one additional axes showing maximal two entries of the stats object. Each entry in this list is a list consisting of three entries: key, label and color. info can be None. In this case no additional axes is plotted.

• show_vlines – If True, show vertical alignment grid lines on plot at positions of the major x-tick marks.

• show_traces – If True, plot the individual traces in the stream in an additional set of axes below the plot of the stacked trace. If False, info will also be set to None and the only thing plotted is the stacked trace.

rf.imaging.plot_stations(inventory, label_stations=True, ax=None, crs=None, **kwargs)

Plot stations.

Parameters:

• inventory – station inventory

• label_stations – weather to label stations

• ax – geoaxes (default None: new ax will be created)

• crs – coordinate reference system for new geoaxis, (default: None, then AzimuthalEquidistant projection with appropriate center is used.)

• **kwargs – other kwargs are passed to ax.scatter() call

profile Module

Functions for receiver function profile calculation.

rf.profile.get_profile_boxes(latlon0, azimuth, bins, width=2000)

Create 2D boxes for usage in profile() function.

Parameters:

• latlon0 (tuple) – coordinates of starting point of profile

• azimuth – azimuth of profile direction

• bins (tuple) – Edges of the distance bins in km (e.g. (0, 10, 20, 30))

• width – width of the boxes in km (default: large)

Returns:

List of box dicts. Each box has the entries ‘poly’ (shapely polygon with lonlat corners), ‘length’ (length in km), ‘pos’ (midpoint of box in km from starting coordinates), ‘latlon’ (midpoint of box as coordinates)

rf.profile.profile(stream, boxes, crs=None)

Stack traces in stream by piercing point coordinates in defined boxes.

Parameters:

• stream – stream with pre-calculated piercing point coordinates

• boxes – boxes created with get_profile_boxes()

• crs – cartopy projection (default: AzimuthalEquidistant)

Returns:

profile stream

simple_model Module

Simple move out and piercing point calculation.

class rf.simple_model.SimpleModel(z, vp, vs, n=None)

Bases: object

Simple 1D velocity model for move out and piercing point calculation.

Parameters:

• z – depths in km

• vp – P wave velocities at provided depths in km/s

• vs – S wave velocities at provided depths in km/s

• n – number of support points between provided depths

All arguments can be of type numpy.ndarray or list.

calculate_delay_times(slowness, phase='PS')

Calculate delay times between direct wave and converted phase.

Parameters:

• slowness – ray parameter in s/deg

• phase – Converted phase or multiple (e.g. Ps, Pppp)

Returns:

delay times at different depths

calculate_vertical_slowness(slowness, phase='PS')

Calculate vertical slowness of P and S wave.

Parameters:

• slowness – slowness in s/deg

• phase – Weather to calculate only P, only S or both vertical slownesses

Returns:

vertical slowness of P wave, vertical slowness of S wave at different depths (z attribute of model instance)

moveout(stream, phase='Ps', ref=6.4)

In-place moveout correction to reference slowness.

Parameters:

• stream – stream with stats attributes onset and slowness.

• phase – ‘Ps’, ‘Sp’, ‘Ppss’ or other multiples

• ref – reference slowness (ray parameter) in s/deg

ppoint(stats, depth, phase='S')

Calculate latitude and longitude of piercing point.

Piercing point coordinates and depth are saved in the pp_latitude, pp_longitude and pp_depth entries of the stats object or dictionary.

Parameters:

• stats – Stats object or dictionary with entries slowness, back_azimuth, station_latitude and station_longitude

• depth – depth of interface in km

• phase – ‘P’ for piercing point of P wave, ‘S’ for piercing point of S wave. Multiples are possible, too.

Returns:

latitude and longitude of piercing point

ppoint_distance(depth, slowness, phase='S')

Calculate horizontal distance between piercing point and station.

Parameters:

• depth – depth of interface in km

• slowness – ray parameter in s/deg

• phase – ‘P’ or ‘S’ for P wave or S wave. Multiples are possible.

Returns:

horizontal distance in km

stretch_delay_times(slowness, phase='Ps', ref=6.4)

Stretch delay times of provided slowness to reference slowness.

First, calculate delay times (time between the direct wave and the converted phase or multiples at different depths) for the provided slowness and reference slowness. Secondly, stretch the the delay times of provided slowness to reference slowness.

Parameters:

• slowness – ray parameter in s/deg

• phase – ‘Ps’, ‘Sp’ or multiples

• ref – reference ray parameter in s/deg

Returns:

original delay times, delay times stretched to reference slowness

rf.simple_model.load_model(fname='iasp91')

Load model from file.

Parameters:

fname – path to model file or ‘iasp91’

Returns:

SimpleModel instance

The model file should have 4 columns with depth, vp, vs, n. The model file for iasp91 starts like this:

#IASP91 velocity model
#depth  vp    vs   n
 0.00  5.800 3.360 0
 0.00  5.800 3.360 0
10.00  5.800 3.360 4

harmonics Module

Harmonic decomposition

rf.harmonics.decomp_one(baz, azim=0)

Build jacobian for harmonic decomposition for one RF component.

Parameters:

• baz – array of back azimuths for RFs in degrees

• azim – azimuth along which to decompose the RFs. This can be used with some kind of optimization minimize components

Returns:

jacobian matrix as np.ndarray, 5xN where N is time series length

rf.harmonics.decomp_two(baz, azim=0, scalars=(3, 0.3))

Build jacobian for harmonic decomposition with two components of RFs.

Parameters:

• baz – array of back azimuths for RFs in degrees

• azim – azimuth along which to decompose the RFs. This can be used with some kind of optimization minimize components

• scalars – scalar multipliers for R/Q and T constant terms. Park and Levin (2016) recommend 3 and 0.3 to bias toward radial component conversions, especially for datasets with patchy back-azimuthal coverage. (1,1) would weight Q/R and T equally in the regression.

Returns:

jacobian matrix as np.ndarray, 10xN where N is the time series length

rf.harmonics.harmonics(stream, components='R', azim=0, method='time', **kwargs)

Perform harmonic decomposition of PRFs.

This implements the method described in Bianchi et al 2010, doi:10.1029/2009JB007061 (equations in supplemental material) and in Park and Levin 2016, doi:10.1093/gji/ggw323 (equations 44/45 and 47).

The harmonic components are returned in an array with rows for constant, cos, sin, cos2, and sin2 terms. The stats dictionaries of the traces in the input stream must have a ‘back_azimuth’ entry and an ‘event_time’ entry.

Parameters:

• stream – RFStream including components to be decomposed

• components – names of components to use in decomposition; can be R, RT, Q, QT, or T (for PRFs)

• azim – azimuth along which to decompose the RFs (default: 0)

• method – domain for harmonic decomposition. Options are ‘time’ -> time domain (Bianchi et al. 2010) ‘freq’ -> frequency domain (Park and Levin 2016)

• **kwargs – other kwargs are passed to underlying functions for building jacobians

Returns:

harmonic components in a Stream

Note

The number of traces returned depends on the number of components. Two-component decomposition returns 10 harmonics: 5 modeled, 5 unmodeled. One-component decomposition returns only 5 modeled harmonics for that single component.

util Module

Utility functions and classes for receiver function calculation.

rf.util.DEG2KM = 111.2

Conversion factor from degrees epicentral distance to km

class rf.util.IterMultipleComponents(stream, key=None, number_components=None)

Bases: object

Return iterable to iterate over associated components of a stream.

Parameters:

• stream – Stream with different, possibly many traces. It is split into substreams with the same seed id (only last character i.e. component may vary)

• key (str or None) – Additionally, the stream is grouped by the values of the given stats entry to differentiate between e.g. different events (for example key=’starttime’, key=’onset’)

• number_components (int, tuple of ints or None) – Only iterate through substreams with matching number of components.

rf.util.direct_geodetic(latlon, azi, dist)

Solve direct geodetic problem with geographiclib.

Parameters:

• latlon (tuple) – coordinates of first point

• azi – azimuth of direction

• dist – distance in km

Returns:

coordinates (lat, lon) of second point on a WGS84 globe

rf.util.iter_event_data(events, inventory, get_waveforms, phase='P', request_window=None, pad=10, pbar=None, **kwargs)

Return iterator yielding three component streams per station and event.

Parameters:

• events – list of events or Catalog instance

• inventory – Inventory instance with station and channel information

• get_waveforms – Function returning the data. It has to take the arguments network, station, location, channel, starttime, endtime.

• phase – Considered phase, e.g. ‘P’, ‘S’, ‘PP’

• request_window (tuple (start, end)) – requested time window around the onset of the phase

• pad (float) – add specified time in seconds to request window and trim afterwards again

• pbar – tqdm instance for displaying a progressbar

• kwargs – all other kwargs are passed to rfstats()

Returns:

three component streams with raw data

Example usage with progressbar:

from tqdm import tqdm
from rf.util import iter_event_data
with tqdm() as t:
    for stream3c in iter_event_data(*args, pbar=t):
        do_something(stream3c)

rf.util.iter_event_metadata(events, inventory, pbar=None)

Return iterator yielding metadata per station and event.

Parameters:

• events – list of events or Catalog instance

• inventory – Inventory instance with station and channel information

• pbar – tqdm instance for displaying a progressbar

rf.util.minimal_example_Srf()

Return S receiver functions calculated from example data.

rf.util.minimal_example_rf()

Return receiver functions calculated from the data returned by read_rf().

batch Module

rf: receiver function calculation - batch command line utility

class rf.batch.ConfigJSONDecoder(*, object_hook=None, parse_float=None, parse_int=None, parse_constant=None, strict=True, object_pairs_hook=None)

Bases: JSONDecoder

Strip lines from comments.

decode(s)

Return the Python representation of s (a str instance containing a JSON document).

exception rf.batch.ParseError

Bases: Exception

rf.batch.init_data(data, client_options=None, plugin=None)

Return appropriate get_waveforms function.

See example configuration file for a description of the options.

rf.batch.iter_event_processed_data(events, inventory, pin, format, yield_traces=False, pbar=None)

Iterator yielding streams or traces which are read from disc.

rf.batch.load_func(modulename, funcname)

Load and return function from Python module.

rf.batch.run(command, conf=None, tutorial=False, **kw)

Create example configuration file and tutorial or load config.

After that call run_commands.

rf.batch.run_cli(args=None)

Command line interface of rf.

After parsing call run.

rf.batch.run_commands(command, commands=(), events=None, inventory=None, objects=None, get_waveforms=None, data=None, plugin=None, phase=None, moveout_phase=None, path_in=None, path_out=None, format='Q', newformat=None, **kw)

Load files, apply commands and write result files.

rf.batch.tqdm()

rf.batch.write(stream, root, format, type=None)

Write stream to one or more files depending on format.