Python SciPy signal.kaiserord用法及代碼示例

本文簡要介紹 python 語言中 scipy.signal.kaiserord 的用法。

用法:

scipy.signal.kaiserord(ripple, width)#

確定 Kaiser 窗口方法的過濾器窗口參數。

此函數返回的參數通常用於使用窗口方法創建有限脈衝響應濾波器，使用 firwin 或 firwin2 。

參數 ：：

ripple： 浮點數

濾波器頻率響應幅度與所需濾波器幅度的偏差(以 dB 為單位)的上限(不包括任何過渡間隔中的頻率)。也就是說，如果 w 是表示為奈奎斯特頻率分數的頻率，A(w) 是濾波器的實際頻率響應，D(w) 是所需的頻率響應，則設計要求是：

abs(A(w) - D(w))) < 10**(-ripple/20)

對於 0 <= w <= 1 並且 w 不在轉換間隔中。

width： 浮點數

過渡區域的寬度，標準化以便 1 對應於 pi 弧度/樣本。也就是說，頻率表示為奈奎斯特頻率的一部分。

返回 ：：

numtaps： int

Kaiser 窗口的長度。

beta： 浮點數

Kaiser 窗口的 beta 參數。

注意：

有幾種方法可以獲取 Kaiser 窗口：

• signal.windows.kaiser(numtaps, beta, sym=True)

• signal.get_window(beta, numtaps)

• signal.get_window(('kaiser', beta), numtaps)

使用由 Kaiser 發現的經驗方程。

參考：

Oppenheim, Schafer，“Discrete-Time 信號處理”，第 475-476 頁。

例子：

我們將使用 Kaiser 窗方法為以 1000 Hz 采樣的信號設計低通 FIR 濾波器。

我們希望在阻帶中至少有 65 dB 的抑製，而在通帶中，增益的變化應該不超過 0.5%。

我們想要一個 175 Hz 的截止頻率，以及 24 Hz 的通帶和阻帶之間的過渡。也就是說，在[0, 163]頻段，增益變化不超過0.5%，在[187, 500]頻段，信號衰減至少65 dB。

>>> import numpy as np
>>> from scipy.signal import kaiserord, firwin, freqz
>>> import matplotlib.pyplot as plt
>>> fs = 1000.0
>>> cutoff = 175
>>> width = 24

Kaiser 方法隻接受一個參數來控製通帶紋波和阻帶抑製，因此我們使用兩者中限製性更強的一個。在這種情況下，通帶紋波為 0.005，即 46.02 dB，因此我們將使用 65 dB 作為設計參數。

使用 kaiserord 確定濾波器的長度和 Kaiser 窗口的參數。

>>> numtaps, beta = kaiserord(65, width/(0.5*fs))
>>> numtaps
167
>>> beta
6.20426

使用 firwin 創建 FIR 濾波器。

>>> taps = firwin(numtaps, cutoff, window=('kaiser', beta),
...               scale=False, fs=fs)

計算濾波器的頻率響應。 w 是頻率數組，h 是對應的頻率響應複數數組。

>>> w, h = freqz(taps, worN=8000)
>>> w *= 0.5*fs/np.pi  # Convert w to Hz.

計算濾波器響應幅度與理想低通濾波器幅度的偏差。過渡區域中的值設置為 nan ，因此它們不會出現在圖中。

>>> ideal = w < cutoff  # The "ideal" frequency response.
>>> deviation = np.abs(np.abs(h) - ideal)
>>> deviation[(w > cutoff - 0.5*width) & (w < cutoff + 0.5*width)] = np.nan

繪製偏差。仔細觀察阻帶的左端會發現，在第一個波瓣中，65 dB 衰減的要求被違反了約 0.125 dB。這對於 Kaiser 窗口方法來說並不少見。

>>> plt.plot(w, 20*np.log10(np.abs(deviation)))
>>> plt.xlim(0, 0.5*fs)
>>> plt.ylim(-90, -60)
>>> plt.grid(alpha=0.25)
>>> plt.axhline(-65, color='r', ls='--', alpha=0.3)
>>> plt.xlabel('Frequency (Hz)')
>>> plt.ylabel('Deviation from ideal (dB)')
>>> plt.title('Lowpass Filter Frequency Response')
>>> plt.show()