% Illustration of the conventional RLS algorithm
% using an adaptive equalization experiment
% 
%
% Programmed by: Dimitris Manolakis, 1999
%
%-----------------------------------------------------------
% Copyright 2000, by Dimitris G. Manolakis, Vinay K. Ingle,
% and Stephen M. Kogon.  For use with the book
% "Statistical and Adaptive Signal Processing"
% McGraw-Hill Higher Education.
%-----------------------------------------------------------


close all;
W=2.9;
Nexp=10;
N=1000;
Nmc=1; % Number of ensemble realizations
M=11;
lambda=0.99;
varv=0.001;
h=zeros(3,1);
er=zeros(N,Nmc);

h(1)=0.5*(1+cos(2*pi*(1-2)/W));
h(2)=0.5*(1+cos(2*pi*(2-2)/W));
h(3)=0.5*(1+cos(2*pi*(3-2)/W));

% Learning curves

hc=[0 h(1) h(2) h(3)]';
n0=7;

t=(1:N)';

for i=1:Nmc

y=sign(rand(N,1)-0.5);
v=sqrt(varv)*randn(N,1);
x=filter(hc,1,y)+v;
x=[zeros(M-1,1);x];
yd=zeros(N+M-1,1); e=yd;
yd(n0+M-1:N+M-1)=y(1:N-n0+1);

%------------------------------------------------
% CRLS Algorithm
%------------------------------------------------

% Initialization

lambda=0.98;
P=(10^-3)*eye(M,M);
c=zeros(M,1);g=c;glambda=g;

% Steady-state

for n=M:M+N-1
     xn=flipud(x(n-M+1:n));
     glambda=P*xn;
     alphal=lambda+conj(glambda')*xn;
     g=glambda/lambda;a(n)=1-conj(g')*xn;
     P=(P-g*conj(glambda'))/lambda;
     P=(P+P')/2;
     e(n)=yd(n)-conj(c')*xn;
     c=c+g*conj(e(n));
end

%------------------------------------------------

end

eplot=e(M:M+N-1).^2;
subplot(2,1,1), plot(t,abs(eplot))
ylabel('|e(n)|^2');
xlabel('n');
subplot(2,1,2), plot(t,a(M:M+N-1));
ylabel('\alpha(n)');
xlabel('n');