5.5. correlate¶
Correlation functions for multi-channel cross-correlation of seismic data.
Various routines used mostly for testing, including links to a compiled routine using FFTW, a Numpy fft routine which uses bottleneck for normalisation and a compiled time-domain routine. These have varying levels of efficiency, both in terms of overall speed, and in memory usage. The time-domain is the most memory efficient but slowest routine (although fastest for small cases of less than a few hundred correlations), the Numpy routine is fast, but memory inefficient due to a need to store large double-precision arrays for normalisation. The fftw compiled routine is fastest and more memory efficient than the Numpy routine.
copyright: | EQcorrscan developers. |
---|---|
license: | GNU Lesser General Public License, Version 3 (https://www.gnu.org/copyleft/lesser.html) |
5.5.1. Classes & Functions¶
fftw_multi_normxcorr |
Use a C loop rather than a Python loop - in some cases this will be fast. |
fftw_normxcorr |
Normalised cross-correlation using the fftw library. |
numpy_normxcorr |
Compute the normalized cross-correlation using numpy and bottleneck. |
time_multi_normxcorr |
Compute cross-correlations in the time-domain using C routine. |
get_array_xcorr |
Get an normalized cross correlation function that takes arrays as inputs. |
get_stream_xcorr |
Return a function for performing normalized cross correlation on lists of streams. |
register_array_xcorr |
Decorator for registering correlation functions. |
5.6. Selecting a correlation function¶
EQcorrscan strives to use sensible default algorithms for calculating correlation values, however, you may want to change how correlations are caclulated to be more advantageous to your specific needs.
There are currently 3 different correlations functions currently included in EQcorrscanq:
eqcorrscan.utils.correlate.numpy_normxcorr()
known as “numpy”eqcorrscan.utils.correlate.time_multi_normxcorr()
known as “time_domain”eqcorrscan.utils.correlate.fftw_normxcorr()
known as “fftw”
Number 3 is the default.
5.7. Switching which correlation function is used¶
You can switch between correlation functions using the xcorr_func paramater included in:
eqcorrscan.core.match_filter.match_filter()
,
eqcorrscan.core.match_filter.Tribe.detect()
,
eqcorrscan.core.match_filter.Template.detect()
by
- passing a string (eg “numpy”, “time_domain”, or “fftw”)
or
- passing a function
for example:
import obspy
from eqcorrscan.utils.correlate import numpy_normxcorr, set_xcorr
from eqcorrscan.core.match_filter import match_filter
# generate some toy templates and stream
random = np.random.RandomState(42)
template = obspy.read()
stream = obspy.read()
for num, tr in enumerate(stream): # iter stream and embed templates
data = tr.data
tr.data = random.randn(6000) * 5
tr.data[100: 100 + len(data)] = data
# do correlation using numpy rather than fftw
match_filter(['1'], [template], stream, .5, 'absolute', 1, False,
xcorr_func='numpy)
# do correlation using a custom function
def custom_normxcorr(templates, stream, pads, *args, **kwargs):
# Just to keep example short call other xcorr function
print('calling custom xcorr function')
return numpy_normxcorr(templates, stream, pads, *args, **kwargs)
match_filter(['1'], [template], stream, .5, 'absolute', 1, False,
xcorr_func=custom_normxcorr)
# prints "calling custom xcorr function
You can also use the set_xcorr object (eqcorrscan.utils.correlate.set_xcorr) to change which correlation function is used:
# this can also be done using the set_xcorr function / context manager
# to change the default xcorr function for all functions that use it
with set_xcorr(custom_normxcorr):
match_filter(['1'], [template], stream, .5, 'absolute', 1, False)
# prints "calling custom xcorr function"
# you can also permanently change the xcorr function (until your python
# kernel is restarted) by calling set_xcorr
set_xcorr(custom_normxcorr)
match_filter(['1'], [template], stream, .5, 'absolute', 1, False)
# prints "calling custom xcorr function
set_xcorr.revert() # change it back