使用 Scipy 的低通 Chebyshev I 型滤波器答案

【问题标题】：Low-pass Chebyshev type-I filter with Scipy使用 Scipy 的低通 Chebyshev I 型滤波器
【发布时间】：2018-03-28 19:56:40
【问题描述】：

我正在阅读一篇论文，试图重现论文的结果。在本文中，他们对原始数据使用低通 Chebyshev I 型滤波器。他们给出了这些参数。

采样频率 = 32Hz, Fcut=0.25Hz, Apass = 0.001dB, Astop = -100dB, Fstop = 2Hz, Order of the filter = 5。我找到了一些资料来帮助我理解这些参数

但是当我查看 scipy.signal.cheby1 时。这个函数需要的参数不同。

cheby1(N, rp, Wn, btype='low', analog=False, output='ba')

这里N：过滤器的顺序； btype：过滤器的类型，在我的例子中，它是“低通”； analog=False，因为数据是采样的，所以是数字的；输出：指定输出的类型。但我不确定 rp, Wn。

在文档中，它说：

rp : 浮动通带中允许低于单位增益的最大纹波。以分贝为单位指定，为正数。

Wn : array_like 给出临界频率的标量或长度为 2 的序列。对于 I 型滤波器，这是过渡带中增益首先下降到 -rp 以下的点。对于数字滤波器，Wn 从 0 归一化到 1，其中 1 是奈奎斯特频率，pi 弧度/样本。（因此，Wn 以半周期/样本为单位。）对于模拟滤波器，Wn 是角频率（例如 rad/s）。

根据这个问题： How To apply a filter to a signal in python

我知道如何使用过滤器。但我不知道如何创建具有与上述相同参数的过滤器。我不知道如何转换这些参数并将它们提供给 Scipy 中的函数。

【问题讨论】：

标签： numpy filter scipy signal-processing data-processing

【解决方案1】：

看看the wikipedia page on the Type I Chebyshev filter。请注意，您的绘图说明了一般过滤器的特性。然而，低通 I 型切比雪夫滤波器在阻带中没有波纹。

您有三个用于设计 I 型切比雪夫滤波器的可用参数：滤波器阶数、纹波因子和截止频率。这是scipy.signal.cheby1的前三个参数：

cheby1 的第一个参数是过滤器的顺序。
第二个参数 rp 对应于维基百科页面中的 δ，显然就是您所说的 Apass。

第三个参数是wn，截止频率表示为奈奎斯特频率的一部分。在你的情况下，你可以写类似

fs = 32      # Sample rate (Hz)
fcut = 0.25  # Desired filter cutoff frequency (Hz)

# Cutoff frequency relative to the Nyquist
wn = fcut / (0.5*fs)

一旦选择了这三个参数，所有其他特征（例如过渡带宽度、Astop、Fstop 等）被确定。因此，您给出的规范“采样频率 = 32Hz，Fcut=0.25Hz，Apass = 0.001dB，Astop = -100dB，Fstop = 2Hz，滤波器阶数 = 5”似乎与 I 型不兼容切比雪夫滤波器。特别是，我在 2 Hz 时获得了大约 -78 dB 的增益。（如果将阶数增加到 6，则 2 Hz 处的增益约为 -103。）

这是一个完整的脚本，后面是它生成的情节。该图仅显示通带，但您可以更改 xlim 和 ylim 函数的参数以查看更多信息。

import numpy as np
from scipy.signal import cheby1, freqz
import matplotlib.pyplot as plt


# Sampling parameters
fs = 32  # Hz

# Desired filter parameters
order = 5
Apass = 0.001  # dB
fcut = 0.25    # Hz

# Normalized frequency argument for cheby1
wn = fcut / (0.5*fs)

b, a = cheby1(order, Apass, wn)

w, h = freqz(b, a, worN=8000)

plt.figure(1)
plt.plot(0.5*fs*w/np.pi, 20*np.log10(np.abs(h)))
plt.axvline(fcut, color='r', alpha=0.2)
plt.plot([0, fcut], [-Apass, -Apass], color='r', alpha=0.2)
plt.xlim(0, 0.3)
plt.xlabel('Frequency (Hz)')
plt.ylim(-5*Apass, Apass)
plt.ylabel('Gain (dB)')
plt.grid()
plt.title("Chebyshev Type I Lowpass Filter")
plt.tight_layout()

plt.show()

【讨论】：