Python SciPy signal.sosfreqz用法及代码示例

本文简要介绍 python 语言中 scipy.signal.sosfreqz 的用法。

用法:

scipy.signal.sosfreqz(sos, worN=512, whole=False, fs=6.283185307179586)#

以 SOS 格式计算数字滤波器的频率响应。

给定 sos，一个形状为 (n, 6) 的数字滤波器二阶部分的数组，计算系统函数的频率响应：

B0(z)   B1(z)         B{n-1}(z)
H(z) = ----- * ----- * ... * ---------
       A0(z)   A1(z)         A{n-1}(z)

对于 z = exp(omega*1j)，其中 B{k}(z) 和 A{k}(z) 是 k-th 二阶部分的传递函数的分子和分母。

参数 ：：

sos： array_like

二阶滤波器系数数组，必须具有形状 (n_sections, 6) 。每行对应一个二阶部分，前三列提供分子系数，后三列提供分母系数。

worN： {无，int，数组}，可选

如果是单个整数，则以那么多频率进行计算(默认为 N=512)。对 FFT 计算使用快速的数字可以导致更快的计算(请参阅 freqz 的注释)。

如果是 数组，计算给定频率的响应(必须是 1-D)。这些与 fs 的单位相同。

whole： 布尔型，可选

通常，频率的计算范围是从 0 到奈奎斯特频率 fs/2(upper-half of unit-circle)。如果整体为真，则计算从 0 到 fs 的频率。

fs： 浮点数，可选

数字系统的采样频率。默认为 2*pi 弧度/样本(所以 w 是从 0 到 pi)。

返回 ：：

w： ndarray

计算 h 的频率，单位与 fs 相同。默认情况下，w 被归一化为范围 [0, pi)(弧度/样本)。

h： ndarray

频率响应，作为复数。

注意：

例子：

设计一个 SOS 格式的 15th-order 带通滤波器。

>>> from scipy import signal
>>> import numpy as np
>>> sos = signal.ellip(15, 0.5, 60, (0.2, 0.4), btype='bandpass',
...                    output='sos')

计算从 DC 到 Nyquist 的 1500 个点的频率响应。

>>> w, h = signal.sosfreqz(sos, worN=1500)

绘制响应。

>>> import matplotlib.pyplot as plt
>>> plt.subplot(2, 1, 1)
>>> db = 20*np.log10(np.maximum(np.abs(h), 1e-5))
>>> plt.plot(w/np.pi, db)
>>> plt.ylim(-75, 5)
>>> plt.grid(True)
>>> plt.yticks([0, -20, -40, -60])
>>> plt.ylabel('Gain [dB]')
>>> plt.title('Frequency Response')
>>> plt.subplot(2, 1, 2)
>>> plt.plot(w/np.pi, np.angle(h))
>>> plt.grid(True)
>>> plt.yticks([-np.pi, -0.5*np.pi, 0, 0.5*np.pi, np.pi],
...            [r'$-\pi$', r'$-\pi/2$', '0', r'$\pi/2$', r'$\pi$'])
>>> plt.ylabel('Phase [rad]')
>>> plt.xlabel('Normalized frequency (1.0 = Nyquist)')
>>> plt.show()

如果将相同的滤波器实现为单个传递函数，则数值误差会破坏频率响应：

>>> b, a = signal.ellip(15, 0.5, 60, (0.2, 0.4), btype='bandpass',
...                    output='ba')
>>> w, h = signal.freqz(b, a, worN=1500)
>>> plt.subplot(2, 1, 1)
>>> db = 20*np.log10(np.maximum(np.abs(h), 1e-5))
>>> plt.plot(w/np.pi, db)
>>> plt.ylim(-75, 5)
>>> plt.grid(True)
>>> plt.yticks([0, -20, -40, -60])
>>> plt.ylabel('Gain [dB]')
>>> plt.title('Frequency Response')
>>> plt.subplot(2, 1, 2)
>>> plt.plot(w/np.pi, np.angle(h))
>>> plt.grid(True)
>>> plt.yticks([-np.pi, -0.5*np.pi, 0, 0.5*np.pi, np.pi],
...            [r'$-\pi$', r'$-\pi/2$', '0', r'$\pi/2$', r'$\pi$'])
>>> plt.ylabel('Phase [rad]')
>>> plt.xlabel('Normalized frequency (1.0 = Nyquist)')
>>> plt.show()