Matlab - Algoritmo RLS

clear all

 clf

 N=6; % Number of antennas

 SNR=20; % SIgnal to noise ratio

 mu=.01; % LMS step

 niter=50000; % Number of iterations

 Npro=10; % Number of promedies

 lambda = 0.99; % Forgetting factor RLS

 delta = 0.001; % Initialization of P(0)

 % Scenario

 a=25/180*pi; % Desired signal

 b1=125/180*pi; % Interferring signal 1

 b2=-90/180*pi; % Interferring signal 2

 g=1; % Gain of our desired signal

 s2n=10.^(-SNR/10);

 % Definition of signals

 sa=1/sqrt(N)*exp(j*pi*(0:N-1)*sin(a))';

 sb1=1/sqrt(N)*exp(j*pi*(0:N-1)*sin(b1))';

 sb2=1/sqrt(N)*exp(j*pi*(0:N-1)*sin(b2))';

 % Correlation matrices

 Rs=sa*sa';

 Rint=sb1*sb1'+sb2*sb2'+s2n*eye(N);

 Rx=Rs+Rint;

 % Projection and blocking matrices

 A=eye(N);

 A(:,1)=sa;

 % Gram-Schmidt algorithm

 Aort=zeros(size(A));

 Aort(:,1)=A(:,1)/sqrt(real((A(:,1)'*A(:,1))));

 for k1=2:N

    for k2=1:(k1-1)

     Aort(:,k1)=Aort(:,k1)+(Aort(:,k2)'*A(:,k1))*Aort(:,k2);

    end

    Aort(:,k1)=A(:,k1)-Aort(:,k1);

    Aort(:,k1)=Aort(:,k1)/sqrt(real((Aort(:,k1)'*Aort(:,k1))));

 end

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 C=Aort(:,1); % Restrictions subspace

 B=Aort(:,2:end); % Restrictions' orthogonal subspace

 P_C=C*inv(C'*C)*C';

 P_C_ort=eye(N)-P_C;

 % Optimal solution

 wopt=inv(Rx)*C*inv(C'*inv(Rx)*C)*g;

 wq=sa;

 e_gslc_LMS=zeros(Npro,niter); % Error sequence LMS

 e_gslc_RLS=zeros(Npro,niter); % Error sequence RLS

 % GSLC algorithm

 for kk=1:Npro

    Xa = kron(sa,ones(1,niter)).*kron(ones(N,1),(sign(randn(1,niter))+1i*sign((randn(1,niter)))/sqrt(2)));

    Xb1 = kron(sb1,ones(1,niter)).*kron(ones(N,1),(sign(randn(1,niter))+1i*sign((randn(1,niter)))/sqrt(2)));

    Xb2 = kron(sb2,ones(1,niter)).*kron(ones(N,1),(sign(randn(1,niter))+1i*sign((randn(1,niter)))/sqrt(2)));

    X = Xa + Xb1 + Xb2 + (randn(size(Xa)) + 1i*randn(size(Xa)))/sqrt(2)*sqrt(s2n);

    w_gslc_LMS = zeros(N,niter);

    w_gslc_RLS = zeros(N,niter);

    xa = B'*X;

    wa_LMS = ones(N-1,niter);

    wa_RLS = ones(N-1,niter);

    P = inv(delta)*eye(N-1);

    for k = 2:niter

       y_LMS = X(:,k)'*w_gslc_LMS(:,k-1);

       y_RLS = X(:,k)'*w_gslc_RLS(:,k-1);

       % LMS

       wa_LMS(:,k) = wa_LMS(:,k-1) + mu*xa(:,k)*(X(:,k)'*wq - xa(:,k)'*wa_LMS(:,k-1));

       % RLS

       z = P*xa(:,k);

       alfa = X(:,k)'*wq - xa(:,k)'*wa_RLS(:,k-1);

       ge = z/(lambda + xa(:,k)'*z);

       P = (P - ge*z')/lambda;

       wa_RLS(:,k) = wa_RLS(:,k-1) + alfa*ge; 

       w_gslc_LMS(:,k) = wq - B*wa_LMS(:,k);

       w_gslc_RLS(:,k) = wq - B*wa_RLS(:,k);

    end

    e_gslc_LMS(kk,:) = sum(abs(w_gslc_LMS - kron(wopt,ones(1,niter))).^2,1)/N;

    e_gslc_RLS(kk,:) = sum(abs(w_gslc_RLS - kron(wopt,ones(1,niter))).^2,1)/N;

 end

 figure(1)

 plot(10*log10(mean(e_gslc_LMS)),'b'), hold on

 plot(10*log10(mean(e_gslc_RLS)),'r')

 grid

 xlabel('Iterations')

 ylabel('dB')

 title('GSLC algorithm - Learning curve of LMS/RLS algorithms')