LMS,NLMS,RLS算法程序
%%%%%%%%%%%%%%%%%%%% LMS程序 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
ntr=9000;
s=sign(2*rand(1,ntr)-1);%生成双极性码
b = cos(pi/5)*[0.2 0.7 0.9];
% b= [0.005 0.0092 -0.0242 0.8542 -0.2182 0.0492 -0.0162];
a = 1;
% Transmission channel filter传输通道滤波器
channel = dfilt.df2t(b,a);
%信道输出信号x(n)
x=filter(b,a,s);
%增加信道噪声
snr=30;%信噪比
u=awgn(x,snr,'measured');% Add white Gaussian noise.即滤波器的输入向量
%Trained Equalization
M=11; %滤波器的抽头数
% d is training signal is a shifted version of the original transmitted
% d是训练信号
% signal s[n].
d = [zeros(1,M) s(1:ntr-M+1)];
mu=0.03; %收敛因子u
w=zeros(1,M);%滤波器系数初始化
%LMS算法
for i=M:ntr
y(i)=w*u(i:-1:i-M+1)';
e(i)=d(i)-y(i);
w=w+2*mu*e(i)*u(i:-1:i-M+1);
end
figure(1)
subplot(4,1,1);stem(u(5000:5500),'r');title('均衡前信号u(n)');axis([1 50 -2 2]);
subplot(4,1,2);stem(y(5000:5500),'r');title('均衡后输出信号y(n)');axis([1 60 -2 2]);
subplot(4,1,3);stem(d(5000:5500),'r');title('期望输出d(n)');axis([1 50 -2 2]);
subplot(4,1,4);stem(e(5000:5500),'r');title('e(n)');axis([1 50 -0.5 0.5]);
figure (2)
%subplot(2,1,1)
plot(abs(e).^2);axis([-100 6000 -0.2 2]);
% semilogy(abs(e).^2);
xlabel('迭代数量');
ylabel('均方误差')
title('LMS算法 SNR=30 u=0.03 M=11');
%%%%%%%%%%%%%%%%% NLMS程序 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all
close all
ntr=9000;
s=sign(2*rand(1,ntr)-1);%生成双极性码
b = cos(pi/5)*[0.2 0.7 0.9];
a = 1;
channel = dfilt.df2t(b,a);
snr=30;%信噪比
x=filter(b,a,s);
u=awgn(x,snr);%在信号x中加入高斯白噪声
%在系统中插入一个均衡器，抵消系统的码间串扰
%NLMS算法
% d is training signal is a shifted version of the original transmitted
% signal s[n].
d = [zeros(1,10) s(1:ntr-10)];%期望响应
%Trained Equalization
M=9; %滤波器的抽头数
c=0.001;
mu=0.2; %收敛因子
w=zeros(1,M);%滤波器系数初始化
for i=M:ntr
y(i)=w*u(i:-1:i-M+1)';
e(i)=d(i)-y(i);
w=w+2.*mu.*e(i).*u(i:-1:i-M+1)/(c+u(i:-1:i-M+1)*u(i:-1:i-M+1)');
end
%训练误差
figure(1)
subplot(3,1,1);stem(u(5000:5500));title('均衡前信号u(n)');axis([1 60 -1.5 1.5]);%由u(n)可看出系统中已发生码间串扰
subplot(3,1,2);stem(d(5000:5500));title('期待输出信号d（n）');axis([1 60 -1.5 1.5]);
subplot(3,1,3);stem(y(5000:5500));title('均衡后信号y（n）');axis([1 60 -1.5 1.5]);
figure(2)
plot(abs(e).^2);axis([-100 6000 -0.2 2.5]);
xlabel('迭代数量');%Number of iterations
ylabel('均方误差')
title('NLMS算法 SNR=30 u=0.2 M=9')
%%%%%%%%% RLS程序 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clc
clear all
close all hidden
%% 初始化参数设置（用户可根据需要修改此处参数）
nFrame = 1; % 仿真发送的帧数（可作为蒙特卡洛仿真次数）
ntr = 8000; % 每帧的发送
的数据点数
step = 0.001; % 学习步长
snr=30;
%dB = 25; % 信噪比
e = zeros(1,ntr); % 瞬时计算误差
E = zeros(1,ntr); % 统计均方误差（为平滑用）
%h = [0.18,0.3,1,0.18]; %ISI信道参数
%h = h/norm(h);
M = 11; % 均衡器阶数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
d=sign(2*rand(1,ntr)-1);
LP = ntr-M % 移掉几个数据避免负或零下标
X = zeros(M+1,LP ); % 组织均衡器输入矩阵
b = cos(pi/5)*[0.2 0.7 0.9];
a = 1;
channel = dfilt.df2t(b,a)
%snr=30;
x=filter(b,a,d);
u=awgn(x,snr);
% 循环开始（迭代次数=nFrame*nPacket）
%for kk = 1:nFrame % 此处修改蒙特卡洛次数
W = zeros(M+1,1);
for i=1:LP
X(:,i)=u(i+M:-1:i).';
end
%% RLS 迭代算法%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%初始化RLS相关参数%
lambda = 0.9; %遗忘因子,最大为1
beta = std(u)^2; %信号功率
P = 1/beta*eye(M+1,M+1); %初始化相关矩阵
% 开始均衡器参数调整
for i = 1:LP
K = (lambda^(-1)*P*X(:,i))/(1+lambda^(-1)*(X(:,i))'*P*X(:,i)); %计算增益矢量
e(i) = W'*X(:,i)-d(i); %误差
W = W - conj(K)*e(i); %权值更新
% P = lambda^(-1)*P-lambda^(-1)*K*(X(:,i))'*P; %相关矩阵更新
E(i) = E(i)+e(i)^2; %均方误差值
end
% 计算误码
y = W'*X; % 均衡后信号
figure(1)
subplot(3,1,1)
stem(u(300:350),'r');axis([1 50 -2 2]);
title('u(n)均衡前信号')
subplot(3,1,2)
stem(d(300:350),'r');axis([1 50 -2 2]);
title('d(n)期望信号')
subplot(3,1,3)
stem(y(300:350),'r');axis([1 50 -2 2]);
title('y(n)均衡后信号')
figure(2)
plot(abs(e).^2);axis([-200 8000 -0.2 3]);
title('误差平方（收敛曲线）')