scipy.signal.

csd#

scipy.signal.csd(x, y, fs=1.0, window='hann', nperseg=None, noverlap=None, nfft=None, detrend='constant', return_onesided=True, scaling='density', axis=-1, average='mean')[原始碼]#

使用 Welch’s method 估計交叉功率譜密度 Pxy。

參數：

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 是 str 或 tuple，則設定為 256；如果 window 是 array_like，則設定為視窗的長度。
noverlap：int，選用: 區段之間重疊的點數。如果 None，noverlap = nperseg // 2。預設值為 None。
nfftint，選用: 使用的 FFT 長度，如果需要零填充 FFT。如果 None，則 FFT 長度為 nperseg。預設值為 None。
detrendstr 或 function 或 False，選用: 指定如何對每個區段進行去趨勢化。如果 detrend 是字串，則會將其作為 type 引數傳遞至 detrend 函數。如果它是函數，則會接受一個區段並傳回一個去趨勢化的區段。如果 detrend 是 False，則不進行去趨勢化。預設值為 ‘constant’。
return_onesidedbool，選用: 如果 True，則傳回實數資料的單邊頻譜。如果 False，則傳回雙邊頻譜。預設值為 True，但對於複數資料，始終傳回雙邊頻譜。
scaling{ ‘density’, ‘spectrum’ }，選用: 選擇計算交叉頻譜密度 (‘density’)，其中 Pxy 的單位為 V**2/Hz，以及計算交叉頻譜 (‘spectrum’)，其中 Pxy 的單位為 V**2，如果 x 和 y 以 V 為單位測量，而 fs 以 Hz 為單位測量。預設值為 ‘density’
axisint，選用: 計算兩個輸入的 CSD 所沿著的軸；預設值是最後一個軸（即 axis=-1）。
average{ ‘mean’, ‘median’ }，選用: 平均週期圖時使用的方法。如果頻譜是複數，則實部和虛部的平均值會分開計算。預設值為 ‘mean’。

在版本 1.2.0 中新增。

傳回：

fndarray: 樣本頻率的陣列。
Pxyndarray: x,y 的交叉頻譜密度或交叉功率譜。

另請參閱

periodogram: 簡單、可選修改的週期圖
lombscargle: 用於不均勻取樣資料的 Lomb-Scargle 週期圖
welch: Welch’s 方法的功率譜密度。[等同於 csd(x,x)]
coherence: Welch’s 方法的平方同調性。

註解

依照慣例，Pxy 是使用 X 的共軛 FFT 乘以 Y 的 FFT 來計算。

如果輸入序列的長度不同，則較短的序列將會進行零填充以匹配。

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

請參閱 SciPy 使用者指南的頻譜分析章節，以了解頻譜密度和（振幅）頻譜的縮放比例。

在版本 0.16.0 中新增。

參考文獻

[1]

P. Welch，“The use of the fast Fourier transform for the estimation of power spectra: A method based on time averaging over short, modified periodograms”，IEEE Trans. Audio Electroacoust. vol. 15, pp. 70-73, 1967。

[2]

Rabiner, Lawrence R., and B. Gold。“Theory and Application of Digital Signal Processing” Prentice-Hall, pp. 414-419, 1975

範例

>>> 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, Pxy = signal.csd(x, y, fs, nperseg=1024)
>>> plt.semilogy(f, np.abs(Pxy))
>>> plt.xlabel('frequency [Hz]')
>>> plt.ylabel('CSD [V**2/Hz]')
>>> plt.show()