如何在 Python 中为采样信号编写低通滤波器?

信息处理 过滤器设计 低通滤波器 Python
2021-12-31 03:07:23

我有一些信号每 1 ns(1e-9 秒)采样一次,并且有 1e4 个点。我需要从这个信号中滤除高频。假设我需要过滤高于 10 MHz 的频率。我希望对于低于截止频率的频率信号将保持不变。这意味着对于低于截止频率的频率,滤波器的增益将为 1。我希望能够指定过滤顺序。我的意思是,一阶滤波器在截止频率后有 20 分贝/十倍频斜率(功率滚降),二阶滤波器在截止频率后有 40 分贝/十倍频斜率,依此类推。代码的高性能很重要。

4个回答

使用黄油 函数设计的滤波器的频率响应为:

巴特沃斯滤波器响应

但是没有理由将滤波器限制为恒定单调滤波器设计。如果您希望在阻带中有更高的衰减和更陡峭的过渡带,则存在其他选项。有关使用iirdesing指定过滤器的更多信息,请参阅如黄油设计的频率响应图所示,截止频率(-3dB 点)远未​​达到目标。这可以通过在滤波前进行下采样来缓解(设计功能将很难使用如此窄的滤波器,2% 的带宽)。让我们看看使用指定的截止值过滤原始采样率。

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

from scipy.signal import fir_filter_design as ffd
from scipy.signal import filter_design as ifd

# setup some of the required parameters
Fs = 1e9           # sample-rate defined in the question, down-sampled

# remez (fir) design arguements
Fpass = 10e6       # passband edge
Fstop = 11.1e6     # stopband edge, transition band 100kHz
Wp = Fpass/(Fs)    # pass normalized frequency
Ws = Fstop/(Fs)    # stop normalized frequency

# iirdesign agruements
Wip = (Fpass)/(Fs/2)
Wis = (Fstop+1e6)/(Fs/2)
Rp = 1             # passband ripple
As = 42            # stopband attenuation

# Create a FIR filter, the remez function takes a list of 
# "bands" and the amplitude for each band.
taps = 4096
br = ffd.remez(taps, [0, Wp, Ws, .5], [1,0], maxiter=10000) 

# The iirdesign takes passband, stopband, passband ripple, 
# and stop attenuation.
bc, ac = ifd.iirdesign(Wip, Wis, Rp, As, ftype='ellip')  
bb, ab = ifd.iirdesign(Wip, Wis, Rp, As, ftype='cheby2') 

原始采样率过滤器

如前所述,因为我们试图过滤如此小百分比的带宽,所以滤波器不会有一个尖锐的截止。在这种情况下,使用低通滤波器,我们可以降低带宽以获得更好看的滤波器。python/scipy.signal resample 函数可用于减少带宽。

请注意,重采样函数将执行过滤以防止混叠。也可以执行预过滤(以减少混叠),在这种情况下,我们可以简单地重新采样 100 并完成,但问题是关于创建过滤器。对于这个例子,我们将下采样 25 并创建一个新的过滤器

R = 25;            # how much to down sample by
Fsr = Fs/25.       # down-sampled sample rate
xs = signal.resample(x, len(x)/25.)

如果我们更新 FIR 滤波器的设计参数,新的响应是。

# Down sampled version, create new filter and plot spectrum
R = 25.             # how much to down sample by
Fsr = Fs/R          # down-sampled sample rate
Fstop = 11.1e6      # modified stopband
Wp = Fpass/(Fsr)    # pass normalized frequency
Ws = Fstop/(Fsr)    # stop normalized frequency
taps = 256
br = ffd.remez(taps, [0, Wp, Ws, .5], [1,0], maxiter=10000) 

下采样滤波器响应

对下采样数据进行操作的滤波器具有更好的响应。使用 FIR 滤波器的另一个好处是您将获得线性相位响应。

这行得通吗?

from __future__ import division
from scipy.signal import butter, lfilter

fs = 1E9 # 1 ns -> 1 GHz
cutoff = 10E6 # 10 MHz
B, A = butter(1, cutoff / (fs / 2), btype='low') # 1st order Butterworth low-pass
filtered_signal = lfilter(B, A, signal, axis=0)

不过,您是对的,文档不是很完整。它看起来像是butter一个包装器iirfilter有更好的文档记录

N : int 过滤器的顺序。Wn : array_like 给出临界频率的标量或长度为 2 的序列。

不过,这些东西大部分是从 matlab 克隆的,所以你也可以查看他们的文档

归一化截止频率 Wn 必须是介于 0 和 1 之间的数字,其中 1 对应于奈奎斯特频率,每个样本 π 弧度。

更新:

我为这些功能添加了文档:) Github 让它变得简单。

不确定您的应用程序是什么,但您可能想查看 Gnuradio:http ://gnuradio.org/doc/doxygen/classgr__firdes.html

信号处理块是用 C++ 编写的(尽管 Gnuradio 流程图是用 Python 编写的),但您确实说过高性能很重要。

这个 FIR 滤波器效果很好。注意它两次应用滤波器,“前进”和“后退”,以补偿信号偏移(filtfilt功能不起作用,不知道为什么):

def firfilt(interval, freq, sampling_rate):
    nfreq = freq/(0.5*sampling_rate)
    taps =  sampling_rate + 1
    a = 1
    b = scipy.signal.firwin(taps, cutoff=nfreq)
    firstpass = scipy.signal.lfilter(b, a, interval)
    secondpass = scipy.signal.lfilter(b, a, firstpass[::-1])[::-1]
    return secondpass

THIS是过滤器设计和使用的重要资源,我从哪里获取此代码,以及从哪里可以获取带通和高通滤波器示例