将过滤实现为乘法或卷积在数值上更稳定吗？

我正在编写一个程序，以使用五阶Butterworth滤波器离线过滤20,000个采样信号。我可以访问FFT实现。似乎有两种选择可以实现滤波：

在时域中将信号与脉冲响应进行卷积，或者将信号与频率中的脉冲响应相乘域和反变换结果

在理论上FT情况下，这些方法是相同的。但是，我认为使用DFT在现实生活中进行操作是不同的。其中一种方法在数值上更稳定吗？还有其他我应该注意的问题吗？计算数量并不重要。

对于大多数信号长度，FFT方法将更快地进行计算。使用时域卷积，只有较短的长度更快。

#1 楼

使用卷积，您不会遇到任何稳定性问题，因为没有递归过滤，因此您不会累积任何错误。换句话说，系统全为零，没有极点。我曾经听说过，但没有为我自己证实，基于FFT的卷积比时域卷积具有更低的误差，仅仅是因为它具有O（n log n）个算术运算而不是O（n ^ 2）。
就我所知，通常情况下，巴特沃思过滤器是作为递归（IIR）过滤器实现的，所以这是一个不同的话题。 IIR滤波器具有极点和零点，因此实际上可能存在稳定性问题。另外，对于IIR滤波器，不是基于FFT的方法不是一种选择，但是从好的方面来说，IIR滤波器的阶数往往很低。

关于IIR滤波器的稳定性问题，它们倾向于高阶问题-我会抛出一个数字，说大约6阶正推动它。取而代之的是，它们通常被实现为级联双二阶（二阶滤波器部分）。对于您的5阶滤波器，将其编写为z域传递函数（它将是5阶有理函数），然后将其分解为5个极点和5个零。收集复共轭，您将拥有两个双二阶和一个一阶滤波器。通常，随着极点越来越靠近单位圆，稳定性问题趋于出现。

IIR滤波器中还可能存在噪声和极限周期问题，因此存在不同的滤波器拓扑结构（即直接滤波器）。形式I，直接形式II）具有不同的数值属性，但我不会想过这一点-仅使用双精度，几乎可以肯定就足够了。

$ \ begingroup $
它仅适用于FIR滤波器是什么意思？我认为IIR滤波器无论如何都必须进行采样。通常会在时域中实现IIR滤波器来避免这种情况吗？
$ \ endgroup $
–安德烈亚斯（Andreas）
2011年11月7日，下午5:27

$ \ begingroup $
据我所知，IIR滤波器始终在时域中实现。通过诸如y（n）= b0 * x（n）+ b1 * x（n-1）+ b2 * x的差分方程来定义IIR滤波器（例如，此处为二阶滤波器或“二阶”）。（n-2）-a1 * y（n-1）-a2 * y（n-2）。请注意，这是先前输入样本（x值）和先前输出样本（y值）的组合。 FIR滤波器仅取决于过去的输入，因此可以接受高效的频域实现。 IIR滤波器没有，但无论如何还是非常有效的，因为IIR滤波器往往要低很多。
$ \ endgroup $
– schnarf
2011-12-7 22:03

$ \ begingroup $
IIR滤波器趋向于低阶的原因是，与FIR滤波器相比，极点（先前输出样本的反馈）使滤波器变得陡峭，系数很小。当我说低得多时，典型的IIR滤波器可能是二阶的（5个系数），而用于同一任务的典型FIR滤波器可能具有数千个系数。
$ \ endgroup $
– schnarf
2011年12月7日在22:09

#2 楼

在几乎所有情况下，您的最佳选择既不是卷积也不是FFT，而只是直接应用IIR滤波器（使用例如sosfilt（）函数）。就CPU和内存消耗而言，这将大大提高效率。

是否产生数值差异取决于特定的过滤器。如果极点非常非常靠近单位圆，则可能会出现一些差异。即使有一些技巧也可以提供帮助。请勿使用传递函数表示形式和filter（），而应将极点和零与sosfilt（）一起使用。这是一个差异的示例。

filter（）在约44.1kHz的15Hz截止频率时变坏。对于sosfilt（），截止频率可以很好地低于1/100 Hz @ 44.1kHz，没有任何问题。

如果您遇到稳定性问题，FFT也无济于事。由于您的滤波器是IIR滤波器，因此脉冲响应是无限的，因此必须首先将其截断。在这些非常低的频率下，脉冲响应变得很长，以至于FFT也变得不切实际。

例如，如果您想要在44.1 kHz处的1/100 Hz截止频率，并且想要100 dB脉冲响应中的动态范围，则大约需要2500万个采样！在44.1 kHz的频率下将近10分钟，比原始信号长很多倍

$ \ begingroup $
这并不能真正回答有关数值问题的问题，但是我并不知道过滤器存在的问题-谢谢！我的高通截止频率为0.5 Hz @ 250 Hz。过滤器出现问题的原因是什么？我自己编写实现。
$ \ endgroup $
–安德烈亚斯（Andreas）
2011年11月8日，晚上7:52

#3 楼

您为什么认为情况会有所不同？理论概念应转化为实际应用，唯一的区别是浮点问题，这是我们无法逃脱的。您可以使用MATLAB中的一个简单示例轻松验证：

x=randn(5,1);
y=randn(5,1);
X=fft(x,length(x)+length(y)-1);
Y=fft(y,length(x)+length(y)-1);

z1=conv(x,y);z2=ifft(X.*Y);
z1-z2

ans =

   1.0e-15 *

   -0.4441
   -0.6661
         0
   -0.2220
    0.8882
   -0.2220
         0
   -0.4441
    0.8882

如您所见，误差大约是机器精度。应该没有任何理由不使用FFT方法。如Endolith所述，使用FFT方法进行滤波/卷积/互相关等操作更为普遍，除了非常小的样本（如本例中所示）外，它要快得多。并不是说永远都不会进行时域处理...所有这些都归结为应用程序，需求和约束。

$ \ begingroup $
我认为最初的问题是钻探基于FFT的滤波固有的浮点问题，而不是在时域直接实现滤波器。例如，如果滤波器的确很长，或者FFT的实现不正确，这可能是定点信号处理的真正问题。对于双精度浮点数中长度为5的序列，您绝对不会看到任何影响。
$ \ endgroup $
–Jason R
2011年11月7日15:14

$ \ begingroup $
@JasonR如果在上述示例中将序列的长度扩展为1e6，则错误仍然具有机器精度。您提到的错误主要是由于不良的滤波器设计或不良的FFT实施造成的。如果这些都很好，我不明白为什么时域的卷积应该给出与频域不同的答案。
$ \ endgroup $
–乳香
2011年11月7日15:41

$ \ begingroup $
不应根据您在哪个域中进行计算而给出不同的答案。我唯一的要点是两种方法之间的实际数学运算差异很大，因此，取决于每种可用的实现，您可能会看到明显的差异。对于双精度，假设您具有良好的实现，除非在极端情况下，否则我不会期望有任何区别。
$ \ endgroup $
–Jason R
2011年11月7日15:45

编程黑洞网

将过滤实现为乘法或卷积在数值上更稳定吗？

评论

#1 楼

评论

#2 楼

评论

#3 楼

评论