2020/12/21

DSP 相關 Correlation

DSP 相關 Correlation

包含交互相關 cross-correlation 以及 自相關 autocorrelation 及其應用。

correlation 相關就是兩個隨機變數之間,是否具有關聯性。相關係數可用來測量兩個隨機變數集合的相關性,通常介於 -1 (負相關) ~ +1 (正相關) 之間。

交互相關 cross correlation

給定兩個函數 \(x(t), h(t)\),連續時間的 cross correlation 定義為:

\(y(t) = \int_{-∞}^{∞}x^*(τ)h(t+τ)dτ\) 或

\(y(t) = x(t) ⊗ h(t)\)

其中 \(x^*(τ)\) 是共軛複數 complex conjugate 。但數位訊號的輸入資料都是實數,所以 cross correlation 跟卷積運算類似,只是把 \(h(t-τ)\) 換成 \(h(t+τ)\)

數位訊號 x[n], h[n] 的交互相關定義為:

\( y[n] = \sum_{k=-∞}^{∞}x^*[k] \cdot h[n+k] \) 或

\( y[n]=x[n] ⊗ h[n] \)

h[n] 可視為第二個數位訊號,用來進行訊號比對

因輸入訊號的數量 N=7,脈衝響應數量為 M = 5,full cross-correlation 結果長度為 \(M+N-1 = 11\)

計算前先將 \(x[n]\) 兩邊補上 \(M-1=4\) 個 0,稱為 Zero-Padding,但不像卷積運算一樣,不將 h(n) 旋轉 180 度

x h
0 0 0 0 1 2 4 3 2 1 1 0 0 0 0 1 2 3 1 1
1 2 3 1 1
\(y[0] = 0 \cdot 1 + 0 \cdot 2 + 0 \cdot 3 + 0 \cdot 1 + 1 \cdot 1 = 1 \)
x h
0 0 0 0 1 2 4 3 2 1 1 0 0 0 0 1 2 3 1 1
1 2 3 1 1
\(y[1] = 0 \cdot 1 + 0 \cdot 2 + 0 \cdot 3 + 1 \cdot 1 + 2 \cdot 1 = 3 \)

最後卷積運算結果為

n 0 1 2 3 4 5 6 7 8 9 10
y[n] 1 3 9 15 22 22 18 11 7 3 1

實際運算時,會擷取卷積運算的部分結果

n 0 1 2 3 4 5 6
x[n] 1 2 4 3 2 1 1
y[n] 9 15 22 22 18 11 7

結果在 n=2, 3 的位置,數值最大,也就是 x[n] 與 h[n] 在這個位置的相關性最高,波形最相似。

cross correlation 經典的應用是訊號比對 Signal Matching,例如:音樂搜尋系統,聲紋比對系統

numpy 有提供 full, same 兩種 cross correlation 運算的結果

import numpy as np

x = np.array( [ 1, 2, 4, 3, 2, 1, 1 ] )
h = np.array( [ 1, 2, 3, 1, 1 ] )
y = np.correlate( x, h, 'full' )
y1 = np.correlate( x, h, 'same' )

print( "x =", x )
print( "h =", h )
print( "Full Correlation y =", y )
print( "Correlation y =", y1 )

執行結果

$ python correlation.py
x = [1 2 4 3 2 1 1]
h = [1 2 3 1 1]
Full Correlation y = [ 1  3  9 15 22 22 18 11  7  3  1]
Correlation y = [ 9 15 22 22 18 11  7]

如果 h[n] 左右對稱且輸入的數位訊號是實數,則 cross correlation 的結果跟卷積運算一樣。另外 cross correlation 跟卷積一樣,滿足交換率、結合率、分配率。

自相關 autocorrelation

給定兩個函數 \(x(t), h(t)\),連續時間的 autocorrelation 定義為:

\(y(t) = \int_{-∞}^{∞}x^*(τ)x(t+τ)dτ\) 其中 τ 稱為延遲 lag

其中 \(x^*(τ)\) 是共軛複數 complex conjugate ,自相關 autocorrelation 是延遲 τ 的函數。但數位訊號的輸入資料都是實數,所以 autocorrelation 可解釋為:訊號與訊號本身的延遲訊號,進行積分運算的結果。

數位訊號的卷積定義為:

\( R[l] = \sum_{k=-∞}^{∞}x^*[k] \cdot x[n+l] \) 其中 l 是延遲 lag

autocorrelation 是延遲 l 的函數,且 l 為正整數。可用來偵測數位訊號是否具有週期性。通常週期性的數位訊號,自關性較強。 noise 因為每個 sample 在時間軸上都是獨立的,所以不具有 autocorrelation

因輸入訊號的數量 N=5,脈衝響應數量為 M = 5,autocorrelation 結果長度為 5

計算前先將 \(x[n]\) 後面補上 \(N-1=4\) 個 0

x h
1 2 1 2 1 0 0 0 0 1 2 1 2 1
1 2 1 2 1
\(R[0] = 1 \cdot 1 + 2 \cdot 2 + 1 \cdot 1 + 2 \cdot 2 + 1 \cdot 1 = 11 \)
x h
1 2 1 2 1 0 0 0 0 1 2 1 2 1
1 2 1 2 1
\(y[1] = 2 \cdot 1 + 1 \cdot 2 + 2 \cdot 1 + 1 \cdot 2 + 0 \cdot 1 = 8 \)

最後 autocorrelation 運算結果為

n 0 1 2 3 4
R[n] 11 8 6 4 1

cross correlation 運算結果的後半部分,就是 autocorrelation 的結果

import numpy as np

def autocorr( x ):
    # cross correlation 運算結果的後半部分,就是 autocorrelation
    R = np.correlate( x, x, 'full' )
    return R[ int( R.size / 2 ) : ]

def main( ):
    x = np.array( [ 1, 2, 1, 2, 1 ] )
    R = autocorr( x )
    print( "x =", x )
    print( "Autocorrelation =", R )

main( )

執行結果

$ python autocorrelation.py
x = [1 2 1 2 1]
Autocorrelation = [11  8  6  4  1]

autocorrelation 的應用

週期性訊號的自相關性較強,雜訊不具有自相關性。自相關適合用來偵測訊號的週期性,即使週期性的訊號混入了雜訊,一樣能找出自相關性。

import numpy as np
import numpy.random as random
import matplotlib.pyplot as plt

def autocorr( x ):
    R = np.correlate( x, x, 'full' )
    return R[ int( R.size / 2 ) : ]

def main( ):
    t = np.linspace( 0, 1,  100, endpoint = False )     # 定義時間陣列
    x = 10 * np.cos( 2 * np.pi * 5 * t )                # 原始訊號
    noise = random.uniform( -2, 2, 100 )                # 雜訊(均勻分佈)  
    y = x + noise

    auto_corr1 = autocorr( x )
    auto_corr2 = autocorr( noise )
    auto_corr3 = autocorr( y )

    plt.figure( 1 )                                     
    plt.subplot( 121 )
    plt.plot( t, x )
    plt.xlabel( 't (second)' )
    plt.ylabel( 'Amplitude' )

    plt.subplot( 122 )                                  
    plt.plot( auto_corr1 )
    plt.xlabel( 'Lag' )
    plt.ylabel( 'Auto Correlation' )

    plt.figure( 2 )
    plt.subplot( 121 )
    plt.plot( t, noise )
    plt.xlabel( 't (second)' )
    plt.ylabel( 'Amplitude' )

    plt.subplot( 122 )
    plt.plot( auto_corr2 )
    plt.xlabel( 'Lag' )
    plt.ylabel( 'Auto Correlation' )

    plt.figure( 3 ) 
    plt.subplot( 121 )
    plt.plot( t, y )
    plt.xlabel( 't (second)' )
    plt.ylabel( 'Amplitude' )

    plt.subplot( 122 )
    plt.plot( auto_corr3 )
    plt.xlabel( 'Lag' )
    plt.ylabel( 'Auto Correlation' )

    plt.show( )
    
main( )

左圖 x[n] 是原始訊號,右邊是 autocorrelation 的結果

雜訊,及其 autocorrelation

弦波加上雜訊的原始訊號,右圖是 autocorrelation

References

數位訊號處理:Python程式實作(附範例光碟)(第二版)

沒有留言:

張貼留言