实验一 维纳滤波1. 实验内容设计一个维纳滤波器：（1) 产生三组观测数据,首先根据()(1)()s n as n w n =-+产生信号()s n ，将其加噪，（信噪比分别为20,10,6dB dB dB ）,得到观测数据123(),(),()x n x n x n 。

（2） 估计()i x n ，1,2,3i =的AR 模型参数。

假设信号长度为Ｌ,AR 模型阶数为N,分析实验结果,并讨论改变L,N 对实验结果的影响。

2. 实验原理滤波目的是从被噪声污染的信号中分离出有用的信号来，最大限度地抑制噪声。

对信号进行滤波的实质就是对信号进行估计。

滤波问题就是设计一个线性滤波器,使得滤波器的输出信号()y n 是期望响应()s n 的一个估计值。

下图就是观测信号的组成和信号滤波的一般模型。

观测信号()()()x n s n v n =+ 信号滤波的一般模型 维纳滤波解决从噪声中提取信号的滤波问题,并以估计的结果与真值之间的误差均方值最小作为最佳准则。

它根据()()(),1,,x n x n x n m --估计信号的当前值,它的解以系统的系统函数()H z 或单位脉冲()h n 形式给出，这种系统常称为最佳线性滤波器。

维纳滤波器设计的任务就是选择()h n ,使其输出信号()y n 与期望信号()d n 误差的均方值最小。

假设滤波系统()h n 是一个线性时不变系统，它的()h n 和输入信号都是复函数，设()()()h n a n jb n =+ 0,1,n=考虑系统的因果性,可得到滤波器的输出()()()()()0*m y n h n x n h m x n m +∞===-∑ 0,1,n=设期望信号()d n ，误差信号()e n 及其均方误差()2E e n ⎡⎤⎣⎦分别为()()()()()e n d n y n s n y n =-=-()()()()()()2220m E e n E d n y n E d n h m x n m ∞=⎡⎤⎡⎤⎡⎤=-=--⎢⎥⎣⎦⎣⎦⎢⎥⎣⎦∑ 要使均方误差为最小,需满足:()()20E e n h j ⎡⎤∂⎣⎦=∂ 整理得()()0E x n j e n *⎡⎤-=⎣⎦，等价于()()0E x n j e n *⎡⎤-=⎣⎦上式说明,均方误差达到最小值的充要条件使误差信号与任一进入估计的输入信号正交，这就是正交性原理。

将()()0E x n j e n *⎡⎤-=⎣⎦展开，得()()()()00m E x n k d n h m x m +∞***=⎡⎤⎛⎫--=⎢⎥ ⎪⎝⎭⎣⎦∑整理得 ()()()0dx xx m r k h m r m k +∞*=-=-∑ 0,1,2,k= 等价于()()()()()0dx xx xx m r k h m r k m h k r k +∞==-=*∑ 0,1,2,k=此式称为维纳-霍夫（W ｉen ｅr －Ｈolf)方程。

解此方程可得到最优权系数012,,,h h h ，此式是Wiener 滤波器的一般方程。

定义121M M h h h h ⨯⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥⎣⎦()()()1011xd xd xd xd M r r R r M ⨯⎡⎤⎢⎥⎢⎥=⎢⎥⎢⎥-⎢⎥⎣⎦()()()()()()()()()011102120xx xx xx xxxx xx xx xx xx xx M Mr r r M r r r M R r M r M r ***⨯⎡⎤-⎢⎥-⎢⎥=⎢⎥⎢⎥--⎢⎥⎣⎦则维纳-霍夫方程可写成矩阵形式xd xx R R h =求逆,得1xx xd h R R -=此式表明，已知期望信号与观测数据的互相关函数及观测数据的自相关函数时,可以通过矩阵求逆运算，得到维纳滤波器的最佳解。

3. 实验结果及分析（1)当L=2００,N=6信噪比为20ｄB 的滤波效果信噪比为１０dＢ的滤波效果信噪比为6dB的滤波效果（2)当L＝2０0，N=60信噪比为20ｄB的滤波效果信噪比为1０dＢ的滤波效果信噪比为6dB的滤波效果（３)当L=６00，Ｎ＝6信噪比为20dB的滤波效果信噪比为10dＢ的滤波效果信噪比为６dB的滤波效果实验分析:别取信号长度为200、600，滤波器长度为6、６0,加噪信噪比为20dB、10dB、6dB，组合进行实验。

每组实验得到的最小均方误差统计如下表。

由此表可以看出，信号长度越长，最小均方误差（绝对值）越大，精度越差；在信噪比较大（误差影响较小）的滤波过程中,滤波器长度约长,最小均方误差（绝对值)越小,精度越好。

表1 最小均方误差统计表加噪信号信噪比L=2０0N＝6L=200Ｎ＝6０L=600N＝61()x n20dB －178５．３971-１574.１1０7-6８8８.41032()x n10dB -845.５714 -１727.5260-1021９．33973()x n６ｄB -1213.636７-2１3５.459６－77４3．9３5８对于相同信号和滤波器(这里取L=２00,N=6),信噪比越大，最小距离误差约小;而当信噪比较小时,信号与噪声值接近,导致滤波效果受到影响,最小距离误差变大。

4. 源代码clｅaｒ；clc；％初始化变量L＝200；%信号长度N=6;%滤波器的阶次a=０.96;ｗn=randｎ（L,1）;%wn为用于生成信号的噪声信号,随机生成一个Ｌ*1矩阵,生成矩阵的元素值在％区间（０.０，１．0）之间sn=ｚｅroｓ(L,１）;%sｎ为信号，生成一个Ｌ＊1的零矩阵ｈn=zeros（N,1)；%hｎ为系统单位脉冲响应生成一个N*１的零矩阵rxx=zｅros(N，１）;%rxx为自相关函数,生成一个N＊1的零矩阵rｘd=zeroｓ(Ｎ，1);％ｒｘd为互相关函数,生成一个N＊1的零矩阵ｙn＝zeros(L,1）;%ｙn为输出信号,生成一个L*１的零矩阵xt=ｚerｏs（L＋N,1）；%生成一个（L+Ｎ)＊１的零矩阵gn=zｅｒos（L,1);%gn为yn与sn最小距离误差信号，生成一个L*1的零矩阵％根据给定公式ｓ（ｎ）=aｓ(n-１）+w(n）,生成理想信号for i=2:Ｌsn(ｉ,１)=a＊ｓn(i-1,1)+ｗn（i,1);ｅｎdsn(1，1)=ｗn(1,１);suｂpｌot(2,2，1)；pｌot（ｓn,'r＇）,ａxｉs([0 2０0 -10 10])，ｘlabeｌ('时间'),ｙlaｂe ｌ('幅度'),ｔｉｔle（＇sn')%生成期望信号方差ｃdcd=(ｖaｒ（wn))/(１－a＾2）；%对信号加噪ｘ1=aｗgn(sｎ,20);x2=ａｗgｎ(sn，１0)；x３=awgn(sn,6）；subploｔ(２,2，2)plot（ｘ3,＇g')，axiｓ(［0 ２00 -1０ 10]）,xｌａbeｌ('时间＇）,ｙｌab ｅl(＇幅度')，titｌe('x3')；%生成输入信号与理想信号的互相关函数,此处x1为输入信号，ｓn为期望信号for i=1:N，for ｍ=i:1:L,rxd(i，1)＝rxd(i,1)+x３（m,１）＊ｓn(m-ｉ＋1，1);enｄｅnd％生成输入信号的自相关函数foｒ i=1：N,for ｍ=i：１:L,rxx(i,１)=rxx(i，1)+x３（m，1)*x3(m-ｉ+1，1);ｅｎdend%将自相关函数生成托普勒斯矩阵rxx1=ｔoepｌｉｔｚ(rxx）;%生成逆矩阵iｒｘx=inｖ(rｘx1）;%生成滤波器系数ｈ(n)hn=irxx*rxd;ｆｏr i=1:Lxt(ｉ＋N,１)=ｘ3（i,１);eｎｄ%实现滤波for i＝１：L,for m=１:Ｎ，ｙn(ｉ,1）=ｙn(ｉ,1)＋ｘt(i+N+１-m,1)*ｈn(m,1);eｎｄend％计算最小均方误差信号eｎeｎ=０;eｎ=cｄ-(rxｄ＇)*hn;%生成最小距离误差信号gnｇn＝yｎ-sｎ；%画出滤波后的信号时域图sｕｂpｌｏt(2,2,3);plot（yn),aｘiｓ（[0 2００-1０ 1０]),xｌabｅｌ('时间＇),ｙｌabｅl（'幅度＇),titｌe('yn＇);%画出理想信号与输出信号对比图ｓｕbplot(2,2，4）;ｐlot(sｎ,'r')，aｘis([0 200 -１０ 10]),ｘlaｂｅl('时间'）,ｙｌabel （'幅度＇),title('ｓn与ｙn对比')；hold ｏn;pｌot(yｎ，＇b')，axiｓ(［0 2０0 -１0 10］)；hｏld oｆｆ;%画出最小距离误差信号图ｆigure；plｏt（gｎ）,ａxiｓ([０２0０－2 2]),ｘlaｂｅl('时间'),yｌabel('幅度'),titlｅ(＇gn'）;。