6.1.3. subspace¶
This module contains functions relevant to executing subspace detection for earthquake catalogs.
We recommend that you read Harris’ detailed report on subspace detection theory which can be found here: https://e-reports-ext.llnl.gov/pdf/335299.pdf
- copyright:
EQcorrscan developers.
- license:
GNU Lesser General Public License, Version 3 (https://www.gnu.org/copyleft/lesser.html)
Subspace detection for either single-channel cases, or network cases. This is modelled on the methods described by Harris. This method allows for slightly more variation in detected waveforms than the traditional matched-filter method. In this method, templates are constructed either by using the empirical subspace method, or by computing the basis vectors by singular-value decomposition. Both methods are provided as part of EQcorrscan in the clustering module.
6.1.3.1. Classes¶
- class eqcorrscan.core.subspace.Detector(name=None, sampling_rate=None, multiplex=None, stachans=None, lowcut=None, highcut=None, filt_order=None, data=None, u=None, sigma=None, v=None, dimension=None)[source]¶
Class to serve as the base for subspace detections.
- Parameters:
name (str) – Name of subspace detector, used for book-keeping
sampling_rate (float) – Sampling rate in Hz of original waveforms
multiplex (bool) – Is this detector multiplexed.
stachans (list) – List of tuples of (station, channel) used in detector. If multiplexed, these must be in the order that multiplexing was done.
lowcut (float) – Lowcut filter in Hz
highcut (float) – Highcut filter in Hz
filt_order (int) – Number of corners for filtering
data (numpy.ndarray) – The actual detector
u (numpy.ndarray) – Full rank U matrix of left (input) singular vectors.
sigma (numpy.ndarray) – Full rank vector of singular values.
v (numpy.ndarray) – Full rank right (output) singular vectors.
dimension (int) – Dimension of data.
Warning
Changing between scipy.linalg.svd solvers (obvious changes between scipy version 0.18.x and 0.19.0) result in sign changes in svd results. You should only run a detector created using the same scipy version as you currently run.
Methods
construct
(streams, lowcut, highcut, ...[, ...])Construct a subspace detector from a list of streams, full rank.
detect
(st, threshold, trig_int[, moveout, ...])Detect within continuous data using the subspace method.
energy_capture
([stachans, size, show])Calculate the average percentage energy capture for this subspace.
partition
(dimension)Partition subspace into desired dimension.
plot
([stachans, size, show])Plot the output basis vectors for the detector at the given dimension.
read
(filename)Read detector from a file, must be HDF5 format.
write
(filename)Write detector to a file - uses HDF5 file format.
- __init__(name=None, sampling_rate=None, multiplex=None, stachans=None, lowcut=None, highcut=None, filt_order=None, data=None, u=None, sigma=None, v=None, dimension=None)[source]¶
- construct(streams, lowcut, highcut, filt_order, sampling_rate, multiplex, name, align, shift_len=0, reject=0.3, no_missed=True, plot=False)[source]¶
Construct a subspace detector from a list of streams, full rank.
Subspace detector will be full-rank, further functions can be used to select the desired dimensions.
- Parameters:
streams (list) – List of
obspy.core.stream.Stream
to be used to generate the subspace detector. These should be pre-clustered and aligned.lowcut (float) – Lowcut in Hz, can be None to not apply filter
highcut (float) – Highcut in Hz, can be None to not apply filter
filt_order (int) – Number of corners for filter.
sampling_rate (float) – Desired sampling rate in Hz
multiplex (bool) – Whether to multiplex the data or not. Data are multiplexed according to the method of Harris, see the multi function for details.
name (str) – Name of the detector, used for book-keeping.
align (bool) – Whether to align the data or not - needs to be done at some point
shift_len (float) – Maximum shift allowed for alignment in seconds.
reject (float) – Minimum correlation to include traces - only used if align=True.
no_missed (bool) – Reject streams with missed traces, defaults to True. A missing trace from lots of events will reduce the quality of the subspace detector if multiplexed. Only used when multi is set to True.
plot (bool) – Whether to plot the alignment stage or not.
Note
The detector will be normalized such that the data, before computing the singular-value decomposition, will have unit energy. e.g. We divide the amplitudes of the data by the L1 norm of the data.
Warning
EQcorrscan’s alignment will attempt to align over the whole data window given. For long (more than 2s) chunks of data this can give poor results and you might be better off using the
eqcorrscan.utils.stacking.align_traces()
function externally, focusing on a smaller window of data. To do this you would align the data prior to running construct.
- detect(st, threshold, trig_int, moveout=0, min_trig=0, process=True, extract_detections=False, cores=1)[source]¶
Detect within continuous data using the subspace method.
- Parameters:
st (obspy.core.stream.Stream) – Un-processed stream to detect within using the subspace detector.
threshold (float) – Threshold value for detections between 0-1
trig_int (float) – Minimum trigger interval in seconds.
moveout (float) – Maximum allowable moveout window for non-multiplexed, network detection. See note.
min_trig (int) – Minimum number of stations exceeding threshold for non-multiplexed, network detection. See note.
process (bool) – Whether or not to process the stream according to the parameters defined by the detector. Default is True, which will process the data.
extract_detections (bool) – Whether to extract waveforms for each detection or not, if True will return detections and streams.
cores (int) – Number of threads to process data with.
- Returns:
list of
eqcorrscan.core.match_filter.Detection
- Return type:
list
Warning
Subspace is currently in beta, see note in the subspace tutorial for information.
Note
If running in bulk with detectors that all have the same parameters then you can pre-process the data and set process to False. This will speed up this detect function dramatically.
Warning
If the detector and stream are multiplexed then they must contain the same channels and multiplexed in the same order. This is handled internally when process=True, but if running in bulk you must take care.
Note
Non-multiplexed, network detection. When the detector is not multiplexed, but there are multiple channels within the detector, we do not stack the single-channel detection statistics because we do not have a one-size-fits-all solution for computing delays for a subspace detector (if you want to implement one, then please contribute it!). Therefore, these parameters provide a means for declaring a network coincidence trigger using single-channel detection statistics, in a similar fashion to the commonly used network-coincidence trigger with energy detection statistics.
- energy_capture(stachans='all', size=(10, 7), show=False)[source]¶
Calculate the average percentage energy capture for this subspace.
- Returns:
Percentage energy capture
- Return type:
- partition(dimension)[source]¶
Partition subspace into desired dimension.
- Parameters:
dimension (int) – Maximum dimension to use.
- plot(stachans='all', size=(10, 7), show=True, *args, **kwargs)[source]¶
Plot the output basis vectors for the detector at the given dimension.
Corresponds to the first n horizontal vectors of the V matrix.
- Parameters:
stachans (list) – list of tuples of station, channel pairs to plot.
stachans – List of tuples of (station, channel) to use. Can set to ‘all’ to use all the station-channel pairs available. If detector is multiplexed, will just plot that.
size (tuple) – Figure size.
show (bool) – Whether or not to show the figure.
- Returns:
Figure
- Return type:
matplotlib.pyplot.Figure
Note
See
eqcorrscan.utils.plotting.subspace_detector_plot()
for example.
6.1.3.2. Functions¶
- eqcorrscan.core.subspace.read_detector(filename)[source]¶
Read detector from a filename.
- Parameters:
filename (str) – Filename to save the detector to.
- Returns:
Detector object
- Return type:
- eqcorrscan.core.subspace.multi(stream)[source]¶
Internal multiplexer for multiplex_detect.
- Parameters:
stream (obspy.core.stream.Stream) – Stream to multiplex
- Returns:
trace of multiplexed data
- Return type:
Maps a standard multiplexed stream of seismic data to a single traces of multiplexed data as follows:
Input: x = [x1, x2, x3, …] y = [y1, y2, y3, …] z = [z1, z2, z3, …]
Output: xyz = [x1, y1, z1, x2, y2, z2, x3, y3, z3, …]
- eqcorrscan.core.subspace.subspace_detect(detectors, stream, threshold, trig_int, moveout=0, min_trig=1, parallel=True, num_cores=None)[source]¶
Conduct subspace detection with chosen detectors.
- Parameters:
detectors (list) – list of
eqcorrscan.core.subspace.Detector
to be used for detection.stream (obspy.core.stream.Stream) – Stream to detect within.
threshold (float) – Threshold between 0 and 1 for detection, see
Detector.detect()
trig_int (float) – Minimum trigger interval in seconds.
moveout (float) – Maximum allowable moveout window for non-multiplexed, network detection. See note.
min_trig (int) – Minimum number of stations exceeding threshold for non-multiplexed, network detection. See note in
Detector.detect()
.parallel (bool) – Whether to run detectors in parallel in groups.
num_cores (int) – How many cpu cores to use if parallel==True. If set to None (default), will use all available cores.
- Return type:
list
- Returns:
List of
eqcorrscan.core.match_filter.Detection
detections.
Note
This will loop through your detectors using their detect method. If the detectors are multiplexed it will run groups of detectors with the same channels at the same time.