6.2.4.1.10. eqcorrscan.utils.clustering.extract_detections

eqcorrscan.utils.clustering.extract_detections(detections, templates, archive, arc_type, extract_len=90.0, outdir=None, extract_Z=True, additional_stations=[])[source]

Extract waveforms associated with detections

Takes a list of detections for the template, template. Waveforms will be returned as a list of obspy.core.stream.Stream containing segments of extract_len. They will also be saved if outdir is set. The default is unset. The default extract_len is 90 seconds per channel.

Parameters:
  • detections (list) – List of eqcorrscan.core.match_filter.Detection.

  • templates (list) – A list of tuples of the template name and the template Stream used to detect detections.

  • archive (str) – Either name of archive or path to continuous data, see eqcorrscan.utils.archive_read() for details

  • arc_type (str) – Type of archive, either FDSN, day_vols

  • extract_len (float) – Length to extract around the detection (will be equally cut around the detection time) in seconds. Default is 90.0.

  • outdir (str) – Default is None, with None set, no files will be saved, if set each detection will be saved into this directory with files named according to the detection time, NOT than the waveform start time. Detections will be saved into template subdirectories. Files written will be multiplexed miniseed files, the encoding will be chosen automatically and will likely be float.

  • extract_Z (bool) – Set to True to also extract Z channels for detections delays will be the same as horizontal channels, only applies if only horizontal channels were used in the template.

  • additional_stations (list) – List of tuples of (station, channel) to also extract data for using an average delay.

Returns:

list of obspy.core.streams.Stream

Return type:

list

>>> from eqcorrscan.utils.clustering import extract_detections
>>> from eqcorrscan.core.match_filter import Detection
>>> from obspy import read, UTCDateTime
>>> # Get the path to the test data
>>> import eqcorrscan
>>> import os
>>> TEST_PATH = os.path.dirname(eqcorrscan.__file__) + '/tests/test_data'
>>> # Use some dummy detections, you would use real one
>>> detections = [Detection(
...     template_name='temp1', detect_time=UTCDateTime(2012, 3, 26, 9, 15),
...     no_chans=2, chans=['WHYM', 'EORO'], detect_val=2, threshold=1.2,
...     typeofdet='corr', threshold_type='MAD', threshold_input=8.0),
...               Detection(
...     template_name='temp2', detect_time=UTCDateTime(2012, 3, 26, 18, 5),
...     no_chans=2, chans=['WHYM', 'EORO'], detect_val=2, threshold=1.2,
...     typeofdet='corr', threshold_type='MAD', threshold_input=8.0)]
>>> archive = os.path.join(TEST_PATH, 'day_vols')
>>> template_files = [os.path.join(TEST_PATH, 'temp1.ms'),
...                   os.path.join(TEST_PATH, 'temp2.ms')]
>>> templates = [('temp' + str(i), read(filename))
...              for i, filename in enumerate(template_files)]
>>> extracted = extract_detections(detections, templates,
...                                archive=archive, arc_type='day_vols')
>>> print(extracted[0].sort())
2 Trace(s) in Stream:
AF.EORO..SHZ | 2012-03-26T09:14:15.000000Z - 2012-03-26T09:15:45.000000Z | 1.0 Hz, 91 samples
AF.WHYM..SHZ | 2012-03-26T09:14:15.000000Z - 2012-03-26T09:15:45.000000Z | 1.0 Hz, 91 samples
>>> print(extracted[1].sort())
2 Trace(s) in Stream:
AF.EORO..SHZ | 2012-03-26T18:04:15.000000Z - 2012-03-26T18:05:45.000000Z | 1.0 Hz, 91 samples
AF.WHYM..SHZ | 2012-03-26T18:04:15.000000Z - 2012-03-26T18:05:45.000000Z | 1.0 Hz, 91 samples
>>> # Extract from stations not included in the detections
>>> extracted = extract_detections(
...    detections, templates, archive=archive, arc_type='day_vols',
...    additional_stations=[('GOVA', 'SHZ')])
>>> print(extracted[0].sort())
3 Trace(s) in Stream:
AF.EORO..SHZ | 2012-03-26T09:14:15.000000Z - 2012-03-26T09:15:45.000000Z | 1.0 Hz, 91 samples
AF.GOVA..SHZ | 2012-03-26T09:14:15.000000Z - 2012-03-26T09:15:45.000000Z | 1.0 Hz, 91 samples
AF.WHYM..SHZ | 2012-03-26T09:14:15.000000Z - 2012-03-26T09:15:45.000000Z | 1.0 Hz, 91 samples
>>> # The detections can be saved to a file:
>>> extract_detections(detections, templates, archive=archive,
...                    arc_type='day_vols',
...                    additional_stations=[('GOVA', 'SHZ')], outdir='.')