我在我的python程序中使用cython进行相关计算。我有两个音频数据集,我需要知道它们之间的时差。根据开始时间切割第二组,然后在第一组上滑动。有两个for循环:一个滑动集合,而内部循环则计算该点的相关性。此方法效果很好,并且足够准确。
问题在于,使用纯python会花费超过一分钟的时间。使用我的cython代码,大约需要17秒。这仍然太多了。您是否有任何提示可以加快此代码的速度:
import numpy as np cimport numpy as np cimport cython FTYPE = np.float ctypedef np.float_t FTYPE_t @cython.boundscheck(False) def delay(np.ndarray[FTYPE_t, ndim=1] f, np.ndarray[FTYPE_t, ndim=1] g): cdef int size1 = f.shape[0] cdef int size2 = g.shape[0] cdef int max_correlation = 0 cdef int delay = 0 cdef int current_correlation, i, j # Move second data set frame by frame for i in range(0, size1 - size2): current_correlation = 0 # Calculate correlation at that point for j in range(size2): current_correlation += f[<unsigned int>(i+j)] * g[j] # Check if current correlation is highest so far if current_correlation > max_correlation: max_correlation = current_correlation delay = i return delay
编辑: 现在有scipy.signal.fftconvolve一种方法将成为我在下面描述的基于FFT的卷积方法的首选方法。我将保留原始答案以解释速度问题,但实际上使用scipy.signal.fftconvolve。
scipy.signal.fftconvolve
原始答案: 通过将 FFT 问题从O(n ^ 2)转换为O(n log n),使用 FFT 和 卷积定理 可以显着提高速度。这对于像您这样的长数据集特别有用,并且视长度而定,可以使速度提高1000倍甚至更多。这也很容易做到:只需对两个信号进行FFT,对乘积进行逆FFT,然后对其进行逆FFT。numpy.correlate在互相关例程中不使用FFT方法,最好在很小的内核中使用。
numpy.correlate
这是一个例子
from timeit import Timer from numpy import * times = arange(0, 100, .001) xdata = 1.*sin(2*pi*1.*times) + .5*sin(2*pi*1.1*times + 1.) ydata = .5*sin(2*pi*1.1*times) def xcorr(x, y): return correlate(x, y, mode='same') def fftxcorr(x, y): fx, fy = fft.fft(x), fft.fft(y[::-1]) fxfy = fx*fy xy = fft.ifft(fxfy) return xy if __name__ == "__main__": N = 10 t = Timer("xcorr(xdata, ydata)", "from __main__ import xcorr, xdata, ydata") print 'xcorr', t.timeit(number=N)/N t = Timer("fftxcorr(xdata, ydata)", "from __main__ import fftxcorr, xdata, ydata") print 'fftxcorr', t.timeit(number=N)/N
给出每个周期的运行时间(以秒为单位,表示10,000个长波形)
xcorr 34.3761689901 fftxcorr 0.0768054962158
很明显fftxcorr方法要快得多。
如果将结果绘制出来,您会发现它们在零时移附近非常相似。但是请注意,随着距离的增加,xcorr会减少,而fftxcorr不会。这是因为处理波形移位时不重叠的波形部分有点模棱两可。xcorr将其视为零,而FFT将波形视为周期性,但如果出现问题,则可以通过零填充来解决。