SciPy-1-12-中文文档-十二-

SciPy 1.12 中文文档（十二）

原文：docs.scipy.org/doc/scipy-1.12.0/index.html

scipy.signal.cubic

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.cubic.html#scipy.signal.cubic

scipy.signal.cubic(x)

自 SciPy 1.11.0 版本起不推荐使用：scipy.signal.cubic在 SciPy 1.11 中已弃用，并将在 SciPy 1.13 中移除。

具体的等价（对于浮点数组x）是：

>>> from scipy.interpolate import BSpline
>>> out = BSpline.basis_element([-2, -1, 0, 1, 2])(x)
>>> out[(x < -2 | (x > 2)] = 0.0 

三次 B 样条。

这是bspline的一个特例，相当于bspline(x, 3)。

参数：

xarray_like

一个结点向量

返回：

resndarray

三次 B 样条基函数值

另请参阅

bspline

B 样条阶数为 n 的基函数

quadratic

二次 B 样条。

示例

我们可以计算几个阶数的 B 样条基函数：

>>> import numpy as np
>>> from scipy.signal import bspline, cubic, quadratic
>>> bspline(0.0, 1)
1 

>>> knots = [-1.0, 0.0, -1.0]
>>> bspline(knots, 2)
array([0.125, 0.75, 0.125]) 

>>> np.array_equal(bspline(knots, 2), quadratic(knots))
True 

>>> np.array_equal(bspline(knots, 3), cubic(knots))
True 

scipy.signal.quadratic

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.quadratic.html#scipy.signal.quadratic

scipy.signal.quadratic(x)

自 1.11.0 版本起弃用：scipy.signal.quadratic

对于浮点数数组x的确切等效是：

>>> from scipy.interpolate import BSpline
>>> out = BSpline.basis_element([-1.5, -0.5, 0.5, 1.5])(x)
>>> out[(x < -1.5 | (x > 1.5)] = 0.0 

一个二次 B 样条。

这是bspline的特殊情况，等同于bspline(x, 2)。

参数：

xarray_like

结点向量

返回：

resndarray

二次 B 样条基函数值

另请参阅

bspline

阶数为 n 的 B 样条基函数

cubic

一个三次 B 样条。

示例

我们可以计算多个阶次的 B 样条基函数：

>>> import numpy as np
>>> from scipy.signal import bspline, cubic, quadratic
>>> bspline(0.0, 1)
1 

>>> knots = [-1.0, 0.0, -1.0]
>>> bspline(knots, 2)
array([0.125, 0.75, 0.125]) 

>>> np.array_equal(bspline(knots, 2), quadratic(knots))
True 

>>> np.array_equal(bspline(knots, 3), cubic(knots))
True 

scipy.signal.gauss_spline

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.gauss_spline.html#scipy.signal.gauss_spline

scipy.signal.gauss_spline(x, n)

n 阶 B 样条基函数的高斯近似。

参数：

x array_like

结节向量

n整型

样条的阶数。必须为非负数，即 n >= 0

返回值：

res ndarray

B 样条基函数值由均值为零的高斯函数近似。

注意事项：

B 样条基函数可以用均值为零、标准差等于(\sigma=(n+1)/12)的高斯函数很好地近似：

[\frac{1}{\sqrt {2\pi\sigma²}}exp(-\frac{x²}{2\sigma})]

参考文献：

[1]

Bouma H., Vilanova A., Bescos J.O., ter Haar Romeny B.M., Gerritsen F.A. (2007) 基于 B 样条的快速精确高斯导数。在：Sgallari F., Murli A., Paragios N. (eds) 计算机视觉中的尺度空间与变分方法。SSVM 2007. 计算机科学讲座笔记，4485. Springer, Berlin, Heidelberg

[2]

folk.uio.no/inf3330/scripting/doc/python/SciPy/tutorial/old/node24.html

示例

我们可以计算由高斯分布近似的 B 样条基函数：

>>> import numpy as np
>>> from scipy.signal import gauss_spline, bspline
>>> knots = np.array([-1.0, 0.0, -1.0])
>>> gauss_spline(knots, 3)
array([0.15418033, 0.6909883, 0.15418033])  # may vary 

>>> bspline(knots, 3)
array([0.16666667, 0.66666667, 0.16666667])  # may vary 

scipy.signal.cspline1d

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.cspline1d.html#scipy.signal.cspline1d

scipy.signal.cspline1d(signal, lamb=0.0)

计算秩-1 数组的三次样条系数。

假设镜像对称边界条件，找到一维信号的三次样条系数。为了从样条表示中恢复信号，使用长度为 3 的 FIR 窗口 [1.0, 4.0, 1.0]/ 6.0 镜像对称卷积这些系数。

参数：

signal ndarray

一个表示信号样本的秩-1 数组。

lamb float, optional

平滑系数，默认为 0.0。

返回值：

c ndarray

三次样条系数。

参见

cspline1d_eval

在新点集上评估三次样条。

示例

我们可以使用三次样条来滤波信号，以减少并平滑高频噪声：

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from scipy.signal import cspline1d, cspline1d_eval
>>> rng = np.random.default_rng()
>>> sig = np.repeat([0., 1., 0.], 100)
>>> sig += rng.standard_normal(len(sig))*0.05  # add noise
>>> time = np.linspace(0, len(sig))
>>> filtered = cspline1d_eval(cspline1d(sig), time)
>>> plt.plot(sig, label="signal")
>>> plt.plot(time, filtered, label="filtered")
>>> plt.legend()
>>> plt.show() 

scipy.signal.qspline1d

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.qspline1d.html#scipy.signal.qspline1d

scipy.signal.qspline1d(signal, lamb=0.0)

计算秩为 1 的数组的二次样条系数。

参数：

signalndarray

代表信号样本的秩为 1 的数组。

lambfloat, optional

平滑系数（现在必须为零）。

返回：

cndarray

二次样条系数。

另请参见

qspline1d_eval

在新的点集上评估二次样条。

注意事项

假设镜像对称边界条件，为 1-D 信号找到二次样条系数。为了从样条表示中恢复信号，使用长度为 3 的 FIR 窗口[1.0, 6.0, 1.0] / 8.0 镜像对称卷积这些系数。

示例

通过二次样条可以滤波信号，以减少和平滑高频噪声：

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from scipy.signal import qspline1d, qspline1d_eval
>>> rng = np.random.default_rng()
>>> sig = np.repeat([0., 1., 0.], 100)
>>> sig += rng.standard_normal(len(sig))*0.05  # add noise
>>> time = np.linspace(0, len(sig))
>>> filtered = qspline1d_eval(qspline1d(sig), time)
>>> plt.plot(sig, label="signal")
>>> plt.plot(time, filtered, label="filtered")
>>> plt.legend()
>>> plt.show() 

scipy.signal.cspline2d

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.cspline2d.html#scipy.signal.cspline2d

scipy.signal.cspline2d(input, lambda=0.0, precision=-1.0)

二维立方（3 阶）B 样条的系数。

返回二维输入图像的规则间隔输入网格上的三阶 B 样条系数。

参数：

输入ndarray

输入信号。

λfloat

指定传递函数中平滑的程度。

精度float

指定计算无限和以应用镜像对称边界条件所需的精度。

返回：

输出ndarray

过滤后的信号。

示例

示例可以在 tutorial 中找到。

scipy.signal.qspline2d

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.qspline2d.html#scipy.signal.qspline2d

scipy.signal.qspline2d(input, lambda=0.0, precision=-1.0)

2-D 二次（2 阶）B 样条的系数：

返回二维输入图像在规则间隔的输入网格上的二阶 B 样条系数。

参数：

inputndarray

输入信号。

lambdafloat

指定传递函数中的平滑量。

precisionfloat

指定计算应用镜像对称边界条件所需的无限和的精度。

返回：

outputndarray

过滤后的信号。

示例

示例参见教程。

scipy.signal.cspline1d_eval

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.cspline1d_eval.html#scipy.signal.cspline1d_eval

scipy.signal.cspline1d_eval(cj, newx, dx=1.0, x0=0)

在新一组点上评估三次样条曲线。

dx 是旧的采样间距，而 x0 是旧的原点。换句话说，cj 表示样条系数的旧样本点（结点）是等间距点：

oldx = x0 + j*dx j=0…N-1，其中 N=len(cj)

边界使用镜像对称边界条件处理。

参数：

cjndarray

三次样条曲线系数

newxndarray

新一组点。

dxfloat，可选

旧的采样间距，默认值为 1.0。

x0int，可选

旧的原点，默认值为 0。

返回：

resndarray

评估了三次样条曲线点。

另请参见

cspline1d

计算一维数组的三次样条系数。

示例

我们可以使用三次样条滤波器来过滤信号，以减少和平滑高频噪声：

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from scipy.signal import cspline1d, cspline1d_eval
>>> rng = np.random.default_rng()
>>> sig = np.repeat([0., 1., 0.], 100)
>>> sig += rng.standard_normal(len(sig))*0.05  # add noise
>>> time = np.linspace(0, len(sig))
>>> filtered = cspline1d_eval(cspline1d(sig), time)
>>> plt.plot(sig, label="signal")
>>> plt.plot(time, filtered, label="filtered")
>>> plt.legend()
>>> plt.show() 

scipy.signal.qspline1d_eval

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.qspline1d_eval.html#scipy.signal.qspline1d_eval

scipy.signal.qspline1d_eval(cj, newx, dx=1.0, x0=0)

在新的一组点上评估二次样条。

参数：

cjndarray

二次样条系数

newxndarray

新的一组点。

dxfloat，可选

旧的样本间距，默认值为 1.0。

x0int，可选

旧原点，默认值为 0。

返回：

resndarray

评估了二次样条点。

参见

qspline1d

为 rank-1 数组计算二次样条系数。

注意

dx是旧的样本间距，而x0是旧的原点。换句话说，旧样本点（结点点）为均匀间隔的点：

oldx = x0 + j*dx  j=0...N-1, with N=len(cj) 

边界使用镜像对称边界条件处理。

示例

我们可以使用二次样条对信号进行滤波，以减少和平滑高频噪声：

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from scipy.signal import qspline1d, qspline1d_eval
>>> rng = np.random.default_rng()
>>> sig = np.repeat([0., 1., 0.], 100)
>>> sig += rng.standard_normal(len(sig))*0.05  # add noise
>>> time = np.linspace(0, len(sig))
>>> filtered = qspline1d_eval(qspline1d(sig), time)
>>> plt.plot(sig, label="signal")
>>> plt.plot(time, filtered, label="filtered")
>>> plt.legend()
>>> plt.show() 

scipy.signal.spline_filter

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.spline_filter.html#scipy.signal.spline_filter

scipy.signal.spline_filter(Iin, lmbda=5.0)

对一个排名为 2 的数组进行三次样条（立方）平滑滤波。

使用（立方）平滑样条滤波器对输入数据集Iin进行滤波。

参数：

Iinarray_like

输入数据集

lmbdafloat，可选

样条平滑的衰减值，默认为5.0。

返回：

resndarray

过滤后的输入数据

示例

我们可以使用三次 B 样条滤波器来过滤多维信号（例如：2D 图像）：

>>> import numpy as np
>>> from scipy.signal import spline_filter
>>> import matplotlib.pyplot as plt
>>> orig_img = np.eye(20)  # create an image
>>> orig_img[10, :] = 1.0
>>> sp_filter = spline_filter(orig_img, lmbda=0.1)
>>> f, ax = plt.subplots(1, 2, sharex=True)
>>> for ind, data in enumerate([[orig_img, "original image"],
...                             [sp_filter, "spline filter"]]):
...     ax[ind].imshow(data[0], cmap='gray_r')
...     ax[ind].set_title(data[1])
>>> plt.tight_layout()
>>> plt.show() 

scipy.signal.order_filter

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.order_filter.html#scipy.signal.order_filter

scipy.signal.order_filter(a, domain, rank)

在一个 N 维数组上执行顺序滤波。

在输入数组上执行顺序滤波。domain 参数充当以每个像素为中心的蒙版。domain 的非零元素用于选择围绕每个输入像素的元素，并放置在一个列表中。列表被排序，该像素的输出是在排序列表中对应于 rank 的元素。

参数：

andarray

N 维输入数组。

domainarray_like

与a具有相同维数的蒙版数组。每个维度应该有奇数个元素。

rankint

一个非负整数，用于从排序列表中选择元素（0 对应最小元素，1 是下一个最小元素，依此类推）。

返回：

outndarray

与a相同形状的数组中的有序滤波结果。

示例

>>> import numpy as np
>>> from scipy import signal
>>> x = np.arange(25).reshape(5, 5)
>>> domain = np.identity(3)
>>> x
array([[ 0,  1,  2,  3,  4],
 [ 5,  6,  7,  8,  9],
 [10, 11, 12, 13, 14],
 [15, 16, 17, 18, 19],
 [20, 21, 22, 23, 24]])
>>> signal.order_filter(x, domain, 0)
array([[  0.,   0.,   0.,   0.,   0.],
 [  0.,   0.,   1.,   2.,   0.],
 [  0.,   5.,   6.,   7.,   0.],
 [  0.,  10.,  11.,  12.,   0.],
 [  0.,   0.,   0.,   0.,   0.]])
>>> signal.order_filter(x, domain, 2)
array([[  6.,   7.,   8.,   9.,   4.],
 [ 11.,  12.,  13.,  14.,   9.],
 [ 16.,  17.,  18.,  19.,  14.],
 [ 21.,  22.,  23.,  24.,  19.],
 [ 20.,  21.,  22.,  23.,  24.]]) 

scipy.signal.medfilt

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.medfilt.html#scipy.signal.medfilt

scipy.signal.medfilt(volume, kernel_size=None)

对 N 维数组执行中值滤波。

使用由kernel_size 给定的局部窗口大小对输入数组应用中值滤波。数组将自动填充零。

参数：

volume 数组形式

一个 N 维输入数组。

kernel_size 数组形式，可选

标量或长度为 N 的列表，指定每个维度中中值滤波窗口的大小。kernel_size 的元素应为奇数。如果kernel_size 是标量，则在每个维度上使用此标量作为大小。每个维度的默认大小为 3。

返回：

out ndarray

一个与输入大小相同的数组，包含中值滤波后的结果。

警告：

用户警告

如果数组大小在任何维度上小于内核大小

另请参阅

scipy.ndimage.median_filter

scipy.signal.medfilt2d

注意事项

更通用的函数scipy.ndimage.median_filter 具有更有效的中值滤波实现，因此运行速度更快。

对于具有uint8、float32或float64数据类型的二维图像，专用函数scipy.signal.medfilt2d 可能更快。

scipy.signal.medfilt2d

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.medfilt2d.html#scipy.signal.medfilt2d

scipy.signal.medfilt2d(input, kernel_size=3)

对 2 维数组进行中值滤波。

使用由kernel_size（必须为奇数）给定的局部窗口大小对input数组应用中值滤波。数组会自动进行零填充。

参数：

input数组型

一个 2 维输入数组。

kernel_size数组型，可选

标量或长度为 2 的列表，分别指定每个维度中的中值滤波窗口大小。kernel_size的元素应为奇数。如果kernel_size是标量，则在每个维度上使用此标量作为大小。默认为大小为(3, 3)的核。

返回：

out ndarray

与输入大小相同的数组，其中包含中值滤波的结果。

另请参阅

scipy.ndimage.median_filter

注意事项

当输入的数据类型为uint8、float32或float64时，此方法比medfilt更快；对于其他类型，会回退到medfilt。在某些情况下，scipy.ndimage.median_filter可能比此函数更快。

示例

>>> import numpy as np
>>> from scipy import signal
>>> x = np.arange(25).reshape(5, 5)
>>> x
array([[ 0,  1,  2,  3,  4],
 [ 5,  6,  7,  8,  9],
 [10, 11, 12, 13, 14],
 [15, 16, 17, 18, 19],
 [20, 21, 22, 23, 24]]) 

将 i，j 替换为默认 5*5 窗口中的中值

>>> signal.medfilt2d(x, kernel_size=5)
array([[ 0,  0,  2,  0,  0],
 [ 0,  3,  7,  4,  0],
 [ 2,  8, 12,  9,  4],
 [ 0,  8, 12,  9,  0],
 [ 0,  0, 12,  0,  0]]) 

将 i，j 替换为默认 3*3 窗口中的中值

>>> signal.medfilt2d(x)
array([[ 0,  1,  2,  3,  0],
 [ 1,  6,  7,  8,  4],
 [ 6, 11, 12, 13,  9],
 [11, 16, 17, 18, 14],
 [ 0, 16, 17, 18,  0]]) 

将 i，j 替换为默认 5*3 窗口中的中值

>>> signal.medfilt2d(x, kernel_size=[5,3])
array([[ 0,  1,  2,  3,  0],
 [ 0,  6,  7,  8,  3],
 [ 5, 11, 12, 13,  8],
 [ 5, 11, 12, 13,  8],
 [ 0, 11, 12, 13,  0]]) 

将 i，j 替换为默认 3*5 窗口中的中值

>>> signal.medfilt2d(x, kernel_size=[3,5])
array([[ 0,  0,  2,  1,  0],
 [ 1,  5,  7,  6,  3],
 [ 6, 10, 12, 11,  8],
 [11, 15, 17, 16, 13],
 [ 0, 15, 17, 16,  0]]) 

如示例中所示，#内核数量必须是奇数，不能超过原始数组维度

scipy.signal.wiener

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.wiener.html#scipy.signal.wiener

scipy.signal.wiener(im, mysize=None, noise=None)

对 N 维数组执行 Wiener 滤波。

对 N 维数组im应用 Wiener 滤波器。

参数：

imndarray

一个 N 维数组。

mysizeint 或 array_like，可选

一个标量或者一个长度为 N 的列表，其中的元素指定 Wiener 滤波器在每个维度上的窗口大小。mysize 的元素应为奇数。如果 mysize 是标量，则在每个维度上使用此标量作为大小。

noisefloat，可选

用于计算噪声功率。如果为 None，则噪声被估计为输入的局部方差的平均值。

返回：

outndarray

Wiener 滤波后的结果与im具有相同的形状。

注意

此实现类似于 Matlab/Octave 中的 wiener2。更多细节参见[1]

参考文献

[1]

Lim, Jae S., 《二维信号与图像处理》，Englewood Cliffs, NJ, Prentice Hall, 1990, p. 548.

示例

>>> from scipy.datasets import face
>>> from scipy.signal import wiener
>>> import matplotlib.pyplot as plt
>>> import numpy as np
>>> rng = np.random.default_rng()
>>> img = rng.random((40, 40))    #Create a random image
>>> filtered_img = wiener(img, (5, 5))  #Filter the image
>>> f, (plot1, plot2) = plt.subplots(1, 2)
>>> plot1.imshow(img)
>>> plot2.imshow(filtered_img)
>>> plt.show() 

scipy.signal.symiirorder1

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.symiirorder1.html#scipy.signal.symiirorder1

scipy.signal.symiirorder1(input, c0, z1, precision=-1.0)

使用一系列一阶段级联实现具有镜像对称边界条件的平滑 IIR 滤波器。第二个阶段使用了反转序列。这实现了以下传递函数和镜像对称边界条件的系统：

 c0              
H(z) = ---------------------    
        (1-z1/z) (1 - z1 z) 

结果信号将具有镜像对称的边界条件。

参数：

inputndarray

输入信号。

c0, z1scalar

传递函数中的参数。

precision

根据镜像对称输入计算递归滤波器初始条件的精度。

返回：

outputndarray

过滤后的信号。

scipy.signal.symiirorder2

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.symiirorder2.html#scipy.signal.symiirorder2

scipy.signal.symiirorder2(input, r, omega, precision=-1.0)

使用二阶段级联实现具有镜像对称边界条件的平滑 IIR 滤波器。第二阶段使用了反转的序列。这实现了以下传递函数：

 cs²
H(z) = ---------------------------------------
       (1 - a2/z - a3/z²) (1 - a2 z - a3 z² ) 

其中：

a2 = (2 r cos omega)
a3 = - r²
cs = 1 - 2 r cos omega + r² 

参数：

输入ndarray

输入信号。

r, omegafloat

传递函数中的参数。

精度float

指定根据镜像对称输入计算递归滤波器初始条件的精度。

返回：

输出ndarray

过滤后的信号。

scipy.signal.lfilter

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.lfilter.html#scipy.signal.lfilter

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

沿着一维数据使用 IIR 或 FIR 滤波器滤波。

使用数字滤波器滤波数据序列 x。这适用于许多基本数据类型（包括对象类型）。滤波器是标准差分方程的直接 II 转置实现（见注意事项）。

函数 sosfilt（使用 output='sos' 进行滤波器设计）应优先于 lfilter 用于大多数滤波任务，因为二阶段节拍具有较少的数值问题。

参数：

b 数组样

1-D 序列中的分子系数向量。

a 数组样

1-D 序列中的分母系数向量。如果 a[0] 不为 1，则 a 和 b 都将被 a[0] 标准化。

x 数组样

N 维输入数组。

axis 整型，可选

应用线性滤波器的输入数据数组的轴。该滤波器应用于此轴上的每个子数组。默认为 -1。

zi 数组样，可选

滤波器延迟的初始条件。它是长度为 max(len(a), len(b)) - 1 的向量（或者对于 N 维输入是向量数组）。如果 zi 为 None 或未给出，则假定初始休息。详见 lfiltic 获取更多信息。

返回：

y 数组

数字滤波器的输出。

zf 数组，可选

如果 zi 为 None，则不返回，否则 zf 包含最终滤波器延迟值。

另见

lfiltic

构建 lfilter 的初始条件。

lfilter_zi

计算 lfilter 的初始状态（阶跃响应的稳态）。

filtfilt

前向后向滤波器，以获得零相位滤波器。

savgol_filter

Savitzky-Golay 滤波器。

sosfilt

使用级联二阶段节拍滤波数据。

sosfiltfilt

使用二阶段节拍进行前向后向滤波器。

注意事项

该滤波函数实现为直接 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 变换域中的有理传递函数为：

 -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() 

scipy.signal.lfiltic

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.lfiltic.html#scipy.signal.lfiltic

scipy.signal.lfiltic(b, a, y, x=None)

为 lfilter 构造输入和输出向量的初始条件。

给定线性滤波器 (b, a) 和输出 y 以及输入 x 的初始条件，返回 lfilter 使用的状态向量 zi 的初始条件，用于生成输出。

参数：

barray_like

线性滤波器项。

aarray_like

线性滤波器项。

yarray_like

初始条件。

如果 N = len(a) - 1，则 y = {y[-1], y[-2], ..., y[-N]}。

如果 y 太短，会用零填充。

xarray_like，可选

初始条件。

如果 M = len(b) - 1，则 x = {x[-1], x[-2], ..., x[-M]}。

如果没有给出 x，则假设其初始条件为零。

如果 x 太短，会用零填充。

返回：

zindarray

状态向量 zi = {z_0[-1], z_1[-1], ..., z_K-1[-1]}，其中 K = max(M, N)。

另请参见

lfilter，lfilter_zi

scipy.signal.lfilter_zi

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.lfilter_zi.html#scipy.signal.lfilter_zi

scipy.signal.lfilter_zi(b, a)

构造 lfilter 的阶跃响应稳态的初始条件。

为 lfilter 函数计算一个初始状态 zi，对应于阶跃响应的稳态。

此函数的典型用途是设置初始状态，使得滤波器的输出从与待滤波信号的第一个元素相同的值开始。

参数：

b, a array_like (1-D)

IIR 滤波器系数。详见 lfilter。

返回值：

zi 1-D ndarray

滤波器的初始状态。

参见

lfilter，lfiltic，filtfilt

注意事项

具有阶数 m 的线性滤波器具有状态空间表示 (A, B, C, D)，滤波器的输出 y 可以表示为：

z(n+1) = A*z(n) + B*x(n)
y(n)   = C*z(n) + D*x(n) 

其中 z(n) 是长度为 m 的向量，A 的形状为 (m, m)，B 的形状为 (m, 1)，C 的形状为 (1, m)，D 的形状为 (1, 1)（假设 x(n) 是标量）。lfilter_zi 解决：

zi = A*zi + B 

换句话说，它找到了哪个初始条件，使得对全 1 输入的响应是一个常数。

给定滤波器系数 a 和 b，用于线性滤波器的转置直接形式 II 实现的状态空间矩阵，即 scipy.signal.lfilter 使用的实现方式如下：

A = scipy.linalg.companion(a).T
B = b[1:] - a[1:]*b[0] 

假设 a[0] 为 1.0；如果 a[0] 不是 1，a 和 b 首先将被除以 a[0]。

示例

下面的代码创建一个低通 Butterworth 滤波器。然后将该滤波器应用于一个所有值均为 1.0 的数组；输出也全部为 1.0，符合低通滤波器的预期行为。如果未提供 lfilter 的 zi 参数，输出将显示瞬态信号。

>>> from numpy import array, ones
>>> from scipy.signal import lfilter, lfilter_zi, butter
>>> b, a = butter(5, 0.25)
>>> zi = lfilter_zi(b, a)
>>> y, zo = lfilter(b, a, ones(10), zi=zi)
>>> y
array([1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.,  1.]) 

另一个示例：

>>> x = array([0.5, 0.5, 0.5, 0.0, 0.0, 0.0, 0.0])
>>> y, zf = lfilter(b, a, x, zi=zi*x[0])
>>> y
array([ 0.5       ,  0.5       ,  0.5       ,  0.49836039,  0.48610528,
 0.44399389,  0.35505241]) 

注意，lfilter 的 zi 参数是通过 lfilter_zi 计算并缩放为 x[0]。然后输出 y 在输入从 0.5 下降到 0.0 之前没有瞬态信号。

scipy.signal.filtfilt

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.filtfilt.html#scipy.signal.filtfilt

scipy.signal.filtfilt(b, a, x, axis=-1, padtype='odd', padlen=None, method='pad', irlen=None)

对信号应用线性数字滤波器，向前和向后。

此函数对信号应用线性数字滤波器两次，一次向前，一次向后。组合的滤波器具有零相位和原始滤波器两倍的滤波器阶数。

此函数提供处理信号边缘的选项。

函数 sosfiltfilt（和使用 output='sos' 进行滤波器设计）应优先于 filtfilt 用于大多数滤波任务，因为二阶段节省去了更多的数值问题。

参数:

b(N,) array_like

过滤器的分子系数向量。

a(N,) array_like

过滤器的分母系数向量。如果 a[0] 不为 1，则 a 和 b 都将被 a[0] 归一化。

xarray_like

需要过滤的数据数组。

axisint, 可选

要应用滤波器的 x 的轴。默认为 -1。

padtypestr 或 None, 可选

必须是 'odd', 'even', 'constant' 或 None。这决定了要应用滤波器的填充信号的扩展类型。如果 padtype 是 None，则不使用填充。默认值为 'odd'。

padlenint 或 None, 可选

在 x 的两端的 axis 扩展元素的数量。此值必须小于 x.shape[axis] - 1。 padlen=0 表示不填充。默认值为 3 * max(len(a), len(b))。

methodstr, 可选

决定信号边缘处理方法的方法，可以是 “pad” 或 “gust”。当 method 是 “pad” 时，信号被填充；填充的类型由 padtype 和 padlen 决定，irlen 被忽略。当 method 是 “gust” 时，使用 Gustafsson 方法，padtype 和 padlen 被忽略。

irlenint 或 None, 可选

当 method 是 “gust” 时，irlen 指定滤波器的脉冲响应长度。如果 irlen 是 None，则不会忽略脉冲响应的任何部分。对于长信号，指定 irlen 可显著改善滤波器的性能。

返回:

yndarray

输出的滤波后的形状与 x 相同。

另请参阅

sosfiltfilt, lfilter_zi, lfilter, lfiltic, savgol_filter, sosfilt

注意

当 method 为 “pad” 时，函数在给定的轴上以三种方式之一填充数据：奇数、偶数或常数。奇数和偶数扩展在数据端点处具有相应的对称性。常数扩展使用端点处的值延伸数据。在前向和后向传递中，滤波器的初始条件通过使用 lfilter_zi 找到，并通过扩展数据的端点进行缩放。

当 method 为 “gust” 时，使用 Gustafsson 方法 [1]。选择前向和后向传递的初始条件，以便前后向滤波器给出与后前向滤波器相同的结果。

在 scipy 版本 0.16.0 中添加了使用 Gustaffson 方法的选项。

参考文献

[1]

F. Gustaffson，“确定前向-后向滤波中的初始状态”，信号处理交易，Vol. 46，pp. 988-992，1996。

示例

示例将使用 scipy.signal 中的多个函数。

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt 

首先，我们创建一个持续一秒钟的信号，这个信号是两个纯正弦波（频率分别为 5 Hz 和 250 Hz）的和，采样率为 2000 Hz。

>>> t = np.linspace(0, 1.0, 2001)
>>> xlow = np.sin(2 * np.pi * 5 * t)
>>> xhigh = np.sin(2 * np.pi * 250 * t)
>>> x = xlow + xhigh 

现在创建一个低通巴特沃斯滤波器，截止频率为 0.125 倍的奈奎斯特频率，即 125 Hz，并用 filtfilt 应用于 x。结果应该是近似于 xlow，没有相移。

>>> b, a = signal.butter(8, 0.125)
>>> y = signal.filtfilt(b, a, x, padlen=150)
>>> np.abs(y - xlow).max()
9.1086182074789912e-06 

对于这个人工示例，我们得到了一个相当干净的结果，因为奇数扩展是精确的，并且通过适度长的填充，滤波器的瞬态效应在实际数据到达时已经消失。一般来说，边缘处的瞬态效应是不可避免的。

下面的示例演示了选项 method="gust"。

首先，创建一个滤波器。

>>> b, a = signal.ellip(4, 0.01, 120, 0.125)  # Filter to be applied. 

sig 是一个要进行滤波的随机输入信号。

>>> rng = np.random.default_rng()
>>> n = 60
>>> sig = rng.standard_normal(n)**3 + 3*rng.standard_normal(n).cumsum() 

分别对 sig 应用 filtfilt，一次使用 Gustafsson 方法，一次使用填充，并绘制结果进行比较。

>>> fgust = signal.filtfilt(b, a, sig, method="gust")
>>> fpad = signal.filtfilt(b, a, sig, padlen=50)
>>> plt.plot(sig, 'k-', label='input')
>>> plt.plot(fgust, 'b-', linewidth=4, label='gust')
>>> plt.plot(fpad, 'c-', linewidth=1.5, label='pad')
>>> plt.legend(loc='best')
>>> plt.show() 

irlen 参数可用于改善 Gustafsson 方法的性能。

估计滤波器的脉冲响应长度。

>>> z, p, k = signal.tf2zpk(b, a)
>>> eps = 1e-9
>>> r = np.max(np.abs(p))
>>> approx_impulse_len = int(np.ceil(np.log(eps) / np.log(r)))
>>> approx_impulse_len
137 

对较长的信号应用滤波器，有或没有 irlen 参数。y1 和 y2 之间的差异很小。对于长信号，使用 irlen 可显著提高性能。

>>> x = rng.standard_normal(4000)
>>> y1 = signal.filtfilt(b, a, x, method='gust')
>>> y2 = signal.filtfilt(b, a, x, method='gust', irlen=approx_impulse_len)
>>> print(np.max(np.abs(y1 - y2)))
2.875334415008979e-10 

scipy.signal.savgol_filter

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.savgol_filter.html#scipy.signal.savgol_filter

scipy.signal.savgol_filter(x, window_length, polyorder, deriv=0, delta=1.0, axis=-1, mode='interp', cval=0.0)

对数组应用 Savitzky-Golay 滤波器。

这是一个 1-D 滤波器。如果x的维度大于 1，则axis确定应用滤波器的轴。

参数：

xarray_like

要过滤的数据。如果x不是单精度或双精度浮点数组，则在过滤之前将其转换为numpy.float64类型。

window_lengthint

滤波窗口的长度（即系数的数量）。如果mode为‘interp’，window_length必须小于或等于x的大小。

polyorderint

用于拟合样本的多项式的阶数。polyorder必须小于window_length。

derivint，可选

要计算的导数阶数。这必须是非负整数。默认值为 0，表示在不进行微分的情况下过滤数据。

deltafloat，可选

应用过滤器的样本间距。仅在 deriv > 0 时使用。默认值为 1.0。

axisint，可选

要应用过滤器的数组x的轴。默认值为-1。

modestr，可选

必须为‘mirror’、‘constant’、‘nearest’、‘wrap’或‘interp’。这决定了要用于填充信号的填充类型。当mode为‘constant’时，填充值由cval给出。有关‘mirror’、‘constant’、‘wrap’和‘nearest’的更多详细信息，请参阅注释。当选择‘interp’模式（默认情况下）时，不使用扩展。相反，对边缘的最后window_length个值拟合一个polyorder次多项式，并使用此多项式来评估最后window_length // 2个输出值。

cvalscalar，可选

如果mode为‘constant’，则在输入的边缘之外填充的值。默认值为 0.0。

返回：

yndarray，与x相同的形状

过滤后的数据。

另请参阅

savgol_coeffs

注意事项

mode选项的详细信息：

‘mirror’：

以相反顺序重复边缘处的值。不包括最接近边缘的值。

‘nearest’：

扩展包含最接近的输入值。

‘constant’：

扩展包含由cval参数给出的值。

‘wrap’：

扩展包含数组另一端的值。

例如，如果输入为[1, 2, 3, 4, 5, 6, 7, 8]，window_length为 7，则以下显示了各种mode选项的扩展数据（假设cval为 0）：

mode       |   Ext   |         Input          |   Ext
-----------+---------+------------------------+---------
'mirror'   | 4  3  2 | 1  2  3  4  5  6  7  8 | 7  6  5
'nearest'  | 1  1  1 | 1  2  3  4  5  6  7  8 | 8  8  8
'constant' | 0  0  0 | 1  2  3  4  5  6  7  8 | 0  0  0
'wrap'     | 6  7  8 | 1  2  3  4  5  6  7  8 | 1  2  3 

从版本 0.14.0 开始。

示例

>>> import numpy as np
>>> from scipy.signal import savgol_filter
>>> np.set_printoptions(precision=2)  # For compact display.
>>> x = np.array([2, 2, 5, 2, 1, 0, 1, 4, 9]) 

使用窗口长度为 5 和二次多项式进行滤波。对所有其他参数使用默认值。

>>> savgol_filter(x, 5, 2)
array([1.66, 3.17, 3.54, 2.86, 0.66, 0.17, 1\.  , 4\.  , 9\.  ]) 

注意，x 中的最后五个值是抛物线的样本，因此当 mode=’interp’（默认情况）与 polyorder=2 结合使用时，最后三个值保持不变。与 mode=’nearest’ 相比，例如：

>>> savgol_filter(x, 5, 2, mode='nearest')
array([1.74, 3.03, 3.54, 2.86, 0.66, 0.17, 1\.  , 4.6 , 7.97]) 

scipy.signal.deconvolve

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.deconvolve.html#scipy.signal.deconvolve

scipy.signal.deconvolve(signal, divisor)

使用逆滤波器将signal中的divisor去卷积出来。

返回商和余数，使得signal = convolve(divisor, quotient) + remainder

参数：

signal(N,) 数组型

信号数据，通常是记录的信号。

divisor(N,) 数组型

除数数据，通常是应用于原始信号的冲激响应或滤波器

返回：

quotientndarray

商，通常是恢复的原始信号。

remainderndarray

余数

另请参阅

numpy.polydiv

执行多项式除法（相同操作，但也接受 poly1d 对象）

示例

去卷积已经被过滤的信号：

>>> from scipy import signal
>>> original = [0, 1, 0, 0, 1, 1, 0, 0]
>>> impulse_response = [2, 1]
>>> recorded = signal.convolve(impulse_response, original)
>>> recorded
array([0, 2, 1, 0, 2, 3, 1, 0, 0])
>>> recovered, remainder = signal.deconvolve(recorded, impulse_response)
>>> recovered
array([ 0.,  1.,  0.,  0.,  1.,  1.,  0.,  0.]) 

scipy.signal.sosfilt

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.sosfilt.html#scipy.signal.sosfilt

scipy.signal.sosfilt(sos, x, axis=-1, zi=None)

使用级联的二阶段进行数据滤波。

使用数字 IIR 滤波器 sos 过滤数据序列 x。

参数：

sos 数组类型

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

x 数组类型

输入数组的 N 维。

axis 整数，可选

应用线性滤波器的输入数据数组的轴。该滤波器应用于沿此轴的每个子数组。默认为 -1。

zi 数组类型，可选

级联滤波器延迟的初始条件。它是形状为 (n_sections, ..., 2, ...) 的（至少是二维的）向量，其中 ..., 2, ... 表示 x 的形状，但将 x.shape[axis] 替换为 2。如果 zi 为 None 或未给出，则假定初始休息（即全部为零）。请注意，这些初始条件与 lfiltic 或 lfilter_zi 给出的初始条件不同。

返回：

y 数组

数字滤波器的输出。

zf 数组，可选

如果 zi 为 None，则不返回，否则 zf 保存最终的滤波器延迟值。

另请参阅：

zpk2sos, sos2zpk, sosfilt_zi, sosfiltfilt, sosfreqz

注意事项

该滤波器函数实现为直接 II 转置结构的多个二阶滤波器的序列。它旨在减少高阶滤波器的数值精度误差。

0.16.0 版本的新功能。

示例：

使用 lfilter 和 sosfilt 绘制一个 13 阶滤波器的脉冲响应，显示尝试在单个阶段进行 13 阶滤波器时产生的不稳定性（数值误差使一些极点超出单位圆）：

>>> import matplotlib.pyplot as plt
>>> from scipy import signal
>>> b, a = signal.ellip(13, 0.009, 80, 0.05, output='ba')
>>> sos = signal.ellip(13, 0.009, 80, 0.05, output='sos')
>>> x = signal.unit_impulse(700)
>>> y_tf = signal.lfilter(b, a, x)
>>> y_sos = signal.sosfilt(sos, x)
>>> plt.plot(y_tf, 'r', label='TF')
>>> plt.plot(y_sos, 'k', label='SOS')
>>> plt.legend(loc='best')
>>> plt.show() 

scipy.signal.sosfilt_zi

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.sosfilt_zi.html#scipy.signal.sosfilt_zi

scipy.signal.sosfilt_zi(sos)

为阶跃响应稳态的sosfilt构造初始条件。

计算sosfilt函数的初始状态zi，该状态对应于阶跃响应的稳态。

该函数的典型用法是设置初始状态，使滤波器的输出与要滤波信号的第一个元素的值相同。

参数：

sos数组样式

第二阶滤波器系数数组，必须具有形状(n_sections, 6)。参见sosfilt以获取 SOS 滤波器格式规范。

返回：

zi数组

适用于与sosfilt一起使用的初始条件，形状为(n_sections, 2)。

另请参阅

sosfilt, zpk2sos

注释

自 0.16.0 版新增。

示例

对 0 时刻开始的矩形脉冲进行滤波，使用和不使用scipy.signal.sosfilt的zi参数。

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt 

>>> sos = signal.butter(9, 0.125, output='sos')
>>> zi = signal.sosfilt_zi(sos)
>>> x = (np.arange(250) < 100).astype(int)
>>> f1 = signal.sosfilt(sos, x)
>>> f2, zo = signal.sosfilt(sos, x, zi=zi) 

>>> plt.plot(x, 'k--', label='x')
>>> plt.plot(f1, 'b', alpha=0.5, linewidth=2, label='filtered')
>>> plt.plot(f2, 'g', alpha=0.25, linewidth=4, label='filtered with zi')
>>> plt.legend(loc='best')
>>> plt.show() 

scipy.signal.sosfiltfilt

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.sosfiltfilt.html#scipy.signal.sosfiltfilt

scipy.signal.sosfiltfilt(sos, x, axis=-1, padtype='odd', padlen=None)

使用级联二阶节创建前向-后向数字滤波器。

更完整信息，请参见filtfilt方法。

参数：

sosarray_like

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

xarray_like

要进行滤波的数据数组。

axisint，可选

应用滤波器的x的轴。默认为-1。

padtypestr 或 None，可选

必须为'odd'、'even'、'constant'或 None。这决定要用于填充信号的扩展类型，以便应用滤波器。如果padtype为 None，则不使用填充。默认为'odd'。

padlenint 或 None，可选

在应用滤波器之前，沿axis两端延伸x的元素数。该值必须小于x.shape[axis] - 1。padlen=0表示无填充。默认值为：

3 * (2 * len(sos) + 1 - min((sos[:, 2] == 0).sum(),
                            (sos[:, 5] == 0).sum())) 

最后的额外减法试图补偿在原点处的极点和零点（例如，对于奇阶滤波器），以产生与用scipy.signal函数构建的二阶节滤波器的filtfilt相当的padlen估计。

返回：

yndarray

与x具有相同形状的滤波输出。

另请参阅

filtfilt, sosfilt, sosfilt_zi, sosfreqz

注意事项

新版本 0.18.0 中新增。

示例

>>> import numpy as np
>>> from scipy.signal import sosfiltfilt, butter
>>> import matplotlib.pyplot as plt
>>> rng = np.random.default_rng() 

创建一个有趣的信号以进行滤波。

>>> n = 201
>>> t = np.linspace(0, 1, n)
>>> x = 1 + (t < 0.5) - 0.25*t**2 + 0.05*rng.standard_normal(n) 

创建一个低通巴特沃斯滤波器，并用它来滤波x。

>>> sos = butter(4, 0.125, output='sos')
>>> y = sosfiltfilt(sos, x) 

为了比较，使用sosfilt应用一个 8 阶滤波器。滤波器使用x的前四个值的均值进行初始化。

>>> from scipy.signal import sosfilt, sosfilt_zi
>>> sos8 = butter(8, 0.125, output='sos')
>>> zi = x[:4].mean() * sosfilt_zi(sos8)
>>> y2, zo = sosfilt(sos8, x, zi=zi) 

绘制结果。注意y的相位与输入匹配，而y2存在显著的相位延迟。

>>> plt.plot(t, x, alpha=0.5, label='x(t)')
>>> plt.plot(t, y, label='y(t)')
>>> plt.plot(t, y2, label='y2(t)')
>>> plt.legend(framealpha=1, shadow=True)
>>> plt.grid(alpha=0.25)
>>> plt.xlabel('t')
>>> plt.show() 

scipy.signal.hilbert

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.hilbert.html#scipy.signal.hilbert

scipy.signal.hilbert(x, N=None, axis=-1)

使用希尔伯特变换计算解析信号。

默认情况下，转换沿最后一个轴执行。

参数：

xarray_like

信号数据。必须是实数。

Nint，可选

傅里叶分量的数量。默认值为x.shape[axis]

axisint，可选

变换的轴线。默认值为-1。

返回：

xandarray

x的解析信号，沿axis的每个 1-D 数组

注意事项

信号x(t)的解析信号x_a(t)是：

[x_a = F^{-1}(F(x) 2U) = x + i y]

其中F是傅里叶变换，U是单位阶跃函数，y是x的希尔伯特变换。[1]

换句话说，负频谱的负半部分被置零，将实值信号转换为复杂信号。希尔伯特变换信号可以通过np.imag(hilbert(x))获取，原始信号可以通过np.real(hilbert(x))获取。

参考文献

[1]

维基百科，“解析信号”。en.wikipedia.org/wiki/Analytic_signal

[2]

Leon Cohen，“时频分析”，1995. 第二章。

[3]

Alan V. Oppenheim, Ronald W. Schafer. Discrete-Time Signal Processing, Third Edition, 2009. Chapter 12. ISBN 13: 978-1292-02572-8

示例

在这个例子中，我们使用希尔伯特变换来确定调幅信号的幅度包络和即时频率。

>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> from scipy.signal import hilbert, chirp 

>>> duration = 1.0
>>> fs = 400.0
>>> samples = int(fs*duration)
>>> t = np.arange(samples) / fs 

我们创建一个从 20 Hz 到 100 Hz 频率增加并应用幅度调制的啁啾声。

>>> signal = chirp(t, 20.0, t[-1], 100.0)
>>> signal *= (1.0 + 0.5 * np.sin(2.0*np.pi*3.0*t) ) 

幅度包络由解析信号的幅度给出。通过将即时相位相对于时间进行微分，即时频率可以获得。即时相位对应于解析信号的相位角。

>>> analytic_signal = hilbert(signal)
>>> amplitude_envelope = np.abs(analytic_signal)
>>> instantaneous_phase = np.unwrap(np.angle(analytic_signal))
>>> instantaneous_frequency = (np.diff(instantaneous_phase) /
...                            (2.0*np.pi) * fs) 

>>> fig, (ax0, ax1) = plt.subplots(nrows=2)
>>> ax0.plot(t, signal, label='signal')
>>> ax0.plot(t, amplitude_envelope, label='envelope')
>>> ax0.set_xlabel("time in seconds")
>>> ax0.legend()
>>> ax1.plot(t[1:], instantaneous_frequency)
>>> ax1.set_xlabel("time in seconds")
>>> ax1.set_ylim(0.0, 120.0)
>>> fig.tight_layout() 

scipy.signal.hilbert2

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.hilbert2.html#scipy.signal.hilbert2

scipy.signal.hilbert2(x, N=None)

计算 x 的‘2-D’ 解析信号

参数：

xarray_like

二维信号数据。

Nint 或两个 int 的元组，可选

傅里叶分量的数量。默认为 x.shape

返回：

xandarray

沿轴（0,1）取 x 的解析信号。

参考资料

[1]

维基百科，“解析信号”，en.wikipedia.org/wiki/Analytic_signal

scipy.signal.decimate

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.decimate.html#scipy.signal.decimate

scipy.signal.decimate(x, q, n=None, ftype='iir', axis=-1, zero_phase=True)

在应用抗混叠滤波器后对信号进行降采样。

默认情况下，使用阶数为 8 的 Chebyshev I 型滤波器。如果ftype为‘fir’，则使用 30 点 Hamming 窗口的 FIR 滤波器。

参数：

xarray_like

要降采样的信号，作为 N 维数组。

qint

下采样因子。当使用 IIR 下采样时，建议对高于 13 的下采样因子多次调用decimate。

nint，可选

滤波器的阶数（对于‘fir’来说是长度减 1）。对于‘iir’默认为 8，对于‘fir’是下采样因子的 20 倍。

ftypestr {‘iir’，‘fir’}或dlti实例，可选

如果是‘iir’或‘fir’，则指定低通滤波器的类型。如果是dlti对象的实例，则使用该对象在降采样之前进行滤波。

axisint，可选

要降采样的轴。

zero_phasebool，可选

当使用 IIR 滤波器时，通过使用filtfilt而不是lfilter进行滤波，并将输出向后移动滤波器的群延迟来防止相位移动。通常建议使用默认值True，因为通常不希望出现相位移动。

在 0.18.0 版本中新增。

返回：

yndarray

降采样信号。

参见

resample

使用 FFT 方法上下采样。

resample_poly

使用多相滤波和 FIR 滤波器重采样。

注意

zero_phase关键字在 0.18.0 版本中添加。允许使用dlti实例作为ftype在 0.18.0 版本中添加。

示例

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt 

定义波参数。

>>> wave_duration = 3
>>> sample_rate = 100
>>> freq = 2
>>> q = 5 

计算样本数。

>>> samples = wave_duration*sample_rate
>>> samples_decimated = int(samples/q) 

创建余弦波。

>>> x = np.linspace(0, wave_duration, samples, endpoint=False)
>>> y = np.cos(x*np.pi*freq*2) 

降采样余弦波。

>>> ydem = signal.decimate(y, q)
>>> xnew = np.linspace(0, wave_duration, samples_decimated, endpoint=False) 

绘制原始波和降采样波。

>>> plt.plot(x, y, '.-', xnew, ydem, 'o-')
>>> plt.xlabel('Time, Seconds')
>>> plt.legend(['data', 'decimated'], loc='best')
>>> plt.show() 

scipy.signal.detrend

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.detrend.html#scipy.signal.detrend

scipy.signal.detrend(data, axis=-1, type='linear', bp=0, overwrite_data=False)

从数据中去除沿轴的线性趋势。

参数：

data：数组样式

输入数据。

axis：整数，可选

数据去趋势的轴。默认为最后一个轴（-1）。

type：{‘linear’, ‘constant’}，可选

去趋势的类型。如果type == 'linear'（默认），则从data中减去线性最小二乘拟合的结果。如果type == 'constant'，则仅减去data的平均值。

bp：整数数组，可选

断点序列。如果指定，则在data中每个部分之间执行单独的线性拟合。断点被指定为data的索引。当type == 'linear'时，此参数才会生效。

overwrite_data：布尔值，可选

如果为 True，则执行就地去趋势并避免复制。默认为 False。

返回：

ret：ndarray

去趋势后的输入数据。

示例

>>> import numpy as np
>>> from scipy import signal
>>> rng = np.random.default_rng()
>>> npoints = 1000
>>> noise = rng.standard_normal(npoints)
>>> x = 3 + 2*np.linspace(0, 1, npoints) + noise
>>> (signal.detrend(x) - noise).max()
0.06  # random 

scipy.signal.resample

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.resample.html#scipy.signal.resample

scipy.signal.resample(x, num, t=None, axis=0, window=None, domain='time')

使用傅立叶方法沿给定轴将 x 重采样为 num 个样本。

重采样信号从与 x 相同的值开始，但采样间隔为 len(x) / num * (spacing of x)。由于使用了傅立叶方法，信号被假定为周期性的。

参数：

x 数组

要重采样的数据。

num 整数

重采样信号中的样本数。

t 数组，可选

如果给定 t，则假定它是与 x 中信号数据相关联的等间隔采样位置。

axis 整数，可选

被重采样的 x 的轴。默认为 0。

window 数组、可调用对象、字符串、浮点数或元组，可选

指定应用于信号的傅立叶域中的窗口。详情见下文。

domain 字符串，可选

指示输入 x 的域的字符串：time 将输入 x 视为时域（默认），freq 将输入 x 视为频域。

返回：

resampled_x 或 (resampled_x, resampled_t)

要么是重采样后的数组，要么（如果给定了 t）是一个包含重采样后的数组和相应重采样位置的元组。

另请参见

decimate

在应用 FIR 或 IIR 滤波器后对信号进行下采样。

resample_poly

使用多相滤波和 FIR 滤波器进行重采样。

注释

参数 window 控制傅立叶域窗口，在零填充前锐化傅立叶频谱，以减轻对未意图作为带限信号解释的采样信号的响应。

如果 window 是一个函数，则调用它并传入一个指示频率区间的输入向量（即 fftfreq(x.shape[axis])）。

如果 window 是与 x.shape[axis] 长度相同的数组，则假定它是要直接在傅立叶域中应用的窗口（带有直流分量和低频率优先）。

对于任何其他类型的 window，将调用函数 scipy.signal.get_window 来生成窗口。

返回向量的第一个样本与输入向量的第一个样本相同。样本之间的间距从 dx 变为 dx * len(x) / num。

如果 t 不为 None，则仅用于计算重采样位置 resampled_t。

如前所述，resample 使用 FFT 变换，如果输入或输出样本数较大且为质数，则速度可能会非常慢；参见 scipy.fft.fft。

示例

注意，重采样数据的末尾上升以满足下一个周期的第一个样本：

>>> import numpy as np
>>> from scipy import signal 

>>> x = np.linspace(0, 10, 20, endpoint=False)
>>> y = np.cos(-x**2/6.0)
>>> f = signal.resample(y, 100)
>>> xnew = np.linspace(0, 10, 100, endpoint=False) 

>>> import matplotlib.pyplot as plt
>>> plt.plot(x, y, 'go-', xnew, f, '.-', 10, y[0], 'ro')
>>> plt.legend(['data', 'resampled'], loc='best')
>>> plt.show() 

scipy.signal.resample_poly

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.resample_poly.html#scipy.signal.resample_poly

scipy.signal.resample_poly(x, up, down, axis=0, window=('kaiser', 5.0), padtype='constant', cval=None)

使用多相滤波器沿给定轴对x进行重新采样。

信号x通过因子up上采样，然后应用零相位低通 FIR 滤波器，并通过因子down进行下采样。结果的采样率为原始采样率的up / down倍。在滤波步骤期间，默认情况下假设信号边界外的值为零。

参数：

x类数组

要重新采样的数据。

up整数

上采样因子。

down整数

下采样因子。

axis整数，可选

被重新采样的x的轴。默认为 0。

window字符串、元组或类数组，可选

用于设计低通滤波器的期望窗口，或用于使用的 FIR 滤波器系数。详细信息见下文。

padtype字符串，可选

constant, line, mean, median, maximum, minimum 或其他由 scipy.signal.upfirdn 支持的信号扩展模式。更改对边界外值的假设。如果是 constant，假设为 cval（默认为零）。如果是 line，则假设为由第一个和最后一个点定义的线性趋势。mean、median、maximum 和 minimum 的工作方式与 np.pad 中相同，并假设沿轴的数组边界外的值分别为数组的平均值、中位数、最大值或最小值。

新版本 1.4.0 中新增。

cval浮点数，可选

如果padtype='constant'，则使用的值。默认为零。

新版本 1.4.0 中新增。

返回：

resampled_x数组

重新采样后的数组。

另请参阅

decimate

在应用 FIR 或 IIR 滤波器后对信号进行下采样。

resample

使用 FFT 方法上或下采样。

注意

当样本数较大且为质数时，或者当样本数较大且up和down具有较大的最大公约数时，这种多相方法可能比 Fourier 方法更快。所使用的 FIR 滤波器的长度将取决于max(up, down) // gcd(up, down)，并且多相滤波过程中的操作次数将取决于滤波器长度和down（详见scipy.signal.upfirdn）。

参数window指定了 FIR 低通滤波器的设计。

如果window是一个类似数组，则假定它是 FIR 滤波器系数。请注意，FIR 滤波器应用在上采样步骤之后，因此它应设计用于在原始信号的采样频率上比原始频率高up//gcd(up, down)倍。此函数的输出将与此数组相对于中心，因此如果希望得到零相位滤波器（通常情况），最好传递具有奇数样本数的对称滤波器。

对于任何其他类型的窗口，函数scipy.signal.get_window和scipy.signal.firwin被调用以生成适当的滤波器系数。

返回向量的第一个样本与输入向量的第一个样本相同。样本之间的间距从dx变为dx * down / float(up)。

示例

默认情况下，用于 FFT 方法的重采样数据末端上升以满足下一个周期的第一个样本，并且对于多相方法，接近零：

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt 

>>> x = np.linspace(0, 10, 20, endpoint=False)
>>> y = np.cos(-x**2/6.0)
>>> f_fft = signal.resample(y, 100)
>>> f_poly = signal.resample_poly(y, 100, 20)
>>> xnew = np.linspace(0, 10, 100, endpoint=False) 

>>> plt.plot(xnew, f_fft, 'b.-', xnew, f_poly, 'r.-')
>>> plt.plot(x, y, 'ko-')
>>> plt.plot(10, y[0], 'bo', 10, 0., 'ro')  # boundaries
>>> plt.legend(['resample', 'resamp_poly', 'data'], loc='best')
>>> plt.show() 

默认行为可以通过使用padtype选项进行更改：

>>> N = 5
>>> x = np.linspace(0, 1, N, endpoint=False)
>>> y = 2 + x**2 - 1.7*np.sin(x) + .2*np.cos(11*x)
>>> y2 = 1 + x**3 + 0.1*np.sin(x) + .1*np.cos(11*x)
>>> Y = np.stack([y, y2], axis=-1)
>>> up = 4
>>> xr = np.linspace(0, 1, N*up, endpoint=False) 

>>> y2 = signal.resample_poly(Y, up, 1, padtype='constant')
>>> y3 = signal.resample_poly(Y, up, 1, padtype='mean')
>>> y4 = signal.resample_poly(Y, up, 1, padtype='line') 

>>> for i in [0,1]:
...     plt.figure()
...     plt.plot(xr, y4[:,i], 'g.', label='line')
...     plt.plot(xr, y3[:,i], 'y.', label='mean')
...     plt.plot(xr, y2[:,i], 'r.', label='constant')
...     plt.plot(x, Y[:,i], 'k-')
...     plt.legend()
>>> plt.show() 

scipy.signal.upfirdn

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.upfirdn.html#scipy.signal.upfirdn

scipy.signal.upfirdn(h, x, up=1, down=1, axis=-1, mode='constant', cval=0)

上采样、FIR 滤波和下采样。

参数：

harray_like

1-D FIR（有限冲激响应）滤波器系数。

xarray_like

输入信号数组。

upint，可选

采样率上采样。默认为 1。

downint，可选

降采样率。默认为 1。

axisint，可选

应用线性滤波器的输入数据数组的轴。该滤波器应用于沿此轴的每个子数组。默认为 -1。

modestr，可选

要使用的信号扩展模式。集合 {"constant", "symmetric", "reflect", "edge", "wrap"} 对应于 numpy.pad 提供的模式。"smooth" 根据数组末端的最后两个点的斜率进行平滑扩展。"antireflect" 和 "antisymmetric" 是 "reflect" 和 "symmetric" 的反对称版本。模式 “line” 基于沿 axis 定义的线性趋势扩展信号。

新版本 1.4.0 中新增。

cvalfloat，可选

当 mode == "constant" 时使用的常数值。

新版本 1.4.0 中新增。

返回：

yndarray

输出信号数组。除了 axis 外，维度将与 x 相同，axis 的大小将根据 h、up 和 down 参数变化。

注释

该算法是基于 Vaidyanathan 文本第 129 页所示的块图的实现 [1]（图 4.3-8d）。

通过零插入对 P 的因子上采样、长度为 N 的 FIR 滤波和 Q 的因子下采样的直接方法为每个输出样本的复杂度为 O(N*Q)。此处使用的多相实现的复杂度为 O(N/P)。

新版本 0.18 中新增。

参考文献

[1]

P. P. Vaidyanathan，《Multirate Systems and Filter Banks》，Prentice Hall，1993 年。

示例

简单操作：

>>> import numpy as np
>>> from scipy.signal import upfirdn
>>> upfirdn([1, 1, 1], [1, 1, 1])   # FIR filter
array([ 1.,  2.,  3.,  2.,  1.])
>>> upfirdn([1], [1, 2, 3], 3)  # upsampling with zeros insertion
array([ 1.,  0.,  0.,  2.,  0.,  0.,  3.])
>>> upfirdn([1, 1, 1], [1, 2, 3], 3)  # upsampling with sample-and-hold
array([ 1.,  1.,  1.,  2.,  2.,  2.,  3.,  3.,  3.])
>>> upfirdn([.5, 1, .5], [1, 1, 1], 2)  # linear interpolation
array([ 0.5,  1\. ,  1\. ,  1\. ,  1\. ,  1\. ,  0.5])
>>> upfirdn([1], np.arange(10), 1, 3)  # decimation by 3
array([ 0.,  3.,  6.,  9.])
>>> upfirdn([.5, 1, .5], np.arange(10), 2, 3)  # linear interp, rate 2/3
array([ 0\. ,  1\. ,  2.5,  4\. ,  5.5,  7\. ,  8.5]) 

对多个信号应用单个滤波器：

>>> x = np.reshape(np.arange(8), (4, 2))
>>> x
array([[0, 1],
 [2, 3],
 [4, 5],
 [6, 7]]) 

应用于 x 的最后一个维度：

>>> h = [1, 1]
>>> upfirdn(h, x, 2)
array([[ 0.,  0.,  1.,  1.],
 [ 2.,  2.,  3.,  3.],
 [ 4.,  4.,  5.,  5.],
 [ 6.,  6.,  7.,  7.]]) 

应用于 x 的第 0 维度：

>>> upfirdn(h, x, 2, axis=0)
array([[ 0.,  1.],
 [ 0.,  1.],
 [ 2.,  3.],
 [ 2.,  3.],
 [ 4.,  5.],
 [ 4.,  5.],
 [ 6.,  7.],
 [ 6.,  7.]]) 

scipy.signal.bilinear

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.bilinear.html#scipy.signal.bilinear

scipy.signal.bilinear(b, a, fs=1.0)

使用双线性变换从模拟滤波器返回数字 IIR 滤波器。

将一组极点和零点从模拟 s 平面转换到数字 z 平面，使用 Tustin 方法，其中替换s为2*fs*(z-1) / (z+1)，保持频率响应的形状。

参数：

• barray_like

模拟滤波器传递函数的分子。

• aarray_like

模拟滤波器传递函数的分母。

• fsfloat

采样率，作为普通频率（例如赫兹）。此函数中不进行预弯曲。

返回：

• bndarray

转换后的数字滤波器传递函数的分子。

• andarray

转换后的数字滤波器传递函数的分母。

参见

lp2lp, lp2hp, lp2bp, lp2bs

bilinear_zpk

示例

>>> from scipy import signal
>>> import matplotlib.pyplot as plt
>>> import numpy as np 

>>> fs = 100
>>> bf = 2 * np.pi * np.array([7, 13])
>>> filts = signal.lti(*signal.butter(4, bf, btype='bandpass',
...                                   analog=True))
>>> filtz = signal.lti(*signal.bilinear(filts.num, filts.den, fs))
>>> wz, hz = signal.freqz(filtz.num, filtz.den)
>>> ws, hs = signal.freqs(filts.num, filts.den, worN=fs*wz) 

>>> plt.semilogx(wz*fs/(2*np.pi), 20*np.log10(np.abs(hz).clip(1e-15)),
...              label=r'$|H_z(e^{j \omega})|$')
>>> plt.semilogx(wz*fs/(2*np.pi), 20*np.log10(np.abs(hs).clip(1e-15)),
...              label=r'$|H(j \omega)|$')
>>> plt.legend()
>>> plt.xlabel('Frequency [Hz]')
>>> plt.ylabel('Magnitude [dB]')
>>> plt.grid(True) 

scipy.signal.bilinear_zpk

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.bilinear_zpk.html#scipy.signal.bilinear_zpk

scipy.signal.bilinear_zpk(z, p, k, fs)

使用双线性变换从模拟滤波器转换为数字 IIR 滤波器。

将一组模拟 s 平面的极点和零点转换为数字 z 平面，使用 Tustin 方法，用2*fs*(z-1) / (z+1)替换s，保持频率响应的形状。

Parameters:

zarray_like

模拟滤波器传递函数的零点。

parray_like

模拟滤波器传递函数的极点。

kfloat

System gain of the analog filter transfer function.

fsfloat

采样率，作为普通频率（例如赫兹）。此函数中不进行预变形。

返回：

zndarray

Zeros of the transformed digital filter transfer function.

pndarray

Poles of the transformed digital filter transfer function.

kfloat

转换后数字滤波器的系统增益。

See also

lp2lp_zpk, lp2hp_zpk, lp2bp_zpk, lp2bs_zpk

bilinear

Notes

New in version 1.1.0.

Examples

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt 

>>> fs = 100
>>> bf = 2 * np.pi * np.array([7, 13])
>>> filts = signal.lti(*signal.butter(4, bf, btype='bandpass', analog=True,
...                                   output='zpk'))
>>> filtz = signal.lti(*signal.bilinear_zpk(filts.zeros, filts.poles,
...                                         filts.gain, fs))
>>> wz, hz = signal.freqz_zpk(filtz.zeros, filtz.poles, filtz.gain)
>>> ws, hs = signal.freqs_zpk(filts.zeros, filts.poles, filts.gain,
...                           worN=fs*wz)
>>> plt.semilogx(wz*fs/(2*np.pi), 20*np.log10(np.abs(hz).clip(1e-15)),
...              label=r'$|H_z(e^{j \omega})|$')
>>> plt.semilogx(wz*fs/(2*np.pi), 20*np.log10(np.abs(hs).clip(1e-15)),
...              label=r'$|H(j \omega)|$')
>>> plt.legend()
>>> plt.xlabel('Frequency [Hz]')
>>> plt.ylabel('Magnitude [dB]')
>>> plt.grid(True) 

scipy.signal.findfreqs

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.findfreqs.html#scipy.signal.findfreqs

scipy.signal.findfreqs(num, den, N, kind='ba')

找到用于计算模拟滤波器响应的频率数组。

参数:

num, denarray_like, 1-D

滤波器或 LTI 系统传递函数的分子和分母的多项式系数，系数按从高到低的顺序排列。或者传递函数分子和分母的根（即零点和极点）。

Nint

要计算的数组长度。

kindstr {‘ba’, ‘zp’}, 可选

指定分子和分母是否由它们的多项式系数（‘ba’）或它们的根（‘zp’）指定。

返回:

w(N,) ndarray

一个频率的一维数组，对数间隔。

示例

找到跨越滤波器传递函数“有趣部分”的九个频率集合。

H(s) = s / (s² + 8s + 25)

>>> from scipy import signal
>>> signal.findfreqs([1, 0], [1, 8, 25], N=9)
array([  1.00000000e-02,   3.16227766e-02,   1.00000000e-01,
 3.16227766e-01,   1.00000000e+00,   3.16227766e+00,
 1.00000000e+01,   3.16227766e+01,   1.00000000e+02]) 

scipy.signal.firls

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.firls.html#scipy.signal.firls

scipy.signal.firls(numtaps, bands, desired, *, weight=None, nyq=<object object>, fs=None)

使用最小二乘误差最小化的 FIR 滤波器设计。

计算线性相位有限脉冲响应（FIR）滤波器的滤波器系数，其在最小二乘意义上对bands和desired中描述的期望频率响应的最佳逼近（即，在指定的带内加权均方误差的积分最小化）。

参数：

numtaps整数

FIR 滤波器的阶数。numtaps必须为奇数。

bands类数组

一个单调非递减的序列，其中包含 Hz 中的带边。所有元素必须非负且小于或等于nyq给定的奈奎斯特频率。带被指定为频率对，因此，如果使用 1D 数组，则其长度必须为偶数，例如np.array([0, 1, 2, 3, 4, 5])。或者，带可以作为大小为 nx2 的 2D 数组指定，其中 n 是带的数量，例如np.array([[0, 1], [2, 3], [4, 5]])。

desired类数组

与bands大小相同的序列，其中包含每个带的起始点和终点处的期望增益。

weight类数组，可选

在解最小二乘问题时，给每个带区域分配的相对权重。weight的大小必须是bands的一半。

nyq浮点数，可选，已弃用

这是奈奎斯特频率。bands中的每个频率必须介于 0 和nyq（包括）之间。默认为 1。

自 1.0.0 版起已弃用：firls关键字参数nyq已弃用，推荐使用fs，并将在 SciPy 1.14.0 中移除。

fs浮点数，可选

信号的采样频率。bands中的每个频率必须介于 0 和fs/2（包括）之间。默认为 2。

返回：

coeffsndarray

最优（在最小二乘意义上）FIR 滤波器的系数。

另请参见

firwin

firwin2

minimum_phase

remez

注释

此实现遵循[1]中给出的算法。如该文指出，最小二乘设计具有多个优点：

1. 最小二乘意义上的最优。
3. 简单的非迭代方法。
5. 通过解线性方程组获得一般解决方案。
7. 允许使用频率依赖的加权函数。

此函数构造一个 Type I 线性相位 FIR 滤波器，包含满足以下条件的奇数个coeffs，对于(n < numtaps)：

[coeffs(n) = coeffs(numtaps - 1 - n)]

系数的奇数和滤波器的对称性避免了在奈奎斯特频率和 0 频率处可能发生的边界条件（例如，对于 II 型、III 型或 IV 型变体）。

从版本 0.18 开始的新功能。

参考文献

[1]

Ivan Selesnick，最小二乘线性相位 FIR 滤波器设计。OpenStax CNX。2005 年 8 月 9 日。cnx.org/contents/eb1ecb35-03a9-4610-ba87-41cd771c95f2@7

示例

我们希望构建一个带通滤波器。请注意，在我们的阻带和通带之间的频率范围中的行为是未指定的，因此可能会根据我们滤波器的参数而超调：

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt
>>> fig, axs = plt.subplots(2)
>>> fs = 10.0  # Hz
>>> desired = (0, 0, 1, 1, 0, 0)
>>> for bi, bands in enumerate(((0, 1, 2, 3, 4, 5), (0, 1, 2, 4, 4.5, 5))):
...     fir_firls = signal.firls(73, bands, desired, fs=fs)
...     fir_remez = signal.remez(73, bands, desired[::2], fs=fs)
...     fir_firwin2 = signal.firwin2(73, bands, desired, fs=fs)
...     hs = list()
...     ax = axs[bi]
...     for fir in (fir_firls, fir_remez, fir_firwin2):
...         freq, response = signal.freqz(fir)
...         hs.append(ax.semilogy(0.5*fs*freq/np.pi, np.abs(response))[0])
...     for band, gains in zip(zip(bands[::2], bands[1::2]),
...                            zip(desired[::2], desired[1::2])):
...         ax.semilogy(band, np.maximum(gains, 1e-7), 'k--', linewidth=2)
...     if bi == 0:
...         ax.legend(hs, ('firls', 'remez', 'firwin2'),
...                   loc='lower center', frameon=False)
...     else:
...         ax.set_xlabel('Frequency (Hz)')
...     ax.grid(True)
...     ax.set(title='Band-pass %d-%d Hz' % bands[2:4], ylabel='Magnitude')
...
>>> fig.tight_layout()
>>> plt.show() 

scipy.signal.firwin

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.firwin.html#scipy.signal.firwin

scipy.signal.firwin(numtaps, cutoff, *, width=None, window='hamming', pass_zero=True, scale=True, nyq=<object object>, fs=None)

使用窗口方法设计 FIR 滤波器。

此函数计算有限冲激响应滤波器的系数。滤波器将具有线性相位；如果numtaps为奇数则为 Type I，如果numtaps为偶数则为 Type II。

Type II 滤波器在奈奎斯特频率处始终具有零响应，因此如果使用numtaps为偶数且其通带右端在奈奎斯特频率处的情况下调用 firwin，则会引发 ValueError 异常。

参数：

numtaps整数

滤波器的长度（系数数量，即滤波器阶数+1）。如果通带包含奈奎斯特频率，则numtaps必须为奇数。

cutoff浮点数或 1-D 数组

滤波器的截止频率（以与fs相同的单位表示）或截止频率数组（即带边缘）。在后一种情况下，cutoff中的频率应为正且单调增加，在 0 和fs/2之间不应包括值 0 和fs/2。

width浮点数或 None，可选

如果width不为 None，则假定其为过渡区域的大致宽度（以fs的相同单位表示），用于 Kaiser FIR 滤波器设计。在这种情况下，window参数将被忽略。

window字符串或字符串和参数值的元组，可选

所需使用的窗口。有关窗口和所需参数的列表，请参阅scipy.signal.get_window。

pass_zero，可选

如果为 True，则频率为 0 时的增益（即“直流增益”）为 1。如果为 False，则直流增益为 0。也可以是所需滤波器类型的字符串参数（相当于btype在 IIR 设计函数中的参数）。

从版本 1.3.0 开始支持字符串参数。

scale布尔值，可选

设置为 True 以使系数按比例缩放，以便频率响应在某个频率上完全为单位。该频率可以是：

• 如果第一个通带从 0 开始（即 pass_zero 为 True），则直流（DC）为 0。

• fs/2（奈奎斯特频率），如果第一个通带结束于fs/2（即滤波器是单通带高通滤波器）；否则为第一个通带的中心

nyq浮点数，可选，已弃用

这是奈奎斯特频率。cutoff中的每个频率必须介于 0 和nyq之间。默认为 1。

自版本 1.0.0 起不推荐使用：firwin 关键字参数nyq已弃用，推荐使用fs，并将在 SciPy 1.14.0 中移除。

fs浮点数，可选

信号的采样频率。cutoff中的每个频率必须介于 0 和fs/2之间。默认为 2。

返回：

h（numtaps，）ndarray

长度为numtaps的 FIR 滤波器系数。

引发：

ValueError

如果cutoff中的任何值小于等于 0 或大于等于fs/2，如果cutoff的值不是严格单调递增，或者numtaps是偶数但通带包含奈奎斯特频率。

参见

firwin2

firls

minimum_phase

remez

示例

低通从 0 到 f：

>>> from scipy import signal
>>> numtaps = 3
>>> f = 0.1
>>> signal.firwin(numtaps, f)
array([ 0.06799017,  0.86401967,  0.06799017]) 

使用特定的窗口函数：

>>> signal.firwin(numtaps, f, window='nuttall')
array([  3.56607041e-04,   9.99286786e-01,   3.56607041e-04]) 

高通（从 0 到 f）：

>>> signal.firwin(numtaps, f, pass_zero=False)
array([-0.00859313,  0.98281375, -0.00859313]) 

带通：

>>> f1, f2 = 0.1, 0.2
>>> signal.firwin(numtaps, [f1, f2], pass_zero=False)
array([ 0.06301614,  0.88770441,  0.06301614]) 

带阻：

>>> signal.firwin(numtaps, [f1, f2])
array([-0.00801395,  1.0160279 , -0.00801395]) 

多带通（通带为 [0, f1]，[f2, f3] 和 [f4, 1]）：

>>> f3, f4 = 0.3, 0.4
>>> signal.firwin(numtaps, [f1, f2, f3, f4])
array([-0.01376344,  1.02752689, -0.01376344]) 

多带通（通带为 [f1, f2] 和 [f3,f4]）：

>>> signal.firwin(numtaps, [f1, f2, f3, f4], pass_zero=False)
array([ 0.04890915,  0.91284326,  0.04890915]) 

scipy.signal.firwin2

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.firwin2.html#scipy.signal.firwin2

scipy.signal.firwin2(numtaps, freq, gain, *, nfreqs=None, window='hamming', nyq=<object object>, antisymmetric=False, fs=None)

使用窗口方法设计 FIR 滤波器。

根据给定的频率 freq 和相应的增益 gain，此函数构造具有线性相位和（近似）给定频率响应的 FIR 滤波器。

参数：

numtapsint

FIR 滤波器中的 taps 数。numtaps 必须小于 nfreqs。

freqarray_like, 1-D

频率采样点。通常为 0.0 到 1.0，其中 1.0 为奈奎斯特。奈奎斯特频率是 fs 的一半。 freq 中的值必须是非递减的。一个值可以重复一次以实现不连续性。 freq 中的第一个值必须为 0，最后一个值必须为 fs/2。值 0 和 fs/2 不得重复。

gainarray_like

频率采样点处的滤波器增益。根据滤波器类型应用某些增益值的约束条件，请参阅备注以获取详细信息。

nfreqsint, optional

用于构建滤波器的插值网格的大小。为了实现最有效的行为，这应该是一个 2 的幂加 1（例如 129, 257 等）。默认值为大于等于 numtaps 的最小 2 的幂加 1。nfreqs 必须大于 numtaps。

windowstring or (string, float) or float, or None, optional

要使用的窗口函数。默认值为“hamming”。参见 scipy.signal.get_window 获取可能值的完整列表。如果为 None，则不应用窗口函数。

nyqfloat, optional, deprecated

这是奈奎斯特频率。 freq 中的每个频率必须在 0 和 nyq 之间。默认值为 1。

自 1.0.0 版本起已弃用：firwin2 关键字参数 nyq 已弃用，改为使用 fs，将在 SciPy 1.14.0 中删除。

antisymmetricbool, optional

结果脉冲响应是否对称/反对称。更多细节请参见备注。

fsfloat, optional

信号的采样频率。 cutoff 中的每个频率必须在 0 和 fs/2 之间。默认值为 2。

返回：

tapsndarray

FIR 滤波器的滤波器系数，作为长度为 numtaps 的 1-D 数组。

参见

firls

firwin

minimum_phase

remez

备注

从给定的频率和增益集合中，在频率域中构造所需的响应。将逆 FFT 应用于所需的响应以创建相关的卷积核，并返回此卷积核的前 numtaps 系数，按 window 缩放。

FIR 滤波器将具有线性相位。滤波器的类型由 numtaps 的值和 antisymmetric 标志确定。有四种可能的组合：

• 即使 numtaps 为奇数，antisymmetric 为 False，生成类型 I 滤波器。
• 即使 numtaps 为偶数，antisymmetric 为 False，生成类型 II 滤波器。
• 即使 numtaps 为奇数，antisymmetric 为 True，生成类型 III 滤波器。
• 即使 numtaps 为偶数，antisymmetric 为 True，生成类型 IV 滤波器。

除了类型 I 滤波器外，所有滤波器的幅度响应都受以下约束的影响：

• 类型 II – 零频率处为 Nyquist 频率。
• 类型 III – 零频率和 Nyquist 频率处为零。
• 类型 IV – 零频率处为零。

新版本为 0.9.0。

参考文献

[1]

Oppenheim, A. V. 和 Schafer, R. W.，“Discrete-Time Signal Processing”，Prentice-Hall，Englewood Cliffs，New Jersey（1989）。（例如，参见第 7.4 节。）

[2]

Smith, Steven W.，“The Scientist and Engineer’s Guide to Digital Signal Processing”，第十七章。www.dspguide.com/ch17/1.htm

示例

一个低通 FIR 滤波器，其响应在 [0.0, 0.5] 上为 1，并且在 [0.5, 1.0] 上从 1 线性减少到 0：

>>> from scipy import signal
>>> taps = signal.firwin2(150, [0.0, 0.5, 1.0], [1.0, 1.0, 0.0])
>>> print(taps[72:78])
[-0.02286961 -0.06362756  0.57310236  0.57310236 -0.06362756 -0.02286961] 

scipy.signal.freqs

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.freqs.html#scipy.signal.freqs

scipy.signal.freqs(b, a, worN=200, plot=None)

计算模拟滤波器的频率响应。

给定模拟滤波器的 M 阶分子b和 N 阶分母a，计算其频率响应：

 b[0]*(jw)**M + b[1]*(jw)**(M-1) + ... + b[M]
H(w) = ----------------------------------------------
        a[0]*(jw)**N + a[1]*(jw)**(N-1) + ... + a[N] 

参数：

数组b

线性滤波器的分子。

数组a

线性滤波器的分母。

worN，可选

如果为 None，则在响应曲线的有趣部分周围的 200 个频率上计算（由极点零点位置决定）。如果是单个整数，则在那么多频率上计算。否则，计算在给定的角频率（例如，rad/s）处给出的响应worn。

plot可调用函数，可选

接受两个参数的可调用函数。如果给定，则将返回参数w和h传递给 plot。用于在freqs内部绘制频率响应。

返回：

数组w

计算h的角频率。

数组h

频率响应。

另请参见

freqz

计算数字滤波器的频率响应。

注意事项

使用 Matplotlib 的“plot”函数作为plot的可调用函数会产生意外的结果，这会绘制复数传递函数的实部，而不是幅度。尝试lambda w, h: plot(w, abs(h))。

示例

>>> from scipy.signal import freqs, iirfilter
>>> import numpy as np 

>>> b, a = iirfilter(4, [1, 10], 1, 60, analog=True, ftype='cheby1') 

>>> w, h = freqs(b, a, worN=np.logspace(-1, 2, 1000)) 

>>> import matplotlib.pyplot as plt
>>> plt.semilogx(w, 20 * np.log10(abs(h)))
>>> plt.xlabel('Frequency')
>>> plt.ylabel('Amplitude response [dB]')
>>> plt.grid(True)
>>> plt.show() 

scipy.signal.freqs_zpk

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.freqs_zpk.html#scipy.signal.freqs_zpk

scipy.signal.freqs_zpk(z, p, k, worN=200)

计算模拟滤波器的频率响应。

给定滤波器的零点z，极点p和增益k，计算其频率响应：

 (jw-z[0]) * (jw-z[1]) * ... * (jw-z[-1])
H(w) = k * ----------------------------------------
           (jw-p[0]) * (jw-p[1]) * ... * (jw-p[-1]) 

参数：

zarray_like

线性滤波器的零点

parray_like

线性滤波器的极点

kscalar

线性滤波器的增益

worN, 可选

如果为 None，则在响应曲线的有趣部分周围的 200 个频率处计算（由极点位置确定）。 如果为单个整数，则计算该数量的频率。 否则，计算在worN给定的角频率（例如，rad/s）处的响应。

返回：

wndarray

在计算h时使用的角频率。

hndarray

频率响应。

另请参见

freqs

计算 TF 形式模拟滤波器的频率响应

freqz

计算 TF 形式数字滤波器的频率响应

freqz_zpk

计算 ZPK 形式数字滤波器的频率响应

注意

新版 0.19.0 中新增。

示例

>>> import numpy as np
>>> from scipy.signal import freqs_zpk, iirfilter 

>>> z, p, k = iirfilter(4, [1, 10], 1, 60, analog=True, ftype='cheby1',
...                     output='zpk') 

>>> w, h = freqs_zpk(z, p, k, worN=np.logspace(-1, 2, 1000)) 

>>> import matplotlib.pyplot as plt
>>> plt.semilogx(w, 20 * np.log10(abs(h)))
>>> plt.xlabel('Frequency')
>>> plt.ylabel('Amplitude response [dB]')
>>> plt.grid(True)
>>> plt.show() 

scipy.signal.freqz

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.freqz.html#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 

参数：

b：array_like

线性滤波器的分子。如果b的维度大于 1，则假定系数存储在第一维度中，并且b.shape[1:]、a.shape[1:]和频率数组的形状必须兼容广播。

a：array_like

线性滤波器的分母。如果b的维度大于 1，则假定系数存储在第一维度中，并且b.shape[1:]、a.shape[1:]和频率数组的形状必须兼容广播。

worN：{None, int, array_like}，可选

如果是单个整数，则在那么多的频率上进行计算（默认值为 N=512）。这是以下便利的替代方法：

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

使用快速 FFT 计算的数字可以导致更快的计算（见备注）。

如果是 array_like，则在给定的频率上计算响应。这些频率与fs的单位相同。

whole：bool，可选

通常，频率从 0 到 Nyquist 频率 fs/2（单位圆的上半部分）计算。如果whole为 True，则从 0 到 fs 计算频率。如果 worN 是 array_like，则忽略。

plot：callable

一个接受两个参数的可调用函数。如果提供了，返回参数w和h将传递给绘图函数。用于在freqz内绘制频率响应。

fs：float，可选

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

1.2.0 版的新功能。

include_nyquist：bool，可选

如果whole为 False 且worN为整数，则将include_nyquist设置为 True 将包括最后一个频率（Nyquist 频率），否则将被忽略。

1.5.0 版的新功能。

返回：

w：ndarray

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

h：ndarray

频率响应，作为复数。

另请参阅

[freqz_zpk

sosfreqz

备注

当使用 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
>>> b = signal.firwin(80, 0.5, window=('kaiser', 8))
>>> w, h = signal.freqz(b) 

>>> import matplotlib.pyplot as plt
>>> fig, ax1 = plt.subplots()
>>> ax1.set_title('Digital filter frequency response') 

>>> ax1.plot(w, 20 * np.log10(abs(h)), 'b')
>>> ax1.set_ylabel('Amplitude [dB]', color='b')
>>> ax1.set_xlabel('Frequency [rad/sample]') 

>>> ax2 = ax1.twinx()
>>> angles = np.unwrap(np.angle(h))
>>> ax2.plot(w, angles, 'g')
>>> ax2.set_ylabel('Angle (radians)', color='g')
>>> ax2.grid(True)
>>> ax2.axis('tight')
>>> plt.show() 

广播示例

假设我们有两个 FIR 滤波器，它们的系数存储在形状为（2, 25）的数组的行中。为了演示，我们将使用随机数据：

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

为了一次调用freqz计算这两个滤波器的频率响应，我们必须传入b.T，因为freqz期望第一个轴包含系数。然后我们必须通过在长度为 1 的虚拟维度上扩展形状来允许与频率数组进行广播。也就是说，我们传入b.T[..., np.newaxis]，其形状为（25, 2, 1）：

>>> 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 维的。为了使其与频率广播兼容，在调用freqz时，我们通过在虚拟维度上扩展它来实现：

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

scipy.signal.freqz_zpk

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.freqz_zpk.html#scipy.signal.freqz_zpk

scipy.signal.freqz_zpk(z, p, k, worN=512, whole=False, fs=6.283185307179586)

计算 ZPK 形式数字滤波器的频率响应。

给定数字滤波器的零点、极点和增益，计算其频率响应：

(H(z)=k \prod_i (z - Z[i]) / \prod_j (z - P[j]))

其中(k)为增益，(Z)为零点，(P)为极点。

参数：

zarray_like

线性滤波器的零点

parray_like

线性滤波器的极点

k标量

线性滤波器的增益

worN，可选

如果是单个整数，则在该数量的频率上计算（默认值为 N=512）。

如果是 array_like，则计算给定频率处的响应。这些频率与fs具有相同的单位。

whole布尔值，可选

通常，频率从 0 到 Nyquist 频率 fs/2（单位圆的上半部分）计算。如果whole为 True，则从 0 到 fs 计算频率。如果w为 array_like，则忽略。

fs浮点数，可选

数字系统的采样频率。默认为 2pi 弧度/样本（因此w*从 0 到 pi）。

新版本为 1.2.0。

返回：

wndarray

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

hndarray

作为复数的频率响应。

另请参见

freqs

计算 TF 形式模拟滤波器的频率响应

freqs_zpk

计算 ZPK 形式模拟滤波器的频率响应

freqz

计算 TF 形式数字滤波器的频率响应

笔记

新版本为 0.19.0。

示例

在采样率为 1000 Hz 的系统中，设计一个 4 阶数字 Butterworth 滤波器，截止频率为 100 Hz，并绘制其频率响应：

>>> import numpy as np
>>> from scipy import signal
>>> z, p, k = signal.butter(4, 100, output='zpk', fs=1000)
>>> w, h = signal.freqz_zpk(z, p, k, fs=1000) 

>>> import matplotlib.pyplot as plt
>>> fig = plt.figure()
>>> ax1 = fig.add_subplot(1, 1, 1)
>>> ax1.set_title('Digital filter frequency response') 

>>> ax1.plot(w, 20 * np.log10(abs(h)), 'b')
>>> ax1.set_ylabel('Amplitude [dB]', color='b')
>>> ax1.set_xlabel('Frequency [Hz]')
>>> ax1.grid(True) 

>>> ax2 = ax1.twinx()
>>> angles = np.unwrap(np.angle(h))
>>> ax2.plot(w, angles, 'g')
>>> ax2.set_ylabel('Angle [radians]', color='g') 

>>> plt.axis('tight')
>>> plt.show() 

scipy.signal.sosfreqz

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.sosfreqz.html#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 个二阶段传递函数的分子和分母。

参数：

sos 类似数组

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

worN，可选

如果是单个整数，则计算那么多频率（默认为 N=512）。使用 FFT 计算快速的数字可以导致更快的计算（见freqz的注意事项）。

如果是类似数组，则在给定频率处计算响应（必须是 1-D）。这些单位与fs相同。

whole 布尔值，可选

通常，频率从 0 到 Nyquist 频率 fs/2 计算（单位圆的上半部分）。如果whole为 True，则从 0 到 fs 计算频率。

fs 浮点数，可选

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

新版本 1.2.0 中。

返回：

w ndarray

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

h ndarray

频率响应，作为复数。

另请参见

freqz，sosfilt

注意

新版本 0.19.0 中。

示例

设计一个 15 阶带通滤波器的 SOS 格式。

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

在 1500 点处从 DC 到 Nyquist 计算频率响应。

>>> 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() 

scipy.signal.gammatone

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.gammatone.html#scipy.signal.gammatone

scipy.signal.gammatone(freq, ftype, order=None, numtaps=None, fs=None)

Gammatone 滤波器设计。

此函数计算 FIR 或 IIR Gammatone 数字滤波器的系数 [1]。

参数：

freqfloat

滤波器的中心频率（与 fs 相同的单位表示）。

ftype

函数生成的滤波器类型。如果是 ‘fir’，函数将生成一个 N 阶 FIR Gammatone 滤波器。如果是 ‘iir’，函数将生成一个 8 阶数字 IIR 滤波器，建模为 4 阶 Gammatone 滤波器。

orderint, optional

滤波器的阶数。仅在 ftype='fir' 时使用。默认为 4，用于模拟人类听觉系统。必须介于 0 和 24 之间。

numtapsint, optional

滤波器的长度。仅在 ftype='fir' 时使用。默认为 fs*0.015（如果 fs 大于 1000），15（如果 fs 小于或等于 1000）。

fsfloat, optional

信号的采样频率。freq 必须介于 0 和 fs/2 之间。默认为 2。

返回：

b, andarray, ndarray

滤波器的分子 (b) 和分母 (a) 多项式。

Raises:

ValueError

如果 freq 小于或等于 0 或大于或等于 fs/2，如果 ftype 不是 ‘fir’ 或 ‘iir’，如果 order 在 ftype='fir' 时小于或等于 0 或大于 24

参见

firwin

iirfilter

参考文献

[1]

Slaney, Malcolm, “An Efficient Implementation of the Patterson-Holdsworth Auditory Filter Bank”, Apple Computer Technical Report 35, 1993, pp.3-8, 34-39.

示例

以 440 Hz 为中心的 16 采样 4 阶 FIR Gammatone 滤波器

>>> from scipy import signal
>>> signal.gammatone(440, 'fir', numtaps=16, fs=16000)
(array([ 0.00000000e+00,  2.22196719e-07,  1.64942101e-06,  4.99298227e-06,
 1.01993969e-05,  1.63125770e-05,  2.14648940e-05,  2.29947263e-05,
 1.76776931e-05,  2.04980537e-06, -2.72062858e-05, -7.28455299e-05,
 -1.36651076e-04, -2.19066855e-04, -3.18905076e-04, -4.33156712e-04]),
 [1.0]) 

以 440 Hz 为中心的 IIR Gammatone 滤波器

>>> import matplotlib.pyplot as plt
>>> import numpy as np 

>>> b, a = signal.gammatone(440, 'iir', fs=16000)
>>> w, h = signal.freqz(b, a)
>>> plt.plot(w / ((2 * np.pi) / 16000), 20 * np.log10(abs(h)))
>>> plt.xscale('log')
>>> plt.title('Gammatone filter frequency response')
>>> plt.xlabel('Frequency')
>>> plt.ylabel('Amplitude [dB]')
>>> plt.margins(0, 0.1)
>>> plt.grid(which='both', axis='both')
>>> plt.axvline(440, color='green') # cutoff frequency
>>> plt.show() 

scipy.signal.group_delay

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.group_delay.html#scipy.signal.group_delay

scipy.signal.group_delay(system, w=512, whole=False, fs=6.283185307179586)

计算数字滤波器的组延迟。

组延迟测量信号各谱成分幅度包络被滤波器延迟多少个样本。形式上定义为连续（展开）相位的导数：

 d        jw
D(w) = - -- arg H(e)
         dw 

参数：

system数组对（b, a）

滤波器传输函数的分子和分母系数。

w，可选

如果是单个整数，则在那么多的频率上进行计算（默认为 N=512）。

如果是数组形式，则计算给定频率下的延迟。这些频率与fs单位相同。

whole布尔值，可选

通常，频率从 0 到奈奎斯特频率 fs/2（单位圆的上半部分）计算。如果whole为 True，则从 0 到 fs 计算频率。如果 w 是数组形式，则忽略。

fs浮点数，可选

数字系统的采样频率。默认为 2*pi 弧度/样本（所以 w 在 0 到 pi 之间）。

新版功能于版本 1.2.0 中添加。

返回：

wndarray

计算组延迟的频率，单位与fs相同。默认情况下，w被归一化到范围 0, pi)（弧度/样本）。

gdndarray

组延迟。

另请参阅

[freqz

数字滤波器的频率响应

注释

MATLAB 中的类似函数称为grpdelay。

如果数字系统的传输函数(H(z))在单位圆上有零点或极点，则在相应频率下的组延迟是未定义的。当出现这种情况时，会发出警告，并将组延迟设置为这些频率上的 0。

关于组延迟的数值计算的详细信息，请参考[1]。

新版功能于版本 0.16.0 中添加。

参考文献

[1]

Richard G. Lyons，《理解数字信号处理，第 3 版》，第 830 页。

示例

>>> from scipy import signal
>>> b, a = signal.iirdesign(0.1, 0.3, 5, 50, ftype='cheby1')
>>> w, gd = signal.group_delay((b, a)) 

>>> import matplotlib.pyplot as plt
>>> plt.title('Digital filter group delay')
>>> plt.plot(w, gd)
>>> plt.ylabel('Group delay [samples]')
>>> plt.xlabel('Frequency [rad/sample]')
>>> plt.show() 

scipy.signal.iirdesign

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.iirdesign.html#scipy.signal.iirdesign

scipy.signal.iirdesign(wp, ws, gpass, gstop, analog=False, ftype='ellip', output='ba', fs=None)

完整的 IIR 数字和模拟滤波器设计。

根据给定的基本类型的通带和阻带频率及增益构造模拟或数字 IIR 滤波器的最小阶数。以分子、分母（‘ba’）、极点-零点（‘zpk’）或二阶段（‘sos’）形式返回输出。

参数：

wp, wsfloat 或 array like, 形状 (2,)

通带和阻带边缘频率。可能的取值为标量（适用于低通和高通滤波器）或范围（适用于带通和带阻滤波器）。对于数字滤波器，这些频率与fs（采样频率）的单位相同。默认情况下，fs是每个样本的 2 个半周期，因此这些频率被归一化为 0 到 1，其中 1 是奈奎斯特频率。例如：

• 低通：wp = 0.2，ws = 0.3
• 高通：wp = 0.3，ws = 0.2
• 带通：wp = [0.2, 0.5]，ws = [0.1, 0.6]
• 带阻：wp = [0.1, 0.6]，ws = [0.2, 0.5]

对于模拟滤波器，wp 和 ws 是角频率（例如，rad/s）。注意，对于带通和带阻滤波器，通带必须严格位于阻带内，反之亦然。

gpassfloat

通带中的最大损失（dB）。

gstopfloat

在阻带中的最小衰减（dB）。

analogbool, 可选

当为 True 时，返回模拟滤波器，否则返回数字滤波器。

ftypestr, 可选

要设计的 IIR 滤波器类型：

• Butterworth：‘butter’
• Chebyshev I：‘cheby1’
• Chebyshev II：‘cheby2’
• Cauer/elliptic：‘ellip’

output, 可选

输出的滤波器形式：

• 推荐的二阶段形式：‘sos’
• 分子/分母（默认）：‘ba’
• 极点-零点：‘zpk’

一般推荐使用二阶段形式（‘sos’），因为推断分子/分母形式（‘ba’）的系数会受到数值不稳定性的影响。出于向后兼容性的考虑，默认形式是分子/分母形式（‘ba’），其中‘b’ 和 ‘a’ 分别是系数的常用名称。

注意：有时使用二阶段形式（‘sos’）会伴随额外的计算成本：因此建议对数据密集型应用进行探索，也要考虑分子/分母形式（‘ba’）。

fsfloat, 可选

数字系统的采样频率。

新版本 1.2.0 中的新增功能。

返回：

b, andarray, ndarray

IIR 滤波器的分子（b）和分母（a）多项式。仅在output='ba'时返回。

z, p, kndarray, ndarray, float

IIR 滤波器传递函数的零点、极点和系统增益。仅在output='zpk'时返回。

sosndarray

IIR 滤波器的二阶段表示。仅在output='sos'时返回。

另请参见

butter

使用阶数和临界点设计滤波器

cheby1, cheby2, ellip, bessel

buttord

根据通带和阻带规格找到阶数和关键点

cheb1ord, cheb2ord, ellipord

iirfilter

使用阶数和关键频率进行一般滤波器设计

笔记

'sos'输出参数是在 0.16.0 版本中添加的。

示例

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt
>>> import matplotlib.ticker 

>>> wp = 0.2
>>> ws = 0.3
>>> gpass = 1
>>> gstop = 40 

>>> system = signal.iirdesign(wp, ws, gpass, gstop)
>>> w, h = signal.freqz(*system) 

>>> fig, ax1 = plt.subplots()
>>> ax1.set_title('Digital filter frequency response')
>>> ax1.plot(w, 20 * np.log10(abs(h)), 'b')
>>> ax1.set_ylabel('Amplitude [dB]', color='b')
>>> ax1.set_xlabel('Frequency [rad/sample]')
>>> ax1.grid(True)
>>> ax1.set_ylim([-120, 20])
>>> ax2 = ax1.twinx()
>>> angles = np.unwrap(np.angle(h))
>>> ax2.plot(w, angles, 'g')
>>> ax2.set_ylabel('Angle (radians)', color='g')
>>> ax2.grid(True)
>>> ax2.axis('tight')
>>> ax2.set_ylim([-6, 1])
>>> nticks = 8
>>> ax1.yaxis.set_major_locator(matplotlib.ticker.LinearLocator(nticks))
>>> ax2.yaxis.set_major_locator(matplotlib.ticker.LinearLocator(nticks)) 

scipy.signal.iirfilter

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.iirfilter.html#scipy.signal.iirfilter

scipy.signal.iirfilter(N, Wn, rp=None, rs=None, btype='band', analog=False, ftype='butter', output='ba', fs=None)

给定阶数和临界点，设计 IIR 数字和模拟滤波器。

设计一个 N 阶数字或模拟滤波器，并返回滤波器系数。

参数：

N 整数

滤波器的阶数。

Wn 类似数组

给出临界频率的标量或长度为 2 的序列。

对于数字滤波器，Wn的单位与fs相同。默认情况下，fs为每个采样 2 个半周期，因此这些值从 0 到 1 进行归一化，其中 1 为奈奎斯特频率。(Wn因此是半周期/每个样本。)

对于模拟滤波器，Wn 是一个角频率（例如，rad/s）。

当 Wn 是长度为 2 的序列时，Wn[0]必须小于Wn[1]。

rp 浮点数, 可选

对于 Chebyshev 和 elliptic 滤波器，提供通带中的最大波纹。（dB）

rs 浮点数, 可选

对于 Chebyshev 和 elliptic 滤波器，提供阻带中的最小衰减。（dB）

btype {'bandpass', 'lowpass', 'highpass', 'bandstop'}, 可选

滤波器类型。默认为‘bandpass’。

analog 布尔值, 可选

当为 True 时，返回模拟滤波器，否则返回数字滤波器。

ftype 字符串, 可选

要设计的 IIR 滤波器类型：

• Butterworth：‘butter’
• Chebyshev I：‘cheby1’
• Chebyshev II：‘cheby2’
• Cauer/elliptic：‘ellip’
• Bessel/Thomson：‘bessel’

output {'ba', 'zpk', 'sos'}, 可选

输出的过滤器形式：

• 二阶段形式（推荐）：‘sos’
• 默认的分子/分母形式：‘ba’
• 极点-零点形式：‘zpk’

一般推荐使用二阶段形式（‘sos’），因为推断出分子/分母形式（‘ba’）的系数会受到数值不稳定性的影响。出于向后兼容性的考虑，默认形式为分子/分母形式（‘ba’），其中‘b’和‘a’分别指代所用系数的常用名称。

注意：使用二阶段形式（‘sos’）有时会伴随额外的计算成本：因此，建议对于数据密集的用例也研究分子/分母形式（‘ba’）。

fs 浮点数, 可选

数字系统的采样频率。

1.2.0 版中新增。

返回：

b, a 数组, 数组

IIR 滤波器的分子（b）和分母（a）多项式。仅在output='ba'时返回。

z, p, k 数组, 数组, 浮点数

IIR 滤波器传递函数的零点、极点和系统增益。仅在output='zpk'时返回。

sos 数组

IIR 滤波器的二阶段形式表示。仅在output='sos'时返回。

另请参阅

butter

使用阶数和临界点进行滤波器设计。

cheby1, cheby2, ellip, bessel

buttord

从通带和阻带规范中找到阶数和临界点。

cheb1ord, cheb2ord, ellipord

iirdesign

使用通带和阻带规范进行通用滤波器设计。

注意

'sos'输出参数在 0.16.0 版本中被添加。

示例

生成一个从 50 Hz 到 200 Hz 的 17 阶 Chebyshev II 模拟带通滤波器，并绘制频率响应图：

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt 

>>> b, a = signal.iirfilter(17, [2*np.pi*50, 2*np.pi*200], rs=60,
...                         btype='band', analog=True, ftype='cheby2')
>>> w, h = signal.freqs(b, a, 1000)
>>> fig = plt.figure()
>>> ax = fig.add_subplot(1, 1, 1)
>>> ax.semilogx(w / (2*np.pi), 20 * np.log10(np.maximum(abs(h), 1e-5)))
>>> ax.set_title('Chebyshev Type II bandpass frequency response')
>>> ax.set_xlabel('Frequency [Hz]')
>>> ax.set_ylabel('Amplitude [dB]')
>>> ax.axis((10, 1000, -100, 10))
>>> ax.grid(which='both', axis='both')
>>> plt.show() 

在采样率为 2000 Hz 的系统中创建具有相同特性的数字滤波器，并绘制频率响应图。（需要使用二阶段实现以确保这一阶数的滤波器稳定性）：

>>> sos = signal.iirfilter(17, [50, 200], rs=60, btype='band',
...                        analog=False, ftype='cheby2', fs=2000,
...                        output='sos')
>>> w, h = signal.sosfreqz(sos, 2000, fs=2000)
>>> fig = plt.figure()
>>> ax = fig.add_subplot(1, 1, 1)
>>> ax.semilogx(w, 20 * np.log10(np.maximum(abs(h), 1e-5)))
>>> ax.set_title('Chebyshev Type II bandpass frequency response')
>>> ax.set_xlabel('Frequency [Hz]')
>>> ax.set_ylabel('Amplitude [dB]')
>>> ax.axis((10, 1000, -100, 10))
>>> ax.grid(which='both', axis='both')
>>> plt.show() 

scipy.signal.kaiser_atten

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.kaiser_atten.html#scipy.signal.kaiser_atten

scipy.signal.kaiser_atten(numtaps, width)

计算 Kaiser FIR 滤波器的衰减。

给定抽头数N和过渡宽度width，使用 Kaiser 公式计算衰减a，如下所示：

a = 2.285 * (N - 1) * pi * width + 7.95

Parameters:

numtapsint

FIR 滤波器的抽头数量。

widthfloat

滤波器在通带和阻带之间（或一般来说，在任何不连续处）的过渡区域的期望宽度，以奈奎斯特频率的分数形式表示。

Returns:

afloat

波纹的衰减，单位为 dB。

另请参见

kaiserord, kaiser_beta

Examples

假设我们希望使用 Kaiser 窗口方法设计一个 FIR 滤波器，该滤波器有 211 个抽头，并且在采样频率为 480 Hz 的信号中具有 9 Hz 的过渡宽度。以奈奎斯特频率的分数形式表示，宽度为 9/(0.5*480) = 0.0375. 按照以下公式计算近似衰减（以 dB 为单位）：

>>> from scipy.signal import kaiser_atten
>>> kaiser_atten(211, 0.0375)
64.48099630593983 

scipy.signal.kaiser_beta

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.kaiser_beta.html#scipy.signal.kaiser_beta

scipy.signal.kaiser_beta(a)

计算 Kaiser 参数 beta，给定阻带衰减 a。

参数：

afloat

在 dB 中所需的阻带衰减和通带中的最大波动。这应该是一个 正数。

返回值：

betafloat

用于计算 Kaiser 窗口公式中的 beta 参数。

参考文献

Oppenheim, Schafer, “Discrete-Time Signal Processing”, p.475-476.

示例

假设我们想设计一个低通滤波器，在阻带中具有 65 dB 的衰减。通过 kaiser_beta(65) 计算的 Kaiser 窗口参数如下：

>>> from scipy.signal import kaiser_beta
>>> kaiser_beta(65)
6.20426 

scipy.signal.kaiserord

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.kaiserord.html#scipy.signal.kaiserord

scipy.signal.kaiserord(ripple, width)

确定 Kaiser 窗口方法的滤波器窗口参数。

此函数返回的参数通常用于使用窗口法创建有限冲激响应滤波器，可以使用firwin或firwin2。

参数：

ripple浮点数

滤波器的频率响应的幅度与所需滤波器的频率响应的偏差的上限（不包括任何过渡区间中的频率）。也就是说，如果 w 是以 Nyquist 频率的分数表示的频率，则 A(w)是滤波器的实际频率响应，D(w)是期望的频率响应，则设计要求是：

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

对于 0 <= w <= 1 且 w 不在过渡区间内。

width浮点数

过渡区域的宽度，标准化为对应于每个采样π弧度。也就是说，频率表示为 Nyquist 频率的分数。

返回：

numtaps整数

Kaiser 窗口的长度。

beta浮点数

Kaiser 窗口的 beta 参数。

另见

kaiser_beta, kaiser_atten

注意事项

有几种方法可以获取 Kaiser 窗口：

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

• signal.get_window(beta, numtaps)

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

Kaiser 发现的经验方程式被使用。

参考文献

Oppenheim, Schafer, “离散时间信号处理”, pp.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。这对于凯泽窗方法并不罕见。

>>> 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() 

scipy.signal.minimum_phase

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.minimum_phase.html#scipy.signal.minimum_phase

scipy.signal.minimum_phase(h, method='homomorphic', n_fft=None)

将线性相位 FIR 滤波器转换为最小相位

参数：

h数组

线性相位 FIR 滤波器系数。

method

使用的方法：

‘homomorphic’（默认）

此方法[4] [5] 最适用于具有奇数抽头数的滤波器，并且生成的最小相位滤波器的幅度响应近似于原始滤波器幅度响应的平方根。

‘hilbert’

此方法[1]设计用于等波纹滤波器（例如来自remez的滤波器），具有单位或零增益区域。

n_fft整数

FFT 使用的点数。应至少比信号长度大几倍（见注释）。

返回：

h_minimum数组

滤波器的最小相位版本，长度为(length(h) + 1) // 2。

另请参阅

firwin

firwin2

remez

注释

希尔伯特[1]或同态[4] [5]方法都需要选择 FFT 长度以估算滤波器的复合倒谱。

在希尔伯特方法中，偏离理想频谱的epsilon与阻带零点数n_stop和 FFT 长度n_fft有关：

epsilon = 2. * n_stop / n_fft 

例如，有 100 个阻带零点和 FFT 长度为 2048 时，epsilon = 0.0976。如果我们保守地假设阻带零点数比滤波器长度少一个，我们可以将 FFT 长度取为满足epsilon = 0.01的下一个 2 的幂：

n_fft = 2 ** int(np.ceil(np.log2(2 * (len(h) - 1) / 0.01))) 

在n_fft=None时使用的值，此方法对希尔伯特和同态方法都给出了合理的结果。

还存在其他方法创建最小相位滤波器，包括零反转[2]和频谱因子分解[3] [4]。更多信息请参见：

dspguru.com/dsp/howtos/how-to-design-minimum-phase-fir-filters

参考文献

[1] (1,2)

N. Damera-Venkata 和 B. L. Evans，“实数和复数最小相位数字 FIR 滤波器的最优设计”，声学、语音和信号处理，1999 年国际会议记录，凤凰城，AZ，1999，第 1145-1148 卷 3。DOI:10.1109/ICASSP.1999.756179

[2]

X. Chen 和 T. W. Parks，《通过直接因式分解设计最优最小相位 FIR 滤波器》，《信号处理》，第 10 卷，第 4 期，pp. 369-383，1986 年 6 月。

[3]

T. Saramaki，《有限冲激响应滤波器设计》，《数字信号处理手册》，第四章，纽约：Wiley-Interscience，1993 年。

[4] (1,2,3)

J. S. Lim，《信号处理的高级主题》。新泽西州恩格尔伍德克利夫斯：普林斯顿大厅，1988 年。

[5] (1,2)

A. V. Oppenheim, R. W. Schafer, 和 J. R. Buck，《离散时间信号处理》，第二版。新泽西州，上班顶部：普林斯顿大厅，1999 年。

例子

创建一个最优的线性相位滤波器，然后将其转换为最小相位：

>>> import numpy as np
>>> from scipy.signal import remez, minimum_phase, freqz, group_delay
>>> import matplotlib.pyplot as plt
>>> freq = [0, 0.2, 0.3, 1.0]
>>> desired = [1, 0]
>>> h_linear = remez(151, freq, desired, fs=2.) 

将其转换为最小相位：

>>> h_min_hom = minimum_phase(h_linear, method='homomorphic')
>>> h_min_hil = minimum_phase(h_linear, method='hilbert') 

比较这三个滤波器：

>>> fig, axs = plt.subplots(4, figsize=(4, 8))
>>> for h, style, color in zip((h_linear, h_min_hom, h_min_hil),
...                            ('-', '-', '--'), ('k', 'r', 'c')):
...     w, H = freqz(h)
...     w, gd = group_delay((h, 1))
...     w /= np.pi
...     axs[0].plot(h, color=color, linestyle=style)
...     axs[1].plot(w, np.abs(H), color=color, linestyle=style)
...     axs[2].plot(w, 20 * np.log10(np.abs(H)), color=color, linestyle=style)
...     axs[3].plot(w, gd, color=color, linestyle=style)
>>> for ax in axs:
...     ax.grid(True, color='0.5')
...     ax.fill_between(freq[1:3], *ax.get_ylim(), color='#ffeeaa', zorder=1)
>>> axs[0].set(xlim=[0, len(h_linear) - 1], ylabel='Amplitude', xlabel='Samples')
>>> axs[1].legend(['Linear', 'Min-Hom', 'Min-Hil'], title='Phase')
>>> for ax, ylim in zip(axs[1:], ([0, 1.1], [-150, 10], [-60, 60])):
...     ax.set(xlim=[0, 1], ylim=ylim, xlabel='Frequency')
>>> axs[1].set(ylabel='Magnitude')
>>> axs[2].set(ylabel='Magnitude (dB)')
>>> axs[3].set(ylabel='Group delay')
>>> plt.tight_layout() 

scipy.signal.savgol_coeffs

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.savgol_coeffs.html#scipy.signal.savgol_coeffs

scipy.signal.savgol_coeffs(window_length, polyorder, deriv=0, delta=1.0, pos=None, use='conv')

计算 1-D Savitzky-Golay FIR 滤波器的系数。

参数：

window_length：整数

滤波器窗口的长度（即系数的数量）。

polyorder：整数

用于拟合样本的多项式的顺序。polyorder必须小于window_length。

deriv：整数，可选

要计算的导数阶数。这必须是非负整数。默认值为 0，表示在不进行微分的情况下过滤数据。

delta：浮点数，可选

要应用滤波器的样本的间距。仅当 deriv > 0 时使用。

pos：整数或 None，可选

如果 pos 不为 None，则指定窗口内的评估位置。默认值为窗口的中间。

use：字符串，可选

‘conv’或‘dot’。此参数选择系数的顺序。默认值为‘conv’，表示系数按卷积使用的顺序排列。使用‘dot’时，顺序反转，因此通过将系数与数据集点乘来应用滤波器。

返回：

coeffs：1-D ndarray

滤波器系数。

参见

savgol_filter

注意事项

新版本 0.14.0 中引入。

参考文献

A. Savitzky, M. J. E. Golay, 简化最小二乘法的平滑和微分数据处理。分析化学，1964 年，36（8），第 1627-1639 页。罗建文，应奎，白静。2005 年。用于偶数数据的 Savitzky-Golay 平滑和微分滤波器。信号处理。85，7（2005 年 7 月），第 1429-1434 页。

示例

>>> import numpy as np
>>> from scipy.signal import savgol_coeffs
>>> savgol_coeffs(5, 2)
array([-0.08571429,  0.34285714,  0.48571429,  0.34285714, -0.08571429])
>>> savgol_coeffs(5, 2, deriv=1)
array([ 2.00000000e-01,  1.00000000e-01,  2.07548111e-16, -1.00000000e-01,
 -2.00000000e-01]) 

注意，use='dot'仅简单地反转系数。

>>> savgol_coeffs(5, 2, pos=3)
array([ 0.25714286,  0.37142857,  0.34285714,  0.17142857, -0.14285714])
>>> savgol_coeffs(5, 2, pos=3, use='dot')
array([-0.14285714,  0.17142857,  0.34285714,  0.37142857,  0.25714286])
>>> savgol_coeffs(4, 2, pos=3, deriv=1, use='dot')
array([0.45,  -0.85,  -0.65,  1.05]) 

x包含从抛物线 x = t**2 采样的数据，采样点为 t = -1, 0, 1, 2, 3。c保存了在最后一个位置计算导数的系数。当与x点乘时，结果应为 6。

>>> x = np.array([1, 0, 1, 4, 9])
>>> c = savgol_coeffs(5, 2, pos=4, deriv=1, use='dot')
>>> c.dot(x)
6.0 

scipy.signal.remez

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.remez.html#scipy.signal.remez

scipy.signal.remez(numtaps, bands, desired, *, weight=None, Hz=<object object>, type='bandpass', maxiter=25, grid_density=16, fs=None)

使用 Remez 交换算法计算极小值最优滤波器。

使用 Remez 交换算法在指定频段中最小化所需增益与实际增益之间的最大误差，计算有限脉冲响应（FIR）滤波器的滤波器系数。

参数：

numtaps整数

滤波器中所需的阶数。阶数是滤波器中的项数，或者是滤波器阶数加一。

bands数组型

包含带边缘的单调序列。所有元素必须为非负且小于由 fs 给出的采样频率的一半。

desired数组型

包含每个指定频段中所需增益的带的一半大小的序列。

weight数组型，可选

给每个带区域赋予的相对权重。weight 的长度必须是 bands 长度的一半。

Hz标量，可选，已弃用

信号的采样频率（单位：赫兹）。默认为 1。

自 1.0.0 版起弃用：remez关键字参数Hz，将被fs取代，并将在 SciPy 1.14.0 中删除。

type，可选

滤波器类型：

• ‘bandpass’：带通响应。这是默认设置。
• ‘differentiator’：频率比例响应带。
• ‘hilbert’滤波器，具有奇对称性，即类型 III
• （对于偶数阶）或类型 IV（对于奇数阶）线性相位滤波器。

maxiter整数，可选

算法的最大迭代次数。默认为 25。

grid_density整数，可选

网格密度。在remez中使用的密集网格大小为(numtaps + 1) * grid_density。默认为 16。

fs浮点数，可选

信号的采样频率。默认为 1。

返回：

outndarray

一个秩为 1 的数组，包含最优（在最小最大意义上）滤波器的系数。

另见

firls

firwin

firwin2

minimum_phase

参考文献

[1]

J. H. McClellan 和 T. W. Parks，“一种统一的最优 FIR 线性相位数字滤波器设计方法”，IEEE Trans. Circuit Theory, vol. CT-20, pp. 697-701, 1973 年。

[2]

J. H. McClellan, T. W. Parks 和 L. R. Rabiner，“用于设计最优 FIR 线性相位数字滤波器的计算机程序”，IEEE Trans. Audio Electroacoust., vol. AU-21, pp. 506-525, 1973 年。

示例

在这些示例中，remez 用于设计低通、高通、带通和带阻滤波器。定义每个滤波器的参数包括滤波器阶数、频带边界、边界过渡宽度、每个频带中的期望增益以及采样频率。

在所有示例中，我们将使用 22050 Hz 的采样频率。在每个示例中，每个频带中的期望增益为 0（阻带）或 1（通带）。

freqz 用于计算每个滤波器的频率响应，下面定义的实用函数 plot_response 用于绘制响应。

>>> import numpy as np
>>> from scipy import signal
>>> import matplotlib.pyplot as plt 

>>> fs = 22050   # Sample rate, Hz 

>>> def plot_response(w, h, title):
...     "Utility function to plot response functions"
...     fig = plt.figure()
...     ax = fig.add_subplot(111)
...     ax.plot(w, 20*np.log10(np.abs(h)))
...     ax.set_ylim(-40, 5)
...     ax.grid(True)
...     ax.set_xlabel('Frequency (Hz)')
...     ax.set_ylabel('Gain (dB)')
...     ax.set_title(title) 

第一个示例是一个低通滤波器，截止频率为 8 kHz。滤波器长度为 325，从通带到阻带的过渡宽度为 100 Hz。

>>> cutoff = 8000.0    # Desired cutoff frequency, Hz
>>> trans_width = 100  # Width of transition from pass to stop, Hz
>>> numtaps = 325      # Size of the FIR filter.
>>> taps = signal.remez(numtaps, [0, cutoff, cutoff + trans_width, 0.5*fs],
...                     [1, 0], fs=fs)
>>> w, h = signal.freqz(taps, [1], worN=2000, fs=fs)
>>> plot_response(w, h, "Low-pass Filter")
>>> plt.show() 

此示例展示了一个高通滤波器：

>>> cutoff = 2000.0    # Desired cutoff frequency, Hz
>>> trans_width = 250  # Width of transition from pass to stop, Hz
>>> numtaps = 125      # Size of the FIR filter.
>>> taps = signal.remez(numtaps, [0, cutoff - trans_width, cutoff, 0.5*fs],
...                     [0, 1], fs=fs)
>>> w, h = signal.freqz(taps, [1], worN=2000, fs=fs)
>>> plot_response(w, h, "High-pass Filter")
>>> plt.show() 

此示例展示了一个带通滤波器，通带从 2 kHz 到 5 kHz。过渡宽度为 260 Hz，滤波器长度为 63，比其他示例中的要小：

>>> band = [2000, 5000]  # Desired pass band, Hz
>>> trans_width = 260    # Width of transition from pass to stop, Hz
>>> numtaps = 63         # Size of the FIR filter.
>>> edges = [0, band[0] - trans_width, band[0], band[1],
...          band[1] + trans_width, 0.5*fs]
>>> taps = signal.remez(numtaps, edges, [0, 1, 0], fs=fs)
>>> w, h = signal.freqz(taps, [1], worN=2000, fs=fs)
>>> plot_response(w, h, "Band-pass Filter")
>>> plt.show() 

低阶数导致更高的波动和更缓的过渡。

接下来的示例展示了一个带阻滤波器。

>>> band = [6000, 8000]  # Desired stop band, Hz
>>> trans_width = 200    # Width of transition from pass to stop, Hz
>>> numtaps = 175        # Size of the FIR filter.
>>> edges = [0, band[0] - trans_width, band[0], band[1],
...          band[1] + trans_width, 0.5*fs]
>>> taps = signal.remez(numtaps, edges, [1, 0, 1], fs=fs)
>>> w, h = signal.freqz(taps, [1], worN=2000, fs=fs)
>>> plot_response(w, h, "Band-stop Filter")
>>> plt.show() 

scipy.signal.unique_roots

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.unique_roots.html#scipy.signal.unique_roots

scipy.signal.unique_roots(p, tol=0.001, rtype='min')

从根列表中确定唯一根及其重数。

参数：

parray_like

根的列表。

tolfloat，可选

两个根被认为相等的公差。默认值为 1e-3。有关根分组细节，请参阅备注。

rtype，可选

如果多个根在tol范围内，则如何确定返回的根。

• ‘max’、‘maximum’：选择这些根中的最大值
• ‘min’、‘minimum’：选择这些根中的最小值
• ‘avg’、‘mean’：取这些根的平均值

在找到复根的最小或最大值时，首先比较实部，然后再比较虚部。

返回：

uniquendarray

唯一根的列表。

multiplicityndarray

每个根的重数。

备注

如果我们有根a、b和c，使得a接近b，而b接近c（距离小于tol），则并不一定意味着a接近c。这意味着根分组不是唯一的。在此函数中，我们使用“贪婪”分组，按照输入p中给定的顺序遍历根。

此实用函数不专门用于根，而是可用于需要确定唯一性和重数的任何值序列。有关更通用的程序，请参阅numpy.unique。

示例

>>> from scipy import signal
>>> vals = [0, 1.3, 1.31, 2.8, 1.25, 2.2, 10.3]
>>> uniq, mult = signal.unique_roots(vals, tol=2e-2, rtype='avg') 

检查具有大于 1 的重数的根：

>>> uniq[mult > 1]
array([ 1.305]) 

scipy.signal.residue

原文链接：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.residue.html#scipy.signal.residue

scipy.signal.residue(b, a, tol=0.001, rtype='avg')

计算 b(s) / a(s)的部分分式展开。

如果M是分子b的次数，N是分母a的次数：

 b(s)     b[0] s**(M) + b[1] s**(M-1) + ... + b[M]
H(s) = ------ = ------------------------------------------
        a(s)     a[0] s**(N) + a[1] s**(N-1) + ... + a[N] 

然后，部分分式展开 H(s)定义如下：

 r[0]       r[1]             r[-1]
= -------- + -------- + ... + --------- + k(s)
  (s-p[0])   (s-p[1])         (s-p[-1]) 

如果有重复的根（比tol更接近），则 H(s)的项如下：

 r[i]      r[i+1]              r[i+n-1]
-------- + ----------- + ... + -----------
(s-p[i])  (s-p[i])**2          (s-p[i])**n 

该函数用于正幂次 s 或 z 的多项式，如控制工程中的模拟滤波器或数字滤波器。对于 z 的负幂次（典型的数字信号处理中的数字滤波器），请使用residuez。

有关算法的详细信息，请参阅备注。

参数：

barray_like

分子多项式系数。

aarray_like

分母多项式系数。

tolfloat, optional

两个根被视为相等的距离容忍度。默认为 1e-3。详见unique_roots获取更多详情。

rtype, optional

用于计算代表一组相同根的根的方法。默认为'avg'。详见unique_roots获取更多详情。

返回：

rndarray

对应于极点的残余。对于重复的极点，残余按照幂次分数的升序排列。

pndarray

按幅度升序排列的极点。

kndarray

直接多项式项的系数。

另请参阅

invres, residuez, numpy.poly, unique_roots

注意事项

计算使用“透过减法进行紧缩”的算法 —— 第 1 条[1]。

部分分式展开的形式取决于数学上极点的重数。然而，在数值计算中无法精确确定多项式根的重数。因此，你应该将带有给定tol的residue的结果视为对具有经验确定的重数的计算极点的分部分式展开的结果。如果存在接近的极点，tol的选择可能会显著改变结果。

参考文献

[1]

J. F. Mahoney, B. D. Sivazlian，“部分分式展开：计算方法和效率综述”，《计算与应用数学杂志》，第 9 卷，1983 年。

scipy.signal.residuez

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.residuez.html#scipy.signal.residuez

scipy.signal.residuez(b, a, tol=0.001, rtype='avg')

计算 b(z) / a(z)的部分分数展开。

如果M是分子b的度数，N是分母a的度数：

 b(z)     b[0] + b[1] z**(-1) + ... + b[M] z**(-M)
H(z) = ------ = ------------------------------------------
        a(z)     a[0] + a[1] z**(-1) + ... + a[N] z**(-N) 

那么部分分数展开 H(z)的定义如下：

 r[0]                   r[-1]
= --------------- + ... + ---------------- + k[0] + k[1]z**(-1) ...
  (1-p[0]z**(-1))         (1-p[-1]z**(-1)) 

如果有任何重复的根（比tol更接近），则部分分数展开将包含以下术语：

 r[i]              r[i+1]                    r[i+n-1]
-------------- + ------------------ + ... + ------------------
(1-p[i]z**(-1))  (1-p[i]z**(-1))**2         (1-p[i]z**(-1))**n 

此函数用于负幂次 z 的多项式，例如 DSP 中的数字滤波器。对于正幂次，请使用residue。

有关算法的详细信息，请参阅残差的注释。

参数：

barray_like

分子多项式的系数。

aarray_like

分母多项式系数。

tolfloat，可选

两个根被视为相等的容差。默认为 1e-3。有关更多详情，请参见unique_roots。

rtype，可选

用于计算表示一组相同根的根的方法。默认为'avg'。有关更多详情，请参见unique_roots。

返回：

rndarray

对应于极点的残差。对于重复的极点，残差按升序排列以对应于幂分数。

pndarray

按升序排列的极点。

kndarray

直接多项式项的系数。

另请参阅

invresz, residue, unique_roots

scipy.signal.invres

原文：docs.scipy.org/doc/scipy-1.12.0/reference/generated/scipy.signal.invres.html#scipy.signal.invres

scipy.signal.invres(r, p, k, tol=0.001, rtype='avg')

从分式展开计算 b(s) 和 a(s)。

如果 M 是分子 b 的次数，N 是分母 a 的次数：

 b(s)     b[0] s**(M) + b[1] s**(M-1) + ... + b[M]
H(s) = ------ = ------------------------------------------
        a(s)     a[0] s**(N) + a[1] s**(N-1) + ... + a[N] 

然后，分式展开 H(s) 定义如下：

 r[0]       r[1]             r[-1]
= -------- + -------- + ... + --------- + k(s)
  (s-p[0])   (s-p[1])         (s-p[-1]) 

如果有任何重复的根（比 tol 更接近），那么 H(s) 就会有如下项：

 r[i]      r[i+1]              r[i+n-1]
-------- + ----------- + ... + -----------
(s-p[i])  (s-p[i])**2          (s-p[i])**n 

此函数用于正 s 或 z 的正幂多项式，如控制工程中的模拟滤波器或数字滤波器。对于 z 的负幂（DSP 中的数字滤波器典型情况），请使用 invresz。

参数：

rarray_like

对应于极点的残差。对于重复的极点，残差必须按升幂分数顺序排序。

parray_like

极点。相等的极点必须相邻。

karray_like

直接多项式项的系数。

tolfloat，可选

两个根被认为在它们之间的距离方面相等的容差。默认为 1e-3。有关详细信息，请参见 unique_roots。

rtype，可选

用于表示一组相同根的根的计算方法。默认为 ‘avg’。有关详细信息，请参见 unique_roots。

返回：

bndarray

分子多项式系数。

andarray

分母多项式系数。

另请参见

residue, invresz, unique_roots