FIR维纳滤波器的设计
在信号处理的许多实际应用中，人们往往不能直接获得所需要的有用信号，需要从噪声中提取信号。
比如，在信号传输过程中，由于存在信道噪声等干扰，在接收端观测到的信号必然与原始信号不同。
为了从观测数据中尽可能精确地重现原始信号，而最大成都地抑制噪声，需要设计一种滤波器，其输出尽可能逼近原始信号，成为原始信号的最佳估计。
这种滤波器成为最佳滤波器。
维纳(Wiener)滤波器就是用来解决这样一类问题的一种滤波器。
本文将应用MATLAB并结合实例介绍FIR维纳滤波器的设计方法。
一、维纳滤波的原理
维纳滤波的本质是一种最佳估计问题，采用的是最小均方误差准则。
一个线性系统，其单位样本响应为h(n)，当输入一个随机信号)
(
)
(
)
(n
n
s
n
xυ
+
=其中s(n)表示信号，)
(n
υ表示噪声，则输出y(n)为
∑-=
m m
n
x
m
h
n
y)
(
)
(
)
( (1) 系统是通过y(n)来估计s(n)，因此将其称为s(n)的估计值，用)
(ˆn
s表示，即
)
(ˆ
)
(n
s
n
y=(2)
h(n)
)
(ˆ
)
(n
s
n
y=
)
(
)
(
)
(n
n
s
n
xυ
+
=
图1 维纳滤波器基本框图
图1所示为维纳滤波器的基本框图。
式(1)为一卷积，可以理解为从当前和过去的观察值x(n)，x(n-1)，x(n-2)…x(n-m)，…来估计信号的当前值)
(ˆn
s。
维纳滤波器一般有三种用途。
用当前的和过去的观察值x(n)，x(n-1)，x(n-2)，…来估计当前的信号值)
(ˆ
)
(n
s
n
y=称为滤波；用过去的观察值来估计当前的或将来的信号值)0
)(
(ˆ
)
(≥
+
=N
N
n
s
n
y称为预测；
用全部数据来估计过去的信号值)1)((ˆ)(>-=N N n s
n y 称为平滑。
维拉滤波采用的是最小均方误差准则，以)(n s 与)(ˆn s
分别表示信号的真值与估计值，而用e(n)表示它们之间的误差，即
)(ˆ)()(n s
n s n e -= (3)
那么最小均方误差表示为
)(n ξ=[][]min ))(ˆ)(()(22=-=n s
n s E n e E (4)
为了得到使式(4)最小时的维纳滤波器冲击响应，令)(n ξ对)(m h 的导数等于零，即
0)]()([2])
()
()([2)()(=--=∂∂=∂∂m n x n e E m h n e n e E m h n ξ (5) 由式(5)可得
0)]()([=-m n x n e E (6) 式(6)称为正交方程，表明任何时刻的估计误差都与用于估计的所有数据正交。
将式(1)，式(3)代入式(6)可得
)()()(m R m h R xx m
sx -=∑ττ (7)
式(7)称为维纳-霍夫（Wiener-Hopf ）方程，其中)(τsx R 是)(n s 与)(n x 的互相过函数，)(τxx R 是)(n x 的自相关函数，其表达式如下
)]()([)(ττ+=n x n s E R sx (8) )]()([)(ττ+=n x n x E R xx (9)
维纳滤波器有三种情况，一是FIR 维纳滤波器；二是非因果IIR 维纳滤波器；三是因果的IIR 维纳滤波器。
本文仅讨论FIR 维纳滤波器的求解。
二、FIR 维纳滤波器求解
设维纳滤波器的冲击响应序列长度为N ，则冲击响应矢量为 h=T N h h h )]1()1()
0([-
(10)
滤波器输入数据矢量为
)(n x =T N n x n x n x )]1()1()
([+--
(11)
着滤波器输出为
=)(n y )()()(ˆn x h h n x n s
T T == (12) 这样，维纳-霍夫方程可表示成
P=Rh (13)
其中P=)]()([n x n s E ，R=])()([T n x n x E 。
求解式(13)可得
h=P R 1
-
三、实例设计
结合维拉滤波器的原理，利用MATLAB 设计滤波器。
将随机信号X(n)看成是由典型白噪声序列源W(n)激励一个线性系统产生，
用
一
个
差
分
方程
(
)1
.3
5
x n x
n
x n -
-+
---
来描述。
进
行Z 变换得到
1233()1()()1 1.352 1.3380.6620.240X z H z W z z z z z ----=
=
-+-+，那么均
值为１的高斯白噪声序列W （n ）可以用randn 函数产生，再利用函数X=filter(B,A,W)产生随机信号x(n)。
这里将滤波器的阶数设为101，根据维纳－霍夫方程：
-1[]rx1rx2
opt h =⨯，其中rx1是观测信号的自相关函数，rx2是观测信号和期望信
号的互相关函数。
定义维纳滤波的模型，最后带入fi lter 。
具体MATLAB 程序代码如下：
maxlag=100; N=100; %采样次数为100 x=zeros(N,1); y=zeros(N,1); var=1;
%%%%%%%%%%%%%%%%%%列出状态方程 x(1)=randn(1,1); %令x(-1)=x(-2)=x(-3)=x(-4)=0 x(2)=randn(1,1)+1.352*x(1);
x(3)=randn(1,1)+1.352*x(2)-1.338*x(1);
x(4)=randn(1,1)+1.352*x(3)-1.338*x(2)+0.602*x(1); for n=5:N
x(n)=1.352*x(n-1)-1.338*x(n-2)+0.602*x(n-3)-0.24*x(n-4)+randn(1,1); end;
v=randn(N,1); y=x+v; %z_x 为观测样本值=真值+噪声 %%%%%%%%%%%%%%%%%%%%%%%滤波 x = x'; y = y';
xk_s(1)=y(1); %赋初值 xk_s(2)=y(2); xk_s(3)=y(3); xk_s(4)=y(4);
xk=[y(1);y(2);y(3);y(4)];
%%%%%%%%%%%%%%%%%%计算观测函数 %维纳滤波器的生成
[rx,lags]=xcorr(y,maxlag,'biased');%观测信号的自相关函数
rx1=toeplitz(rx(101:end));%对称化自相关函数矩阵使之成为方阵，滤波器的阶数为101阶 rx2=xcorr(x,y,maxlag,'biased');%观测信号与期望信号的互相关函数 rx2=rx2(101:end); h=inv(rx1)*rx2'; %维纳-霍夫方程 xk_s=filter(h,1,y);%加噪信号通过滤波器后的输出 %%%%%%%%%%%%%%%%%%%%%计算误差 e_x=0; eq_x=0; e_x1=N:1;
%计算滤波的均值,计算滤波误差的均值
for i=1:N
e_x(i)=x(i)-xk_s(i); %误差=真实值-滤波估计值
end %%%%%%%%%%%%%%%%%%%%%%%作图
t=1:N;
figure(1);
plot(t,x,'r',t,y,'g',t,xk_s,'b');
legend('实际轨迹','观测样本','估计轨迹');
figure(2);
plot(e_x, 'r');
legend('平均误差曲线');
运行上述MATLAB程序，得到结果如下。
图2 滤波结果对比
从图2、图3中可以看出，设计能很好地估计原始数据。