%FQRRLS3T Problem 3.4
%
% 'ifile.mat' - input file containing:
% I - members of ensemble
% K - iterations
% s - deterministic part of signal to predict
% sigman - standard deviation of noise in signal to predict
% ND - delayer order
% N - filter order
% lambda - forgetting factor
%
% 'ofile.mat' - output file containing:
% ind - sample indexes
% MSE - mean-square error
clear all % clear memory
load ifile; % read input variables
LD=ND+1; % delayer length
L=N+1; % filter length
LD1=LD+1;
L1=L+1;
L2=L1+1; % auxiliary constants
MSE=zeros(K,1); % prepare to accumulate MSE*I
for i=1:I, % ensemble
D=zeros(LD,1); % initial delayer memory
X=zeros(L,1); % initial memory
n=randn(K,1)*sigman; % noise in signal to predict
xq2=zeros(L,1);
dq2=zeros(L,1);
Qthetaf=sparse(eye(L2));
for l=1:L,
Qthetabp(:,((1:L1)+(l-1)*L1))=sparse(eye(L1));
Qtheta(:,((1:L1)+(l-1)*L1))=sparse(eye(L1));
end
d=s(1)+n(1);
D=[d
D(1:ND)]; % initial values
normef=abs(D(LD));
for k=1:LD, % iterations
d=s(k+1)+n(k+1); % sample of signal to predict
D=[d
D(1:ND)]; % new delay vector
x=D(LD); % new input sample
AUX=[x
lambda^(1/2)*xq2];
for l=1:L,
AUX=Qtheta(:,((1:L1)+(l-1)*L1))*AUX;
end
efq1=AUX(1);
xq2=AUX(2:L1);
aux=normef;
normef=sqrt(lambda*aux^2+efq1^2);
Qthetaf(1,1)=1;
Qthetaf(1,L2)=0;
Qthetaf(L2,1)=0;
Qthetaf(L2,L2)=1;
c=[1
zeros(L,1)];
for l=L:-1:1,
c=Qthetabp(:,((1:L1)+(l-1)*L1))*c;
end
ce=[0
c];
for l=1:L,
ce(1:L1)=Qtheta(:,((1:L1)+(l-1)*L1))*ce(1:L1);
end
ce=Qthetaf*ce;
c=ce(2:L2);
alpha=c(1);
for l=1:L,
aux=alpha;
alpha=sqrt(aux^2+(c(l+1))^2);
Qthetabp(1,1+(l-1)*L1)=1;
Qthetabp(1,l+1+(l-1)*L1)=0;
Qthetabp(l+1,1+(l-1)*L1)=0;
Qthetabp(l+1,l+1+(l-1)*L1)=1;
end
AUX=[1
zeros(L1,1)];
for l=1:L,
AUX(1:L1)=Qtheta(:,((1:L1)+(l-1)*L1))*AUX(1:L1);
end
AUX=Qthetaf*AUX;
re=AUX(2:L2);
for l=1:L,
re=(Qthetabp(:,((1:L1)+(l-1)*L1)))'*re;
end
r=re(2:L1);
gammap=1;
for l=1:L,
aux=gammap;
gammap=sqrt(aux^2-(r(l))^2);
Qtheta(1,1+(l-1)*L1)=1;
Qtheta(1,l+1+(l-1)*L1)=0;
Qtheta(l+1,1+(l-1)*L1)=0;
Qtheta(l+1,l+1+(l-1)*L1)=1;
end
X=[x
X(1:N)]; % new input vector
AUX=[d
lambda^(1/2)*dq2];
for l=1:L,
AUX=Qtheta(:,((1:L1)+(l-1)*L1))*AUX;
end
eq1=AUX(1);
dq2=AUX(2:L1);
ep=eq1/gammap; % error sample
MSE(k+1)=MSE(k+1)+ep^2; % accumulate MSE*I
end
for k=LD1:(K-1), % iterations
d=s(k+1)+n(k+1); % sample of signal to predict
D=[d
D(1:ND)]; % new delay vector
x=D(LD); % new input sample
AUX=[x
lambda^(1/2)*xq2];
for l=1:L,
AUX=Qtheta(:,((1:L1)+(l-1)*L1))*AUX;
end
efq1=AUX(1);
xq2=AUX(2:L1);
aux=normef;
normef=sqrt(lambda*aux^2+efq1^2);
costhetaf=lambda^(1/2)*aux/normef;
sinthetaf=efq1/normef;
Qthetaf(1,1)=costhetaf;
Qthetaf(1,L2)=-sinthetaf;
Qthetaf(L2,1)=sinthetaf;
Qthetaf(L2,L2)=costhetaf;
c=[1
zeros(L,1)];
for l=L:-1:1,
c=Qthetabp(:,((1:L1)+(l-1)*L1))*c;
end
ce=[0
c];
for l=1:L,
ce(1:L1)=Qtheta(:,((1:L1)+(l-1)*L1))*ce(1:L1);
end
ce=Qthetaf*ce;
c=ce(2:L2);
alpha=c(1);
for l=1:L,
aux=alpha;
alpha=sqrt(aux^2+(c(l+1))^2);
costhetabp=aux/alpha;
sinthetabp=-c(l+1)/alpha;
Qthetabp(1,1+(l-1)*L1)=costhetabp;
Qthetabp(1,l+1+(l-1)*L1)=sinthetabp;
Qthetabp(l+1,1+(l-1)*L1)=-sinthetabp;
Qthetabp(l+1,l+1+(l-1)*L1)=costhetabp;
end
AUX=[1
zeros(L1,1)];
for l=1:L,
AUX(1:L1)=Qtheta(:,((1:L1)+(l-1)*L1))*AUX(1:L1);
end
AUX=Qthetaf*AUX;
re=AUX(2:L2);
for l=1:L,
re=(Qthetabp(:,((1:L1)+(l-1)*L1)))'*re;
end
r=re(2:L1);
gammap=1;
for l=1:L,
aux=gammap;
gammap=sqrt(aux^2-(r(l))^2);
costheta=gammap/aux;
sintheta=r(l)/aux;
Qtheta(1,1+(l-1)*L1)=costheta;
Qtheta(1,l+1+(l-1)*L1)=-sintheta;
Qtheta(l+1,1+(l-1)*L1)=sintheta;
Qtheta(l+1,l+1+(l-1)*L1)=costheta;
end
X=[x
X(1:N)]; % new input vector
AUX=[d
lambda^(1/2)*dq2];
for l=1:L,
AUX=Qtheta(:,((1:L1)+(l-1)*L1))*AUX;
end
eq1=AUX(1);
dq2=AUX(2:L1);
ep=eq1/gammap; % error sample
MSE(k+1)=MSE(k+1)+ep^2; % accumulate MSE*I
end
end
ind=0:(K-1); % sample indexes
MSE=MSE/I; % calculate MSE
save ofile ind MSE; % write output variables
SQJHC.rar_LMS MATLAB_nlms
版权申诉
29 浏览量
2022-07-14
03:23:09
上传
评论
收藏 92KB RAR 举报
寒泊
- 粉丝: 76
- 资源: 1万+
最新资源
- 基于SpringBoot+MySQL+Guava+Redis+RabbitMQ实现高并发商品限时秒杀系统.rar
- 22222222222222222
- ADSP-BF533的例程源码60:NORFLASH驱动源码(s29al016)
- Linux centos telnet 离线安装rpm
- ADSP-BF533的例程源码58:汉字字库源码
- newtask3-1newtask3-1newtask3-1newtask3-1
- Vue前端开发手册pdf
- 【STM32F103】STM32F103读取PAC功率传感器
- 超级好用的zotero插件,用以进行pdf文件的翻译
- ADSP-BF533的例程源码57:光标显示源码
资源上传下载、课程学习等过程中有任何疑问或建议,欢迎提出宝贵意见哦~我们会及时处理!
点击此处反馈