python(scipy)中的可逆低通(Butterworth)滤波器?

信息处理 过滤器 Python 低通滤波器 scipy
2022-01-29 10:51:43

问题

scipy.signal(或其他python包)是否内置了可逆低通滤波器?如果是这样,它是什么?如果没有,为什么不呢(反转低通滤波器是否特别困难)?

细化

我需要一个可逆(数字,一阶,具体而言)低通滤波器,例如巴特沃斯滤波器。我需要应用它的倒数来用预补偿模拟信号。信号经过滤波后由数字转换为模拟。(x,y)=(time,voltage)

迄今为止完成的工作

天真地,(据我所知)我可以做到以下几点:

import numpy as np 
import scipy
# define signal here
num_order = 1
cutoff_frequency = 1e4 # Hz 
time_constant = 1/(2*np.pi*cutoff_frequency) 
lpf_b, lpf_a = scipy.signal.butter(num_order, btype='low', Wn=1/time_constant/(2*np.pi*sample_rate/2))
precompensated_signal = scipy.signal.lfilter(lpf_a, lpf_b, signal)

其中一个示例信号被定义为(例如)

sample_rate = 1e6
t_period = 2*1e-3
num_samples = int(round(t_period*sample_rate))
amplitude = 0.23
signal = amplitude*np.ones(num_samples)

这种方法在反转高通滤波器时有效。

在正向应用过滤器 ( filtered_signal)

filtered_signal = scipy.signal.lfilter(lpf_b, lpf_a, signal)

产生预期的信号(通过运行代码检查或查看图片): scipy.signal.butter 正向

但是将其向后 ( precompensated_signal) 应用到预补偿会产生一个振荡信号:

scipy.signal.butter 反向

查看系数,我发现

lpf_a = [ 1.0, -0.93906251]
lpf_b = [0.03046875, 0.03046875]

似乎在 b 中具有两个相同的系数作为第一个也是唯一的系数是使滤波器不可逆的原因。

2016年有一个相关的问题和答案。根据答案,取给定的滤波器系数scipy.signal.butter并修改如下

lpf_b, lpf_a = scipy.signal.butter(num_order, btype='low', Wn=1/time_constant/(2*np.pi*sample_rate/2))
lpf_b_2 = [1 + lpf_a[1]] # note that lpf_b_2[0] = lpf_b[0] + lpf_b[1] 
filtered_signal_2 = scipy.signal.lfilter(lpf_b_2, lpf_a, signal)
precompensated_signal_2 = scipy.signal.lfilter(lpf_a, lpf_b_2, signal)

产生一个过滤器,它在正向 ( filtered_signal_2) 上的行为与原始过滤器相同:

带修正滤波器的预补偿信号(只有一个 b 系数)

并且在向后方向上表现得更好(precompensated_signal_2):

在此处输入图像描述

尽管产生的信号仍然值得怀疑。唯一与幅度不同的系数是第一个(即使改变截止频率,这似乎也是如此),令人担忧的是,第一个数据点的电压变得非常大。对于上述链接答案中的过滤器,此功能似乎是正确的。

编辑

ZR Han建议(谢谢!)移动过滤器的极点。据我了解,可以这样做:

z, p, k = scipy.signal.tf2zpk(lpf_b, lpf_a)
# the zero z = [-1.0] as he suggested
z = [-0.95]
lpf_b, lpf_a = scipy.signal.zpk2tf(z, p, k)

这导致

lpf_b = [0.03046875, 0.02894531]
lpf_a = [ 1.0, -0.93906251]

但不幸的是,这并没有消除所有的振荡,通过绘制预补偿信号 ( ) 可以看出scipy.signal.lfilter(lpf_a, lpf_b, signal) with the modified coeffs

带零的预补偿信号向单位圆移动

也许这是由于 z/my 对什么是“小转变”缺乏正式的选择而造成的。由于预补偿信号被转换为模拟信号,不幸的是,该解决方案不能容忍不良行为。

从这里出发的路线

  1. 解释修改后的过滤器对应的内容
  2. 从 scipy 的答案中找到过滤器
  3. 寻找另一个更适合用于转换为模拟信号的滤波器(平滑和有限)

背景:我在数字信号处理方面没有正式的背景,但我可以学习我需要的东西(如果我指向来源,会更快)。

2个回答

是否有可逆低通滤波器

反转低通滤波器有什么特别困难的地方吗?

是的。数字低通滤波器(在最常见的意义上)在奈奎斯特处为零,这意味着逆在奈奎斯特处具有无限增益并且是不稳定的。

似乎在 b 中具有两个相同的系数作为第一个也是唯一的系数是使滤波器不可逆的原因。

是的。这在奈奎斯特创造了零。

一种过滤器,其行为与正向的原始过滤器相同

它不是。将“b”系数更改为单个系数将使极点远离实轴上的单位圆,从而改变高频行为。信号看起来只是一样,因为它的高频成分非常少

根据您的应用程序,您可能最好使用免费的架子过滤器。低通滤波器会破坏不可恢复的信息。实际上,如果任何滤波器具有大量衰减,则反滤波器将具有大量增益,这在噪声性能方面是有问题的。

背景:我在数字信号处理方面没有正式的背景,但我可以学习我需要的东西(如果我指向来源,会更快)。

数字信号处理的数学相当繁重,对理论基础有所了解确实很有帮助。一个好的开始可能是:

https://ccrma.stanford.edu/~jos/filters/

按照Hilmar的建议使用所谓的搁置过滤器(谢谢!),与endolith链接的biquad_cookbook模块(谢谢!)。该软件包(编写于 2013 年)与 (conda) python 3.8.10 一起运行,开箱即用。

小例子:

import numpy as np
import scipy
import biquad_cookbook # need the file in the working directory

precompensated_signals = []
precompensated_and_filtered_signals = []
for g in [3, 6, 9, 12]:
    shelf_b, shelf_a = biquad_cookbook.shelf(Wn=1/time_constant/(2*np.pi*sample_rate/2), 
        dBgain=g, S=1, btype='high', ftype='half', analog=False, output='ba')
    precompensated_signals.append(scipy.signal.lfilter(shelf_b, shelf_a, signal))
    precompensated_and_filtered_signals.append(scipy.signal.lfilter(lpf_b, lpf_a, 
        scipy.signal.lfilter(shelf_b, shelf_a, signal)))

其中lpf_alpf_b在问题中被定义,以产生平滑、有限、尽管轻微振荡的预补偿信号(precompensated_signals注意我们在前向应用了架式滤波器 - 滤波器直接是预补偿滤波器):

在此处输入图像描述

将这些通过原始低通滤波器(precompensated_and_filtered_signals建模为问题中的一阶巴特沃斯滤波器)会产生以下信号:

在此处输入图像描述

所以,这是一个部分解决方案。例如,在 12 dB 增益下,与未预补偿信号相比,信号在大约十分之一的时间内达到所需幅度。权衡是幅度将由于预补偿信号中的纹波(由滤波器产生)而过冲。