16

Python プログラムで相関計算に cython を使用しています。2 つのオーディオ データ セットがあり、それらの時間差を知る必要があります。2 番目のセットは、開始時間に基づいてカットされ、最初のセットを横切ってスライドされます。2 つの for ループがあります。1 つはセットをスライドさせ、内側のループはその時点で相関を計算します。この方法は非常にうまく機能し、十分に正確です。

問題は、純粋な python では、これに 1 分以上かかることです。私の 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
4

3 に答える 3

37

編集：以下で説明するFFTベースの畳み込みアプローチを実行するための好ましいアプローチがあります
。scipy.signal.fftconvolve速度の問題を説明するために元の回答を残しますが、実際にはを使用しますscipy.signal.fftconvolve。

元の回答：FFTと畳み込み定理
を 使用すると、問題をO（n ^ 2）からO（n log n）に変換することで、劇的な速度の向上が得られます。これは、あなたのような長いデータセットに特に役立ち、長さに応じて1000秒以上の速度向上をもたらすことができます。また、簡単に実行できます。両方の信号をFFTし、乗算し、逆FFTするだけです。相互相関ルーチンで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は波形を周期的なものとして扱いますが、それが問題である場合は、ゼロパディングによって修正できます。

于 2009-08-04T17:40:00.020 に答える
2

この種のことの秘訣は、分割して征服する方法を見つけることです.

現在、すべての位置にスライドし、すべての位置ですべてのポイントをチェックしています - 事実上O ( n ^ 2 ) 操作です。

不一致を判断するには、すべてのポイントのチェックとすべての位置の比較を減らす必要があります。

たとえば、「これでも近いですか?」という短い文を使用できます。最初のいくつかの位置をチェックするフィルター。相関関係があるしきい値を超えている場合は、続行してください。それ以外の場合は、あきらめて先に進みます。

8 を掛ける「8 桁ごとにチェック」を行うことができます。これが低すぎる場合は、スキップして先に進みます。これが十分に高い場合は、すべての値をチェックして、最大値が見つかったかどうかを確認します。

問題は、これらすべての乗算を実行するのに必要な時間です -- ( f[<unsigned int>(i+j)] * g[j]) 実際には、これらすべての積を大きな行列に入力し、合計が最大の行を選択しています。「すべて」の製品を計算する必要はありません。最大金額を確実に見つけられるように、十分な数の製品を用意してください。

最大値を見つける際の問題は、それが最大かどうかを確認するためにすべてを合計する必要があることです。これを最小化問題に変えることができれば、中間結果がしきい値を超えると、積と和の計算を簡単に放棄できます。

（これはうまくいくと思います。私は試していません。）

負の数を扱っていた場合max(g)-g[j]は、最大ではなく最小を探すことになります。最初の位置の相関を計算できます。より大きな値に合計されたものはすべてすぐに停止できます。そのオフセットに対して乗算または加算を行う必要はなく、別のオフセットにシフトします。

于 2009-07-29T12:58:22.203 に答える
2
  • 外部ループから range(size2) を抽出できます
  • ループの代わりに sum() を使用して current_correlation を計算できます
  • 相関と遅延をリストに保存し、 max() を使用して最大のものを取得できます
于 2009-07-29T13:04:13.700 に答える