BPMF.template_search

class BPMF.template_search.Beamformer(data=None, network=None, phases=None, travel_times=None, moveouts_relative_to_first=True)

Bases: object

Class for backprojecting waveform features with beamforming.

backproject(waveform_features, reduce='max', device='cpu', out_of_bounds='strict', num_threads=None)

Backproject the waveform features.

Parameters:

• waveform_features ((n_stations, n_components, n_samples) numpy.ndarray) – Features of the waveform time series used for the backprojection onto the grid of theoretical sources.

• device (string, defaults to 'cpu') – Either ‘cpu’ or ‘gpu’, depending on the available hardware and user’s preferences.

• reduce (string, defaults to 'max') – Either ‘max’ or ‘none’. If ‘max’, returns the maximum beam at every time. If ‘none’, returns the full space-time beam.

• out_of_bounds (string, defaults to 'strict') –

  Either ‘strict’ (default) or ‘flexible’. - ‘strict’: A beam is computed if and only if the moveouts point to a

  valid sample (that is, within the bounds of the data stream) for every channel used in the beam.

  • ’flexbile’: A beam is computed as long as the moveouts point to a valid sample for at least one channel. This option is particularly useful for real time applications where an event might have been recorded at the closest stations but not yet at the more distant ones.

• num_threads (int or None, defaults to None) – Number of threads for CPU parallelization. If None (default), uses one thread per available (visible) CPU.

find_detections(detection_threshold, minimum_interevent_time, n_max_stations=None)

Analyze the composite network response to find detections.

Parameters:

• detection_threshold (scalar or (n_samples,) numpy.ndarray, float) – The number of running MADs taken above the running median to define detections.

• minimum_interevent_time (scalar, float) – The shortest duration, in seconds, allowed between two consecutive detections.

• n_max_stations (integer, default to None) – If not None and if smaller than the total number of stations in the network, only extract the n_max_stations closest stations for each theoretical source.

Returns:

detections – Dictionary with data and metadata of the detected earthquakes.

Return type:

dictionary,

property moveouts

property n_phases

property n_sources

property n_stations

property num_sources

plot_detection(detection, figsize=(20, 20), component_aliases={'E': ['E', '2'], 'N': ['N', '1'], 'Z': ['Z']}, n_stations=None)

Plot a detection and the maximum beam.

Parameters:

• detection (dataset.Event) – A dataset.Event instance of a given detection.

• figsize (tuple, default to (20, 20)) – Widht and height of the figure, in inches.

• component_aliases (Dictionary, optional) – Sometimes, components might be named differently than N, E, Z. This dictionary tells the function which alternative component names can be associated with each “canonical” component. For example, component_aliases[‘N’] = [‘N’, ‘1’] means that the function will also check the ‘1’ component in case the ‘N’ component doesn’t exist.

• n_stations (scalar int or None, default to None) – If not None, is the number of stations to plot. The closest stations will be plotted.

Returns:

fig – The Figure instance produced by this method.

Return type:

plt.Figure

plot_likelihood(likelihood=None, time_index=None, **kwargs)

Plot likelihood (beam) slices at a given time.

plot_maxbeam(ax=None, detection=None, **kwargs)

Plot the composite network response.

Parameters:

• ax (plt.Axes, default to None) – If not None, plot in this axis.

• detection (dataset.Event, default to None) – A dataset.Event instance of a given detection.

• figsize (tuple, default to (20, 7)) – Width and height of the figure, in inches.

Returns:

fig – The Figure instance produced by this method.

Return type:

plt.Figure

remove_baseline(window, attribute='composite')

Remove baseline from network response

return_pd_series(attribute='maxbeam')

Return the network response as a Pandas.Series.

set_data(data)

Attribute data to the class instance.

Parameters:

data (dataset.Data) – Instance of dataset.Data.

set_network(network)

Attribute network to the class instance.

network determines which stations are included in the computation of the beams. All arrays with a station axis are ordered according to network.stations.

Parameters:

network (dataset.Network) – The Network instance with the station network information. network can force the network response to be computed only on a subset of the data stored in data.

set_phases(phases)

Attach phases to the class instance.

Parameters:

phases (list of str) – List of strings representing the phase names to read from the travel time table. phases determines the order in which self.moveouts is built.

set_source_coordinates(source_coords)

Attribute _source_coordinates to the class instance.

Parameters:

source_coords (dictionary) – Dictionary with 3 fields: ‘latitude’, ‘longitude’ and ‘depth’

set_travel_times(travel_times)

Attaches travel_times to the class instance.

Parameters:

travel_times (TravelTimes) – TravelTimes instance.

set_weights(weights_phases=None, weights_sources=None)

Set the weights required by beampower.

weights_phases: (n_stations, n_channels, n_phases) np.ndarray, optional

Weight given to each station and channel for a given phase. For example, horizontal components might be given a small or zero weight for the P-wave stacking.

weights_sources: (n_sources, n_stations) np.ndarray, optional

Source-receiver-specific weights. For example, based on the source-receiver distance.

set_weights_sources(n_max_stations, n_min_stations=0)

Set network-geometry-based weights of each source-receiver pair.

Parameters:

n_max_stations (scalar, int) – Maximum number of stations used at each theoretical source. Only the n_max_stations stations will be set a weight > 0.

smooth_maxbeam(window)

Smooth the network response with a gaussian kernel.

property source_coordinates

property stations

weights_station_density(cutoff_dist=None, lower_percentile=0.0, upper_percentile=100.0)

Compute station weights to balance station density.

Areas of high station density produce stronger network responses than areas with sparse coverage, and thus might prevent detecting earthquakes where stations are scarcer.

Parameters:

• cutoff_dist (scalar float, default to None) – All station pairs (i,j) are attributed a number from a gaussian distribution with standard deviation cutoff_dist, in km. The weight of station i is: w_i = 1/(sum_j(exp(-D_ij**2/cutoff_dist**2))). If None, cutoff_dist is set to the median interstation distance.

• lower_percentile (scalar float, default to 0) – If lower_percentile > 0, the weights are clipped above the lower_percentile-th percentile.

• upper_percentile (scalar float, default to 100) – If upper_percentile < 100, the weights are clipped below the upper_percentile-th percentile.

Returns:

weights_sta_density – The station density weights.

Return type:

(n_stations,) numpy.ndarray

class BPMF.template_search.TravelTimes(tt_filename, tt_folder_path='')

Bases: object

Class for handling travel time tables.

convert_to_samples(sampling_rate, remove_tt_seconds=False)

Creates a new self.travel_times_samp attribute.

Parameters:

• sampling_rate (float) – The sampling rate to use to convert times from seconds to samples.

• remove_tt_seconds (boolean, optional) – If True, delete self.travel_times after conversion. This may be necessary to save memory.

get_travel_times_array(units='seconds', stations=None, phases=None, relative_to_first=False)

Parameters:

• units (str, optional) – Either of ‘seconds’ or ‘samples’. If ‘seconds’, build array from self.travel_times. If ‘samples’, build array from self.travel_times_samp.

• stations (list of str, optional) – List of stations to include in the output array. Default is None, which uses all stations.

• phases (list of str, optional) – List of phases to include in the output array. Default is None, which uses all phases.

• relative_to_first (boolean, optional) – If True, the travel times are given as relative to the earliest phase for each source.

property n_sources

property num_sources

property phases

read(phases, source_indexes=None, read_coords=False, stations=None)

Parameters:

• phases (list of str) – List of the seismic phases to read from self.where.

• source_indexes (array-like, optional) – Array-like with the source indexes to read. Default is None (read the whole travel time table).

• read_coords (boolean, optional) – If True, the source coordinates are read from self.where.

• stations (list of str, optional) – If not None, is a list of station names to read a subset of travel times from self.where. Default is None, that is, all stations are read.

property source_coords

property tts

class BPMF.template_search.WaveformTransform(transform_arr, stations, components, starttime, sampling_rate_hz)

Bases: object

Class for waveform transform.

data_frame_view()

Returns self.transform_arr as a pandas.DataFrame.

The data frame’s columns are the component names and the rows are the station names.

property delta

property duration

get_np_array(stations, components=None, verbose=True)

Arguments go to BPMF.utils.get_np_array.

slice(starttime, duration=None, num_samples=None, stations=None, components=None)

Returns a new WaveformTransform instance at the requested time.

Note that either duration or num_samples must be different from None.

Parameters:

• starttime (obspy.UTCDateTime) – Time at which slicing starts.

• duration (float, optional) – If not None, this is the duration, in seconds, of the slice. Defaults to None, in which case duration is calculated from num_samples and self.sampling_rate.

• num_samples (int, optional) – If not None, this is the duration, in samples, of the slice. Defaults to None, in which case num_samples is calculated from duration and self.sampling_rate.

• stations (array-like, optional) – Names of the stations to use in the slice. Defaults to None, in which case stations is taken to be the same as self.stations.

• components (array-like, optional) – Names of the components to use in the slice. Defaults to None, in which case components is taken to be the same as self.components.

Returns:

new_instance – The new WaveformTransform instance sliced from self.

Return type:

BPMF.template_search.WaveformTransform

property sr

property time

BPMF.template_search.baseline(X, w)

BPMF.template_search.envelope(trace)

Compute the envelope of trace.

The envelope is defined as the modulus of the complex analytical signal (a signal whose Fourier transform only has energy in positive frequencies).

Parameters:

trace ((n_samples) numpy.ndarray, float) – The input time series.

Returns:

envelope – The modulus of the analytical signal of the input traces.

Return type:

(n_samples) numpy.ndarray, float

BPMF.template_search.envelope_parallel(traces)

Compute the envelope of traces.

The envelope is defined as the modulus of the complex analytical signal (a signal whose Fourier transform only has energy in positive frequencies).

Parameters:

traces ((n_stations, n_channels, n_samples) numpy.ndarray, float) – The input time series.

Returns:

envelopes – The moduli of the analytical signal of the input traces.

Return type:

(n_stations, n_channels, n_samples) numpy.ndarray, float

BPMF.template_search.saturated_envelopes(traces, anomaly_threshold=1e-11, max_dynamic_range=100000.0)

Compute the saturated envelopes.

Parameters:

• traces ((n_stations, n_components, n_samples) numpy.ndarray) – Input waveform time series.

• anomaly_threshold (scalar, float, default to 1e-11) – Scalar threshold below which the MAD is suspicious. It should be a very small number if you are working on physical unit seismograms.

• max_dynamic_range (scalar, float, default to 1e5) – Higher cutoff on the standardized envelopes. This mitigates the contamination of the network response by transient, undesired high energy signals such as spikes.

BPMF.template_search.time_dependent_threshold(network_response, window, overlap=0.75, CNR_threshold=-10.0)

Compute a time-dependent detection threshold.

Parameters:

• network_response ((n_samples,) numpy.ndarray, float) – Composite network response on which we calculate the detection threshold.

• window (scalar, integer) – Length of the sliding window, in samples, over which we calculate the running statistics used in the detection threshold.

• overlap (scalar, float, default to 0.75) – Ratio of overlap between two contiguous windows.

• CNR_threshold (scalar, float, default to 10) – Number of running MADs above running median that defines the detection threshold.

Returns:

detection_threshold – Detection threshold on the network response.

Return type:

(n_samples,) numpy.ndarray

BPMF.template_search.time_dependent_threshold_pd(network_response, window)

Calculate a time dependent detection threshold using the rolling function from pandas.

Parameters:

• network_response (numpy array,) – Composite network response on which we calculate the detection threshold.

• window (scalar, integer) – Length of the sliding window, in samples, over which we calculate the running statistics used in the detection threshold.

Returns:

detection_threshold – Detection threshold that will serve to select the well backprojected events.

Return type:

numpy array,