有没有人熟悉 Gustafson 的算法，用于最小化前向反向滤波中的瞬态 [1]？我正在尝试实现它，我的第一个猜测是检查 Matlab 的 filtfilt.m，因为他们引用了这篇论文。在 Matlab 函数中，还求解了一个线性方程组，以找到最小化启动瞬态的初始条件 zi，但参考和代码之间的关系对我来说并不明显。关于最小化的唯一代码行是（nfilt 是系数向量的长度）：
zi = ( eye(nfilt-1) - [-a(2:nfilt), [eye(nfilt-2); zeros(1,nfilt-2)]] ) \...
( b(2:nfilt) - b(1)*a(2:nfilt) );
谁能指出我正确的方向，这些线与古斯塔夫森文章中描述的算法有何关系？
[1] Gustafsson, F.“确定前向后向过滤的初始状态”。IEEE® 信号处理汇刊。卷。44，1996 年 4 月，第 988-992 页。
对于任何感兴趣的人，我偶然发现了一篇描述在 matlab 的 filtfilt.m 中实现的方法的论文。附上该论文的链接。至少据我了解，matlab 的 filtfilt.m 没有实现 Gustafson 算法。
OP 问题中的zi = (...)\(...)行决定了过滤器的初始状态。我相信 Python 的scipy. 根据scipy 文档：
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 的输入的响应为常数的初始条件。
（我的重点）
请注意，filtfilt计算初始状态zi，按第一个样本值对其进行缩放，然后将其传递给filter实际应用它的 （文档）。
此处zi提供了如何使用状态空间表示在过滤器中应用此初始状态的基本示例。