scipy.signal.

同調性#

scipy.signal.coherence(x, y, fs=1.0, window='hann', nperseg=None, noverlap=None, nfft=None, detrend='constant', axis=-1)[source]#

使用 Welch 方法估計離散時間訊號 X 和 Y 的強度平方同調性估計值 Cxy。

Cxy = abs(Pxy)**2/(Pxx*Pyy)，其中 Pxx 和 Pyy 是 X 和 Y 的功率頻譜密度估計值，而 Pxy 是 X 和 Y 的交叉頻譜密度估計值。

參數:

xarray_like: 量測值的時間序列
yarray_like: 量測值的時間序列
fsfloat，選填: x 和 y 時間序列的取樣頻率。預設值為 1.0。
windowstr 或 tuple 或 array_like，選填: 想要使用的窗函數。如果 window 是字串或 tuple，它會傳遞給 get_window 以產生窗函數值，預設為 DFT-even。請參閱 get_window 以取得窗函數列表和所需參數。如果 window 是 array_like，它將直接用作窗函數，且其長度必須為 nperseg。預設為 Hann 窗函數。
npersegint，選填: 每個區段的長度。預設值為 None，但如果 window 是字串或 tuple，則設定為 256；如果 window 是 array_like，則設定為窗函數的長度。
noverlapint，選填: 區段之間重疊的點數。如果為 None，則 noverlap = nperseg // 2。預設值為 None。
nfftint，選填: FFT 的長度，如果需要零填充 FFT。如果為 None，則 FFT 長度為 nperseg。預設值為 None。
detrendstr 或 function 或 False，選填: 指定如何對每個區段進行去趨勢。如果 detrend 是字串，它會作為 type 參數傳遞給 detrend 函數。如果是函數，它會接收一個區段並傳回去趨勢後的區段。如果 detrend 是 False，則不進行去趨勢。預設值為 ‘constant’。
axisint，選填: 計算兩個輸入同調性的軸；預設值為最後一個軸 (即 axis=-1)。

回傳:

fndarray: 取樣頻率的陣列。
Cxyndarray: x 和 y 的強度平方同調性。

另請參閱

periodogram: 簡單、可選修改的週期圖
lombscargle: 用於不均勻取樣資料的 Lomb-Scargle 週期圖
welch: 使用 Welch 方法計算功率頻譜密度。
csd: 使用 Welch 方法計算交叉頻譜密度。

註解

適當的重疊量將取決於窗函數的選擇和您的需求。對於預設的 Hann 窗函數，50% 的重疊是在準確估計訊號功率，同時不過度計算任何資料之間合理的權衡。較窄的窗函數可能需要更大的重疊。

在版本 0.16.0 中新增。

參考文獻

[1]

P. Welch，“使用快速傅立葉轉換估計功率譜：一種基於對短時修改週期圖進行時間平均的方法”，IEEE Trans. Audio Electroacoust. vol. 15, pp. 70-73, 1967。

[2]

Stoica, Petre 和 Randolph Moses，“訊號頻譜分析” Prentice Hall，2005

範例

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt
>>> rng = np.random.default_rng()

產生兩個具有一些共同特徵的測試訊號。

>>> fs = 10e3
>>> N = 1e5
>>> amp = 20
>>> freq = 1234.0
>>> noise_power = 0.001 * fs / 2
>>> time = np.arange(N) / fs
>>> b, a = signal.butter(2, 0.25, 'low')
>>> x = rng.normal(scale=np.sqrt(noise_power), size=time.shape)
>>> y = signal.lfilter(b, a, x)
>>> x += amp*np.sin(2*np.pi*freq*time)
>>> y += rng.normal(scale=0.1*np.sqrt(noise_power), size=time.shape)

計算並繪製同調性。

>>> f, Cxy = signal.coherence(x, y, fs, nperseg=1024)
>>> plt.semilogy(f, Cxy)
>>> plt.xlabel('frequency [Hz]')
>>> plt.ylabel('Coherence')
>>> plt.show()

../../_images/scipy-signal-coherence-1.png