滤波器的因果稳定性

这里简单总结一下滤波器系统的因果性以及稳定性

1.因果性

1.1 因果系统：

系统某时刻的输出只取决于该时刻和该时刻以前的输入
如一个常用的滑动平均滤波器，
$$h(n)=[0.2,0.2,0.2,0.2,0.2]，0⩽n⩽4$$
z变换为
$$\sum _{n=0}^{n=4}h(n)z^{-n}$$
差分方程为
$$y(n)=(x(n)+x(n-1)+x(n-2)+x(n-3)+x(n-4))/5$$
这是一个常见的FIR滤波器，系统当前的输出只与当前的输入和之前的四个输入有关，因此这是一个因果系统，看下这个因果滤波器的效果

% close all
h = [1 1 1 1 1 ]/5;
L = -2:2;
figure,plot(L,h)

t = linspace(-pi,pi,100);
rng default  %initialize random number generator
x = sin(t) + 0.25*rand(size(t));

y = conv(x,h);
% y = filter(b,a,x);
figure,plot(x)
hold on,plot(y(1:length(x))),legend('x','y')

可以看到，输出相对于输入有延迟，这个延迟是可以计算得到的，线性相位FIR滤波器的输出有$(N-1)/2$个抽样的延迟

1.2 非因果系统

非因果系统:系统当期的输出与未来的数据有关
还是上面的滑动滤波器，时间区间变为[-2，2]
$$h(n)=[0.2,0.2,0.2,0.2,0.2]，-2⩽n⩽2$$
则z变换相应变为为
$$\sum _{n=-2}^{n=2}h(n)z^{-n}$$
差分方程变为
$$y(n)=(x(n+2)+x(n+1)+x(n)+x(n-1)+x(n-2))/5$$
看看非因果滤波器的效果

% close all
h = [1 1 1 1 1 ]/5;
L = -2:2; % non-casual FIR
N = length(h);
figure,plot(L,h)

t = linspace(-pi,pi,100);
rng default  %initialize random number generator
x = sin(t) + 0.25*rand(size(t));

y = conv(x,h);
% y = filter(b,a,x);
y = zeros(1,length(x));
% x = [zeros(1,(N-1)/2),x,zeros(1,(N-1)/2)];
for i = (N-1)/2:length(x)-(N-1)/2-1
    % positive index
    n = i+1;
    y(n) = (x(n+2)+x(n+1)+x(n)+x(n-1)+x(n-2))/5;
end
figure,plot(x)
hold on,plot(y(1:length(x))),legend('x','y')

代码中循环部分可以清楚看出计算过程，结果如下
在这里插入图片描述

对照上面因果滤波器的结果发现，非因果的输出相对输出没有延时，这就是一个零相位滤波器。
　通常所说的非因果滤波器物理不可实现，是因为在实时系统中无法预知未来的输入。但是，非因果滤波器也不是没有用处了，如果对延时没有严格要求，可以先缓存一段数据，再用非因果滤波器处理，最后再补偿这段延迟，上面的程序就是这样的过程。
　因此，非因果滤波器物理实现过程就是先当非因果滤波器计算，最后再做延时补偿
$y1(n+4)=(x(n+4)+x(n+3)+x(n+2)+x(n+1)+x(n))/5$
$y(n)=y1(n+(N-1)/2)$
　　说到这里也我们再看一下频域滤波的过程，因为频域处理通常是块处理，需要缓存一段数据，流程可以简单表示如下也

从处理过程上看，每次块处理完后，输出得到的$y(t-N+1)$ 这个点的数据，在计算过程用到了$x(t)$，因此频域滤波器是一个非因果系统，也就是文献中通常所说的频域非因果滤波器，当然，前面也介绍了，实现非因果系统的代价就是输出相对输入会有延时

2.稳定性

稳定性是当系统输入有界时，系统输出也有界。线性时不变（LTI）系统稳定的充要条件是
$$\sum _{n=-\infty }^{n=+\infty }\left|h(n)\right|<+\infty$$
即系统的单位脉冲响应绝对可和，从这个条件上看，FIR系统的脉冲响应是有限长的，因此FIR系统是绝对稳定的。

2.1. 因果系统的稳定性

既然FIR系统一定是稳定的，那就来看下IIR，一个简单的IIR滤波器Z方程如下
　　$$H(z)=\frac{b_0}{1-a_1{z}^{-1}}$$
　　差分方程为
　$$y(n)=b_0\ast x(n)-a_1\ast y(n-1)$$
　将使$H(z)$为0的点称为系统零点，使分母为0的点称为系统极点，
　上式中没有零点，只有一个极点$z=a_1$，前面分析过FIR是稳定的，因此零点不影响系统的稳定性，那我们就只分析极点，
　就先令$b_0=1$，可求得上式中的$H(z)的$逆变换单位冲激响应$$h(n)=a_1^n，0⩽n⩽+\infty$$
　当$a_1>1，如a_1=2时$，系统存在极点$z=2$，可以看到h(n)是呈指数增长的，求和一定是无限大的，而从差分方程上看，系统当前的输出是当前的输入加上两倍的上一刻输出，$y(n)$会越变越大，表示为系统存在一个自激振动的过程，因此系统是不稳定的
　当$a_1<1，如a_1=1/2时$，系统存在极点$z=1/2$
　$$\underset{N\to +\infty }{\lim }\sum _{n=0}^N\left|2^{-n}\right|=2$$
　系统冲击响应是绝对可和的，因此系统是稳定的.
　这也对应着稳定性的一个结论：
　对于一个稳定的因果系统，系统函数的全部极点必须在z平面的单位圆内
　给定分子分母系数$b和a$，matlab可以方便的求出零极点，以上面的$H(z)$为例

b = 1;
a=[1 -0.5];  %注意matlab的z方程分母是1+a1*z^-1+....an*z^-n
Hq = dfilt.df2(b,a);
[zq,pq,kq]=zplane(Hq);
fvtool(b,a)

在打开的fvtool里，可以看到这个关于滤波器的幅度响应、相位响应、零极点分布、稳定性等在内的所有信息。