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

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

用法:

scipy.signal.zpk2sos(z, p, k, pairing=None, *, analog=False)#

从系统的零点、极点和增益返回二阶部分

参数 ：：

z： array_like

传递函数的零点。

p： array_like

传递函数的极点。

k： 浮点数

系统增益。

pairing： {无，‘nearest’, ‘keep_odd’，‘minimal’}，可选

用于将极点和零点对组合成部分的方法。如果模拟为假且配对为无，则配对设置为‘nearest’；如果模拟为真，则配对必须为‘minimal’，如果为无，则设置为该值。

analog： 布尔型，可选

如果为真，系统是模拟的，否则是离散的。

返回 ：：

sos： ndarray

二阶滤波器系数数组，形状为 (n_sections, 6) 。有关 SOS 过滤器格式规范，请参阅 sosfilt 。

注意：

用于将 ZPK 转换为 SOS 格式的算法旨在最大限度地减少由于数值精度问题导致的错误。配对算法试图最小化每个双二次部分的峰值增益。这是通过将极点与最近的零配对来完成的，对于discrete-time 系统，从最接近单位圆的极点开始，对于continuous-time 系统，从最接近虚轴的极点开始。

pairing='minimal' 输出可能不适合 sosfilt ，并且 analog=True 输出将永远不适合 sosfilt 。

算法

pairing='nearest' 、 pairing='keep_odd' 和 pairing='minimal' 算法中的步骤大多是共享的。 'nearest' 算法试图最小化峰值增益，而'keep_odd' 在 odd-order 系统应保留一个部分作为一阶的约束下最小化峰值增益。 'minimal' 类似于 'keep_odd' ，但没有引入额外的极点或零点

算法步骤如下：

作为 pairing='nearest' 、 pairing='keep_odd' 的预处理步骤，根据需要将极点或零点添加到原点以获得相同数量的极点和零点以进行配对。如果pairing == 'nearest' 并且有奇数个极点，则在原点添加一个附加极点和一个零。

然后迭代以下步骤，直到不再有极点或零点：

1. 取最接近单位圆(或虚轴，对于 analog=True )的(下一个剩余的)极点(复数或实数)开始新的滤波器部分。

2. 如果极点是实数并且没有其他剩余的实数极点 [1]，则将最接近的实数零添加到该部分并将其保留为一阶部分。请注意，在这一步之后，我们保证为后续配对迭代留下偶数个实极点、复数极点、实零点和复数零点。

3. 别的：

   1. If the pole is complex and the zero is the only remaining real zero*, then pair the pole with the next closest zero (guaranteed to be complex). This is necessary to ensure that there will be a real zero remaining to eventually create a first-order section (thus keeping the odd order).

   2. Else pair the pole with the closest remaining zero (complex or real).

   3. Proceed to complete the second-order section by adding another pole and zero to the current pole and zero in the section:

      1. If the current pole and zero are both complex, add their conjugates.

      2. Else if the pole is complex and the zero is real, add the conjugate pole and the next closest real zero.

      3. Else if the pole is real and the zero is complex, add the conjugate zero and the real pole closest to those zeros.

      4. Else (we must have a real pole and real zero) add the next real pole closest to the unit circle, and then add the real zero closest to that pole.

例子：

为采样率为 8000 Hz、pass-band 转角频率为 1000 Hz 的系统设计 6 阶 low-pass 椭圆数字滤波器。 pass-band 中的纹波不应超过 0.087 dB，stop-band 中的衰减应至少为 90 dB。

在下面对 ellip 的调用中，我们可以使用 output='sos' ，但在本示例中，我们将使用 output='zpk' ，然后使用 zpk2sos 转换为 SOS 格式：

>>> from scipy import signal
>>> import numpy as np
>>> z, p, k = signal.ellip(6, 0.087, 90, 1000/(0.5*8000), output='zpk')

现在转换为 SOS 格式。

>>> sos = signal.zpk2sos(z, p, k)

各部分分子的系数：

>>> sos[:, :3]
array([[0.0014152 , 0.00248677, 0.0014152 ],
       [1.        , 0.72976874, 1.        ],
       [1.        , 0.17607852, 1.        ]])

因为所有的零点都在单位圆上，所以会出现系数的对称性。

各部分分母的系数：

>>> sos[:, 3:]
array([[ 1.        , -1.32544025,  0.46989976],
       [ 1.        , -1.26118294,  0.62625924],
       [ 1.        , -1.2570723 ,  0.8619958 ]])

下一个示例显示了配对选项的效果。我们有一个具有三个极点和三个零点的系统，因此 SOS 阵列将具有形状 (2, 6)。实际上，这意味着在 SOS 表示的原点处有一个额外的极点和一个额外的零。

>>> z1 = np.array([-1, -0.5-0.5j, -0.5+0.5j])
>>> p1 = np.array([0.75, 0.8+0.1j, 0.8-0.1j])

使用pairing='nearest'(默认)，我们得到

>>> signal.zpk2sos(z1, p1, 1)
array([[ 1.  ,  1.  ,  0.5 ,  1.  , -0.75,  0.  ],
       [ 1.  ,  1.  ,  0.  ,  1.  , -1.6 ,  0.65]])

第一段有零点 {-0.5-0.05j, -0.5+0.5j} 和极点 {0, 0.75}，第二段有零点 {-1, 0} 和极点 {0.8+0.1j, 0.8 -0.1j}。请注意，原点处的额外极点和零点已分配给不同的部分。

使用 pairing='keep_odd' ，我们获得：

>>> signal.zpk2sos(z1, p1, 1, pairing='keep_odd')
array([[ 1.  ,  1.  ,  0.  ,  1.  , -0.75,  0.  ],
       [ 1.  ,  1.  ,  0.5 ,  1.  , -1.6 ,  0.65]])

原点处的额外极点和零在同一部分。第一部分实际上是first-order 部分。

使用 pairing='minimal' ，first-order 部分在原点没有额外的极点和零点：

>>> signal.zpk2sos(z1, p1, 1, pairing='minimal')
array([[ 0.  ,  1.  ,  1.  ,  0.  ,  1.  , -0.75],
       [ 1.  ,  1.  ,  0.5 ,  1.  , -1.6 ,  0.65]])