scipy.signal.

freqz#

scipy.signal.freqz(b, a=1, worN=512, whole=False, plot=None, fs=6.283185307179586, include_nyquist=False)[原始碼]#

計算數位濾波器的頻率響應。

給定數位濾波器的 M 階分子 b 和 N 階分母 a,計算其頻率響應

            jw                 -jw              -jwM
   jw    B(e  )    b[0] + b[1]e    + ... + b[M]e
H(e  ) = ------ = -----------------------------------
            jw                 -jw              -jwN
         A(e  )    a[0] + a[1]e    + ... + a[N]e
參數:
barray_like

線性濾波器的分子。如果 b 的維度大於 1,則假定係數儲存在第一個維度中,且 b.shape[1:]a.shape[1:] 和頻率陣列的形狀必須與廣播相容。

aarray_like

線性濾波器的分母。如果 b 的維度大於 1,則假定係數儲存在第一個維度中,且 b.shape[1:]a.shape[1:] 和頻率陣列的形狀必須與廣播相容。

worN{None, int, array_like}, optional

如果是一個整數,則計算該數量的頻率(預設值為 N=512)。這是以下項目的便捷替代方案

np.linspace(0, fs if whole else fs/2, N, endpoint=include_nyquist)

使用對 FFT 計算快速的數字可以加快計算速度(請參閱「註解」)。

如果是一個 array_like,則計算給定頻率下的響應。這些頻率的單位與 fs 相同。

wholebool, optional

通常,頻率是從 0 計算到奈奎斯特頻率 fs/2(單位圓的上半部)。如果 whole 為 True,則計算從 0 到 fs 的頻率。如果 worN 是 array_like,則忽略此參數。

plotcallable

一個接受兩個引數的可呼叫物件。如果給定,則傳回參數 wh 會傳遞給 plot。適用於在 freqz 內部繪製頻率響應。

fsfloat, optional

數位系統的取樣頻率。預設值為 2*pi 弧度/樣本(因此 w 的範圍為 0 到 pi)。

在 1.2.0 版本中新增。

include_nyquistbool, optional

如果 whole 為 False 且 worN 是一個整數,則將 include_nyquist 設定為 True 將包含最後一個頻率(奈奎斯特頻率),否則會被忽略。

在 1.5.0 版本中新增。

傳回值:
wndarray

h 在其中計算的頻率,單位與 fs 相同。預設情況下,w 會正規化到範圍 [0, pi)(弧度/樣本)。

hndarray

頻率響應,以複數表示。

另請參閱

freqz_zpk
freqz_sos

註解

將 Matplotlib 的 matplotlib.pyplot.plot 函數用作 plot 的可呼叫物件會產生非預期的結果,因為這會繪製複數轉移函數的實部,而不是幅度。請嘗試 lambda w, h: plot(w, np.abs(h))

當符合以下條件時,會使用透過 (R)FFT 的直接計算來計算頻率響應

  1. worN 提供了一個整數值。

  2. worN 可以透過 FFT 快速計算(即,next_fast_len(worN) 等於 worN)。

  3. 分母係數是單個值 (a.shape[0] == 1)。

  4. worN 至少與分子係數一樣長 (worN >= b.shape[0])。

  5. 如果 b.ndim > 1,則 b.shape[-1] == 1

對於長 FIR 濾波器,FFT 方法的誤差可能較低,並且比等效的直接多項式計算快得多。

範例

>>> from scipy import signal
>>> import numpy as np
>>> taps, f_c = 80, 1.0  # number of taps and cut-off frequency
>>> b = signal.firwin(taps, f_c, window=('kaiser', 8), fs=2*np.pi)
>>> w, h = signal.freqz(b)
>>> import matplotlib.pyplot as plt
>>> fig, ax1 = plt.subplots(tight_layout=True)
>>> ax1.set_title(f"Frequency Response of {taps} tap FIR Filter" +
...               f"($f_c={f_c}$ rad/sample)")
>>> ax1.axvline(f_c, color='black', linestyle=':', linewidth=0.8)
>>> ax1.plot(w, 20 * np.log10(abs(h)), 'C0')
>>> ax1.set_ylabel("Amplitude in dB", color='C0')
>>> ax1.set(xlabel="Frequency in rad/sample", xlim=(0, np.pi))
>>> ax2 = ax1.twinx()
>>> phase = np.unwrap(np.angle(h))
>>> ax2.plot(w, phase, 'C1')
>>> ax2.set_ylabel('Phase [rad]', color='C1')
>>> ax2.grid(True)
>>> ax2.axis('tight')
>>> plt.show()
../../_images/scipy-signal-freqz-1_00_00.png

廣播範例

假設我們有兩個 FIR 濾波器,其係數儲存在形狀為 (2, 25) 的陣列的列中。為了演示,我們將使用隨機資料

>>> rng = np.random.default_rng()
>>> b = rng.random((2, 25))

為了透過一次呼叫 freqz 來計算這兩個濾波器的頻率響應,我們必須傳入 b.T,因為 freqz 期望第一個軸保存係數。然後我們必須使用長度為 1 的 trivial 維度來擴展形狀,以允許與頻率陣列進行廣播。也就是說,我們傳入形狀為 (25, 2, 1) 的 b.T[..., np.newaxis]

>>> w, h = signal.freqz(b.T[..., np.newaxis], worN=1024)
>>> w.shape
(1024,)
>>> h.shape
(2, 1024)

現在,假設我們有兩個轉移函數,它們具有相同的分子係數 b = [0.5, 0.5]。兩個分母的係數儲存在 2-D 陣列 a 的第一個維度中

a = [   1      1  ]
    [ -0.25, -0.5 ]
>>> b = np.array([0.5, 0.5])
>>> a = np.array([[1, 1], [-0.25, -0.5]])

只有 a 是超過 1-D 的。為了使其與頻率廣播相容,我們在呼叫 freqz 時使用 trivial 維度擴展它

>>> w, h = signal.freqz(b, a[..., np.newaxis], worN=1024)
>>> w.shape
(1024,)
>>> h.shape
(2, 1024)