我在我的 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
。
原答案:
使用 FFT 和 convolution theorem 将问题从 O(n^2) 转换为 O(n log n),从而显着提高速度。这对于像您这样的长数据集特别有用,并且可以提供 1000 秒或更多的速度增益,具体取决于长度。这也很容易做到:只需对两个信号进行 FFT、乘法和逆 FFT 乘积。 numpy.correlate
在互相关例程中不使用 FFT 方法,更适合用于非常小的内核。
这是一个例子
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 将波形视为周期性的,但如果这是一个问题,可以通过零填充来解决。
关于python - Cython 和 numpy 速度,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/1199972/