现代信号处理1

1、现代信号处理 Assignment4姓名：饶康 班级：研1305 学号：2013020141题目1：考虑如下图权值线性组合器，输入端引入随机信号，其平均功率为；假设信号随机抽样相互独立，取。编程实现：（1） 画出LMS算法性能曲面等值线，要求等值线权值间隔不超过1，标明坐标值、均方误差值和性能最小点位置及最小均方误差值，分别对应初始权值和绘出加权值收敛轨迹，迭代次数不小于100次；（2） 计算和时学习曲线的时间常数，绘出学习曲线并在学习曲线中观测时间常数，与理论计算值比较；（3） 计算和时的失调并比较；（4） 分析比较的大小对自适应滤波的影响。要求：写出实验报告：包括原理、方法和结果，并附源代

2、码（加必要的注释）和仿真数据结果。解：（1）LMS 自适应滤波器是使滤波器的输出信号与期望响应之间的误差的均方值为最小，因此称为最小均方（LMS）自适应滤波器。如题图，LMS自适应滤波器是有线性组合器构成，输入为x(k)和x(k-1)，输出y(k)是这些输入加权后的线性组合，即yk=w0xk+w1xk-1 (1)定义权向量为W=w0 w1T，输入信号向量为X=xk xk-1，题图中dk代表期望响应，定义误差信号为k=dk-yk 2写成向量形式为k=dk-XW (3)误差平方为k2=d2k-2dkXW+WTXTXW (4)上式两边去数学期望后，得均方误差为Ek2=Ed2(k)-2EdkXW+WT

3、EXTXW (5)互相关函数行向量定义为RdxT:RdxT=Ed(k)X (6)和自相关函数矩阵Rxx=EXTX (7)则均方误差可表述为Ek2=Ed2(k)-2RdxTW+WTRxxW (8)这表明，均方误差是权系数向量W的二次函数，它是一个中间向上凹的抛物型曲面，是具有唯一最小值的函数。调节权系数使均方误差为最小，相当于沿抛物形曲面下降找最小值。将(8)式对权系数W求导，得到均方误差函数的梯度为k=-2Rdx+2RxxW (9)另k=0，即可求出最佳权系数向量Wopt=Rxx-1Rdx (10)得到最优权值即可求最小均方误差，带入（8）式可得Ek2min=Ed2(k)-RdxTWopt (

4、11)以上计算最优权值的过程需要知道自相关和互相关，还要对矩阵求逆，还是比较麻烦的，而Widrow and Hoff LMS算法给出了一种用迭代方法求最优权值近似值的方法，权值迭代公式如下Wk+1=Wk+2kXT (12)其中是一个控制收敛速度与稳定性的常数，称之为收敛因子（或增益常数）。有了迭代公式就可以在给定初始权值以及收敛因子的条件下求出最优权值的近似值。以上算法使用matlab实现，初始权值选择W=00，=0.1，进行迭代运算200次，得到权值迭代轨迹，然后改变初始权值W=4-10，=0.05，再次进行迭代并再次得到权值轨迹。接下来计算题中系统的性能曲面，输入为单频正弦波，求其自相关可

5、等效成整周期内在时间上求平均，即Exkxk=1Nk=0N-1(sin2kN)2+Erk2=0.51Exkxk-1=1Nk=0N-1sin2kNsin2(k-1)N+Erkrk-1=0.5cos(/8)再计算出期望信号和 输入信号的互相关值，Edkxk=1Nk=0N-12cos2kNsin2kN=0Edkxk-1=1Nk=0N-12cos2kNsin2(k-1)N=-sin(/8)则可写出相关矩阵如下Rxx=0.510.5cos(/8)0.5cos(/8)0.51Rdx=0-sin(/8) 最优权值为W*=Rxx-1Rdxw0*w1*=0.510.5cos(/8)0.5cos(/8)0.510-

6、sin(/8)=3.784-4.178 接着计算最小均方误差为min=Ed2(k)-RdxTW*=0.4011 那么性能函数为MSE=min+W-W*TRxx(W-W*)用matlab计算该函数并画出性能曲面的等值线，在图上标注好均方误差、性能最小点以及最小均方误差，然后将上述迭代的权值再画入图中得到如图1结果图1 性能曲面等值线图（2）学习曲线的时间常数定义为mse=14 (13)当收敛因子分别为0.05和0.1时计算时间常数，见图3。而学习曲线定义为均方误差MSE对应迭代次数的曲线，采用初始权值为W=00，收敛因子分别为0.1和0.05时计算学习曲线MSE(k)，并分别绘于图2 (a)和图

7、2 (b)中如下图2 (a) =0.05时学习曲线 图2 (b) =0.1时学习曲线从图上可以看出越小学习曲线收敛的越慢，也就是说学习曲线的时间常数越小，这与理论计算结果是一致的。（3）失调定义为M= trRxx (14)其中tr(Rxx)是自相关矩阵的迹，即矩阵主对角线元素的和。分别对应收敛因子为0.05和0.1时计算出失调如图3图3 时间常数和失调的计算（4）的取值会影响自适应滤波的收敛性，的取值不能太大。当满足0<<1max (15)时才会收敛，其中max为自相关矩阵的最大自相关值。而的大小还会影响收敛速度，较大的值会使学习曲线的时间常数更大，从而收敛更快，如图2中(a)图收

8、敛速度比(b)图快。但是过大会形成震动性的过渡特性，就像图2中(a)图的收敛值较小，学习曲线比较平滑，(b)图的值较大，就会出现振动的特性。题目2：利用matlab设计一个二阶自适应噪声对消器，对含高斯白噪声信道干扰的正弦信号进行滤波，信噪比为10dB。画出有噪信号、降噪后信号和噪声收敛结果。解：首先设计噪声对消系统模型如图4y(n)信源-d(n)噪声V0(n)AF+噪声V1(n)s(n)e(n)图4 自适应噪声对消系统模型噪声V0(n)、V1(n)、信源s(n)和y(n)是不相关的，但是噪声V0(n)和V1(n)之间是相关的。自适应滤波系统的误差为=sn+V0n-yn (16)取平方的期望后

9、E2=Esn2+EV0n-yn2 (17)当调节滤波器使均方误差最小时，信号功率Es2将不受影响，即Emin2=Esn2+Emin(V0n-yn)2 (18)所以当均方误差最小时EV0n-yn2也最小，此时便达到减弱或对消噪声的效果，此时的误差就是降噪后的信号。自适应滤波部分使用LMS算法，权值迭代的公式如题1中的迭代公式，经过matlab编程得到自适应噪声对消的效果如图5图5 自适应噪声对消附录一：题目1的matlab源代码%LMS算法close allclear allclc%设置参数，N为采样个数，u为增益常数N=16;u=0.1 0.05;%设置迭代次数kK=400;k=1:K;%计算

10、出的自行关阵Rx=0.51 0.5*cos(pi/8);0.5*cos(pi/8) 0.51;Rdx=0;-sin(pi/8);rk=0+0.1*randn(1,K);%平均功率为0.01的随机信号（选用高斯噪声）dk=2*cos(2*pi*k/N);%期望信号xk=sin(2*pi*k/N)+rk;%输入的带噪信号%设置起始权值wk1(:,1)=0;0;wk2(:,1)=4;-10;%求理论最佳权值wk_l=(Rx-1)*Rdx;%性能曲面等值线的等值e=20 15 10 5 2 1 0.8 0.6 0.5 0.42;syms w0syms w1v=w0;w1-wk_l;for i=1:2

11、%用LMS算法迭代求最佳权值 for j=2:K xkk=xk(j),xk(j-1);%输入信号 yk=xkk*wki(:,j-1);%输出信号 err=dk(j)-yk;%误差 wki(:,j)=wki(:,j-1)+2*u(i)*err*xkk'%权值迭代 endend%性能曲面最小值（最小均方误差）Emin=mean(dk.2)-Rdx'*wk_l;figure(1)for m=1:10 z=Emin+v'*Rx*v-e(m); set(ezplot(z,-2 8 -10 0),'color','blue') hold on;end

12、%画迭代时权值的变化plot(wk1(1,:),wk1(2,:),'r');hold onplot(wk2(1,:),wk2(2,:),'r');hold on%标注最佳权值的位置plot(wk_l(1),wk_l(2),'*');title('性能曲面等值线')%计算学习曲线for i=1:2 wk1(:,1)=0;0; wk2(:,1)=0;0; for j=2:K xk=sin(2*pi*j/N) sin(2*pi*(j-1)/N)+rk(j);%输入信号 yk=xk*wki(:,j-1);%输出信号 err=dk(j)-y

13、k;%误差 wki(:,j)=wki(:,j-1)+2*u(i)*err*xk'%权值迭代 end figure for j=1:K mse(j)=Emin+(wki(:,j)-wk_l)'*Rx*(wki(:,j)-wk_l); end plot(mse) title('LMS算法学习曲线') xlabel('迭代次数k') ylabel('学习曲线MSE') if i=1 legend('u=0.1') else legend('u=0.05') endend%计算时间常数t=1./(4.*u)

14、;%计算失调M=u.*trace(Rx);disp('增益常数为0.1时，时间常数为 ', num2str(t(1),',失调为 ', num2str(M(1)disp('增益常数为0.05时，时间常数为 ', num2str(t(2),',失调为 ', num2str(M(2)附录二：题目2的matlab源代码%用LMS算法设计自适应滤波器close allclear allclcT=10000;t=0:0.001:0.001*(T-1);xs=sin(2*pi*t);xn=randn(1,T);xn=xn-mean(xn);p

15、ower=1/length(xs)*sum(xs.*xs);var=power/10;xn=sqrt(var)/std(xn)*xn;%噪声，作为系统输入x=xs+xn;%带噪声信号e1=sum(xs.2);e2=sum(x-xs).2);snr=10*log10(e1/e2);d=xs+xn;%对于自适应对消器，用xs+xn作为期望信号%设计自适应滤波器N=2;%滤波器阶数w=zeros(N,1);%初始化滤波器权值rx_max=max(eig(xn*xn'); %输入信号相关矩阵的最大特征值u=abs(randn(1)*(1/rx_max); % 收敛因子y=zeros(1,T);for k=2:T y(k)=xn(k) xn(k-1)*w; e(k)=d(k)-y(k); w=w+2*u*e(k)*xn(k) xn(k-1)' % 更新权值ende3=sum(y.2);e4=su