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

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

用法:

scipy.signal.lfilter(b, a, x, axis=-1, zi=None)#

使用 IIR 或 FIR 濾波器沿 one-dimension 過濾數據。

使用數字濾波器過濾數據序列 x。這適用於許多基本數據類型(包括對象類型)。該濾波器是標準差分方程的直接形式 II 轉置實現(見注釋)。

對於大多數過濾任務，函數 sosfilt (以及使用 output='sos' 的濾波器設計)應優先於 lfilter，因為二階部分的數值問題較少。

參數 ：：

b： array_like

一維序列中的分子係數向量。

a： array_like

一維序列中的分母係數向量。如果a[0]不是 1，那麽兩者都是a和b歸一化a[0].

x： array_like

一個 N 維輸入數組。

axis： 整數，可選

沿其應用線性濾波器的輸入數據數組的軸。過濾器沿該軸應用於每個子陣列。默認值為 -1。

zi： 數組，可選

濾波器延遲的初始條件。它是一個長度為的向量(或 N 維輸入的向量數組)max(len(a), len(b)) - 1.如果子為 None 或未給出，則假定初始休息。看lfiltic了解更多信息。

返回 ：：

y： 數組

數字濾波器的輸出。

zf： 數組，可選

如果 zi 為 None，則不返回，否則，zf 保存最終的濾波器延遲值。

注意：

濾波器函數實現為直接 II 轉置結構。這意味著過濾器實現：

a[0]*y[n] = b[0]*x[n] + b[1]*x[n-1] + ... + b[M]*x[n-M]
                      - a[1]*y[n-1] - ... - a[N]*y[n-N]

其中M是分子的次數，N是分母的次數，n是樣本數。它使用以下差分方程實現(假設 M = N)：

a[0]*y[n] = b[0] * x[n]               + d[0][n-1]
  d[0][n] = b[1] * x[n] - a[1] * y[n] + d[1][n-1]
  d[1][n] = b[2] * x[n] - a[2] * y[n] + d[2][n-1]
...
d[N-2][n] = b[N-1]*x[n] - a[N-1]*y[n] + d[N-1][n-1]
d[N-1][n] = b[N] * x[n] - a[N] * y[n]

其中 d 是狀態變量。

z-transform 域中說明此濾波器的有理傳遞函數為：

-1              -M
        b[0] + b[1]z  + ... + b[M] z
Y(z) = -------------------------------- X(z)
                    -1              -N
        a[0] + a[1]z  + ... + a[N] z

例子：

生成要過濾的噪聲信號：

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt
>>> rng = np.random.default_rng()
>>> t = np.linspace(-1, 1, 201)
>>> x = (np.sin(2*np.pi*0.75*t*(1-t) + 2.1) +
...      0.1*np.sin(2*np.pi*1.25*t + 1) +
...      0.18*np.cos(2*np.pi*3.85*t))
>>> xn = x + rng.standard_normal(len(t)) * 0.08

創建一個 3 階低通巴特沃斯濾波器：

>>> b, a = signal.butter(3, 0.05)

將過濾器應用於 xn。使用lfilter_zi選擇過濾器的初始條件：

>>> zi = signal.lfilter_zi(b, a)
>>> z, _ = signal.lfilter(b, a, xn, zi=zi*xn[0])

再次應用過濾器，以與 filtfilt 相同的順序過濾結果：

>>> z2, _ = signal.lfilter(b, a, z, zi=zi*z[0])

使用 filtfilt 應用過濾器：

>>> y = signal.filtfilt(b, a, xn)

繪製原始信號和各種濾波版本：

>>> plt.figure
>>> plt.plot(t, xn, 'b', alpha=0.75)
>>> plt.plot(t, z, 'r--', t, z2, 'r', t, y, 'k')
>>> plt.legend(('noisy signal', 'lfilter, once', 'lfilter, twice',
...             'filtfilt'), loc='best')
>>> plt.grid(True)
>>> plt.show()