5

相互相関しようとしている 2 つのデータ セットがあります。それらは関数に似ているarctanので、信号処理を行う方法を理解するためのモデルとして使用しています。

x = linspace(-15, 15, 2**13)
f1 = arctan(x)
f2 = arctan(x + 2)

ここに画像の説明を入力

私が答える必要がある質問は、緑の信号を青の信号と (ほとんど) 重なるようにするには、どれだけシフトする必要があるかということです。f1と の相互相関関数で最大値を見つけるのと同じくらい簡単だと思いf2、次のアドバイスに大まかに従いました。. これは私が試してきたことです

c = correlate(f1, f2, 'full')
s = arange(1-2**13, 2**13)
dx = 30/2**13
shift = s[c.argmax()]*dx

shift私は多かれ少なかれ正確に 2 に等しいと予想しますが、実際には0.234. これは私には意味がありません。相互相関の最大値の x 座標を見つけました。これは、2 つの信号の最大の重なりがある場所で見つける必要があります。

この種の関数のこの量を計算する方法についてのアイデアはありますか?

編集：実際のデータについては、すべての値が0から1の間になることを追加する必要があります

EDIT EDIT: 次の関数は、実際には私の実際のデータに似ています。

x = linspace(-15, 15, 400)
f1 = (arctan(-x) + pi/2) / pi
f2 = (arctan(-x + 2) + pi/2) / pi

ここに画像の説明を入力

したがって、ここに示す式を使用して: http://paulbourke.net/miscellaneous/correlate/データをパディングして左に 1 を追加し、右に 0 を追加する相互相関関数を作成できます。

def xcorr(x, y);
    mx = x.mean()
    my = y.mean()
    sx = x.std()
    sy = y.std()
    r = zeros(2*len(x))

    for d in range(-len(x), len(x)):
        csum = 0
        for i in range(0, len(x):
            yindx = i - d
            if i - d < 0:
                yval = 1
            elif i - d >= len(x):
                yval = 0
            else:
                yval = y[yindx]
            csum += (x[i] - mx) * (yval - my)
        r[d + len(x)] = csum / (sx * sy)
    return r

この関数を使用すると、できるようになりました

c = xcorr(f1, f2)
s = arange(-400, 400)
dx = 30/400
shift = s[c.argmax()] * dx

これは 2.025 になり、この精度では 2 に限りなく近づきます。ジェイミーが正しかったようです。問題は、numpycorrelateがシグナルのパディングを行う方法にあります。

したがって、明らかに私のxcorr機能は現状では非常に遅いです。問題は、NumPy に同様のことをさせる方法はありますか、それとも C でアルゴリズムを記述し続ける必要があるかということctypesです。

4

2 に答える 2

3

人々が指摘したように、相互相関はデータを超えるパディングによって混乱します。

良いデータを捨てているように感じますが、多くの場合、データセットをトリミングして、相関関係を仮定なしで実行できるようにする方が良い場合がよくあります (少なくとも、実際のデータをパディング用に作成されたデータと相関させる代替手段と比較した場合)。 .

x = linspace(-15, 15, 4000)
f1 = (arctan(-x) + pi/2) / pi
f2 = (arctan(-x + 2) + pi/2) / pi

L4 = int(len(f2)/8)
sf2 = f2[L4:-L4]

c = correlate(f1-mean(f1), sf2-mean(f1), 'same')
print "peak correlation occurs at:", x[argmax(c)]  # -2.02925731433

plt.plot(x, c)
plt.show()

ここに画像の説明を入力

また、ここで xcorr が最良のアプローチであるかどうかはわかりません。さまざまなシフトの y 軸値間の距離を合計し、最小値を取得するだけで、ゼロがどこにあるかなどのすべての問題から逃れることができます。

于 2013-04-14T04:52:50.507 に答える