Cython和Numpy速度 [英] Cython and numpy speed

问题描述

我在我的python程序中使用cython进行相关计算.我有两个音频数据集，我需要知道它们之间的时差.根据开始时间切割第二组，然后在第一组上滑动.有两个for循环:一个滑动集合，而内部循环则计算该点的相关性.这种方法效果很好，而且足够准确.

I'm using cython for a correlation calculation in my python program. I have two audio data sets and I need to know the time difference between them. The second set is cut based on onset times and then slid across the first set. There are two for-loops: one slides the set and the inner loop calculates correlation at that point. This method works very well and it's accurate enough.

问题在于，使用纯python会花费一分钟以上的时间.使用我的cython代码，大约需要17秒.这仍然太多了.您是否有任何有关如何加快此代码的提示:

The problem is that with pure python this takes more than one minute. With my cython code, it takes about 17 seconds. This still is too much. Do you have any hints how to speed-up this code:

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.

There's now scipy.signal.fftconvolve which would be the preferred approach to doing the FFT based convolution approach that I describe below. I'll leave the original answer to explain the speed issue, but in practice use scipy.signal.fftconvolve.

原始答案:
使用 FFT 和 卷积定理 将大大提高速度.这对于像您这样的长数据集特别有用，并且视长度而定，可以使速度提高1000倍甚至更多.这也很容易做到:只需对两个信号进行FFT，并对乘积进行逆FFT，然后对其进行逆FFT. numpy.correlate在互相关例程中不使用FFT方法，最好在很小的内核中使用.

Original answer:
Using FFTs and the convolution theorem will give you dramatic speed gains by converting the problem from O(n^2) to O(n log n). This is particularly useful for long data sets, like yours, and can give speed gains of 1000s or much more, depending on length. It's also easy to do: just FFT both signals, multiply, and inverse FFT the product. numpy.correlate doesn't use the FFT method in the cross-correlation routine and is better used with very small kernels.

这是一个例子

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个长波形为单位，以秒为单位)

Which gives the running times per cycle (in seconds, for a 10,000 long waveform)

xcorr 34.3761689901
fftxcorr 0.0768054962158

很明显，fftxcorr方法要快得多.

It's clear the fftxcorr method is much faster.

如果绘制结果，您会发现它们在零时移附近非常相似.但是请注意，随着距离的增加，xcorr会减少，而fftxcorr不会.这是因为处理波形移位时不重叠的部分波形有点模棱两可. xcorr将其视为零，而FFT将其视为周期性的波形，但是如果出现问题，可以通过零填充来解决.

If you plot out the results, you'll see that they are very similar near zero time shift. Note, though, as you get further away the xcorr will decrease and the fftxcorr won't. This is because it's a bit ambiguous what to do with the parts of the waveform that don't overlap when the waveforms are shifted. xcorr treats it as zero and the FFT treats the waveforms as periodic, but if it's an issue it can be fixed by zero padding.

这篇关于Cython和Numpy速度的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！