5.2.11.1.3. eqcorrscan.utils.pre_processing.shortproc¶
-
eqcorrscan.utils.pre_processing.
shortproc
(st, lowcut, highcut, filt_order, samp_rate, parallel=False, num_cores=False, starttime=None, endtime=None, seisan_chan_names=False, fill_gaps=True, ignore_length=False, ignore_bad_data=False, fft_threads=1)[source]¶ Basic function to bandpass and downsample.
Works in place on data. This is employed to ensure all parts of the data are processed in the same way.
- Parameters
st (obspy.core.stream.Stream) – Stream to process
lowcut (float) – Low cut for bandpass in Hz
highcut (float) – High cut for bandpass in Hz
filt_order (int) – Number of corners for bandpass filter
samp_rate (float) – Sampling rate desired in Hz
parallel (bool) – Set to True to process traces in parallel, for small numbers of traces this is often slower than serial processing, defaults to False
num_cores (int) – Control the number of cores for parallel processing, if set to False then this will use all the cores available.
starttime (obspy.core.utcdatetime.UTCDateTime) – Desired data start time, will trim to this before processing
endtime (obspy.core.utcdatetime.UTCDateTime) – Desired data end time, will trim to this before processing
seisan_chan_names (bool) – Whether channels are named like seisan channels (which are two letters rather than SEED convention of three) - defaults to True.
fill_gaps (bool) – Whether to pad any gaps found with zeros or not.
ignore_length (bool) – Whether to allow data that are less than 80% of the requested length. Defaults to False which will error if short data are found.
ignore_bad_data (bool) – If False (default), errors will be raised if data are excessively gappy or are mostly zeros. If True then no error will be raised, but an empty trace will be returned.
fft_threads (int) – Number of threads to use for pyFFTW FFT in resampling. Note that it is not recommended to use fft_threads > 1 and num_cores > 1.
- Returns
Processed stream
- Return type
Note
If your data contain gaps you should NOT fill those gaps before using the pre-process functions. The pre-process functions will fill the gaps internally prior to processing, process the data, then re-fill the gaps with zeros to ensure correlations are not incorrectly calculated within gaps. If your data have gaps you should pass a merged stream without the fill_value argument (e.g.: st = st.merge()).
Warning
If you intend to use this for processing templates you should consider how resampling will impact your cross-correlations. Minor differences in resampling between day-long files (which you are likely to use for continuous detection) and shorter files will reduce your cross-correlations!
Example, bandpass
>>> from obspy import read >>> from eqcorrscan.utils.pre_processing import shortproc >>> # Get the path to the test data >>> import eqcorrscan >>> import os >>> TEST_PATH = os.path.dirname(eqcorrscan.__file__) + '/tests/test_data' >>> st = read(TEST_PATH + '/WAV/TEST_/2013-09-01-0410-35.DFDPC_024_00') >>> st = shortproc(st=st, lowcut=2, highcut=9, filt_order=3, samp_rate=20, ... parallel=True, num_cores=2) >>> print(st[0]) AF.LABE..SHZ | 2013-09-01T04:10:35.700000Z - 2013-09-01T04:12:05.650000Z | 20.0 Hz, 1800 samples
Example, low-pass
>>> from obspy import read >>> from eqcorrscan.utils.pre_processing import shortproc >>> # Get the path to the test data >>> import eqcorrscan >>> import os >>> TEST_PATH = os.path.dirname(eqcorrscan.__file__) + '/tests/test_data' >>> st = read(TEST_PATH + '/WAV/TEST_/2013-09-01-0410-35.DFDPC_024_00') >>> st = shortproc(st=st, lowcut=None, highcut=9, filt_order=3, ... samp_rate=20) >>> print(st[0]) AF.LABE..SHZ | 2013-09-01T04:10:35.700000Z - 2013-09-01T04:12:05.650000Z | 20.0 Hz, 1800 samples
Example, high-pass
>>> from obspy import read >>> from eqcorrscan.utils.pre_processing import shortproc >>> # Get the path to the test data >>> import eqcorrscan >>> import os >>> TEST_PATH = os.path.dirname(eqcorrscan.__file__) + '/tests/test_data' >>> st = read(TEST_PATH + '/WAV/TEST_/2013-09-01-0410-35.DFDPC_024_00') >>> st = shortproc(st=st, lowcut=2, highcut=None, filt_order=3, ... samp_rate=20) >>> print(st[0]) AF.LABE..SHZ | 2013-09-01T04:10:35.700000Z - 2013-09-01T04:12:05.650000Z | 20.0 Hz, 1800 samples