IIRSER-LMS

clear all;
close all;
k1=1:800;%LMS迭代次数
k2=1:500;%SER迭代次数
xt1=randn(1,length(k1))*0. 2; %白噪声
dk=zeros(1,length(k1));%定义dk的数组
%---------计算未知系统的输出值-----------------------
b1=1.2; %b1,b2的初始值
b2=-0.6;
dk(1)=xt1(1); %先计算dk的前两个值
dk(2)=xt1(2)+dk(1)*b1;
for i=3:length(k1) %迭代800次计算出期望值
dk(i)=xt1(i)+dk(i-1)*b1+dk (i-2)*b2;
end
%%--------------画性能表面的等高线----------------------
a0=1;
b1=-2:0.01:2; %两个权值的取定
b2=-1:0.01:1;
Ks=zeros(length(b2),length (b1));%定义误差的均方值的数组en=zeros(1,length(k1));%定
义误差的数组
yk=zeros(1,length(k1));%定义yk的数组
yk(1)=xt1(1);%先计算yk的前两个值
yk(2)=xt1(2)+yk(1)*b1(2); for k=1:length(b1) %每取定一个不，下面开始遍历b2
for j=1:length(b2)
for
i=3:length(k1) %取800个时间点上的值计算误差的均方值
yk(i)=xt1(i)+yk(i-1)*b1(k)
+yk(i-2)*b2(j);
en(i)=dk(i)-yk(i);
end
Ks(j,k)=1/length(k1)*sum(e n.^2);
end
end
figure;contour(b1,b2,Ks,[ 0.274 0.224 0.172 0.43 0.31 0.28 0.15 0.09 0.06 0.048 ]);%画出性能表面的等高线xlabel('b1');ylabel('b2'); title('IIR-SER');
%%%%%---------LMS-------------
an0=zeros(1,length(k1)); %此处定义a0，b1,b2
bn1=zeros(1,length(k1));
bn2=zeros(1,length(k1));
M=[0.5 0 0;0 0.1 0;0 0
0.05];%M的值
w=zeros(3,length(k1));%定义3*800的矩阵，其中第一行表示权值a0，第二行表示权值b1，第二行表示权值b2
w(:,1)=[0 0 0]';%权的初始值yt=zeros(1,length(k1));
yt(1)=xt1(1);%先计算yt的前两个值
yt(2)=xt1(2)+yt(1)*w(2,2); U=zeros(3,length(k1));%定义3*800的矩阵，其中第一行表示权值xt，第二行表示权值yt，第二行表示权值yt(t-1)
U(:,1)=[xt1(1) 0 0]';%先计算
U的前两个值
U(:,2)=[xt1(2) yt(1) 0]';
for k=3:length(k1) %根据书中的表达式计算权值
U(:,k)=[xt1(k) yt(k-1) yt(k-2)]';
an0(k)=xt1(k)+w(2,k)*an0(k -1)+w(3,k)*an0(k-2);
bn1(k)=yt(k-1)+w(2,k)*bn1( k-1)+w(3,k)*bn1(k-2);
bn2(k)=yt(k-2)+w(2,k)*bn2( k-1)+w(3,k)*bn2(k-2);
yt(k)=xt1(k)+yt(k-1)*w(2,k )+yt(k-2)*w(3,k);
dsj=-2*(dk(k)-yt(k))*[an0(
k) bn1(k) bn2(k)]';
w(:,k+1)=w(:,k)-M*dsj; end
hold on;
plot(w(2,:),w(3,:),'r'); %%%%%---------SER-------------
k1=k2;
xt1=randn(1,length(k1))*0. 2; %白噪声
dk=zeros(1,length(k1));
%计算未知系统的输出值
b1=1.2;%权的初始值
b2=-0.6;
dk(1)=xt1(1);%先计算dk的前两个值
dk(2)=xt1(2)+dk(1)*b1;
for i=3:length(k1)
dk(i)=xt1(i)+dk(i-1)*b1+dk (i-2)*b2;
end
a=0.93;
av=0.51;
an0=zeros(1,length(k1));%
此处定义a0，b1,b2
bn1=zeros(1,length(k1));
bn2=zeros(1,length(k1));
M=[0.05 0 0;0 0.005 0;0 0 0.0025];%M的值
Q=[1 0 0;0 1 0;0 0 1];
w=zeros(3,length(k1));%定义3*800的矩阵，其中第一行表示权值a0，第二行表示权值b1，第二行表示
权值b2
w(:,1)=[0 0 0]';%权的初始值
yt=zeros(1,length(k1));
yt(1)=xt1(1);%先计算yt的前两个值
yt(2)=xt1(2)+yt(1)*w(2,2); U=zeros(3,length(k1));%定义3*500的矩阵，其中第一行表示权值xt，第二行表示权值yt，第二行表示权值yt(t-1)
U(:,1)=[xt1(1) 0 0]';
U(:,2)=[xt1(2) yt(1) 0]';
for k=3:length(k1)%根据书中的表达式计算权值
U(:,k)=[xt1(k) yt(k-1) yt(k-2)]';
S=Q*U(:,k);
r=a+U(:,k)'*S;
Q=1/a*(Q-1/r*S*S');
an0(k)=xt1(k)+w(2,k)*an0(k -1)+w(3,k)*an0(k-2);
bn1(k)=yt(k-1)+w(2,k)*bn1( k-1)+w(3,k)*bn1(k-2);
bn2(k)=yt(k-2)+w(2,k)*bn2( k-1)+w(3,k)*bn2(k-2);
yt(k)=xt1(k)+yt(k-1)*w(2,k )+yt(k-2)*w(3,k);
dsj=-2*(dk(k)-yt(k))*[an0( k) bn1(k) bn2(k)]';
w(:,k+1)=w(:,k)-M*av*(1-a^ (k+1))/(1-a)*Q*dsj;
end
hold on;
plot(w(2,:),w(3,:),'b'); legend('
','IIR-LMS','IIR-SER');。