%multiuser_dectect_rayleigh.m

clear all;
close all;

%gold squence
load('ch200.mat','-mat')
N=31;
step=2000;
K=10;
S1=[1 1 1 1 1 -1 -1 -1 1 1 -1 1 1 1 -1 1 -1 1 -1 -1 -1 -1 1 -1 -1 1 -1 1 1 -1 -1];%a(n)=a(2)+a(-1)
S2=[1 -1 -1 -1 -1 1 -1 1 1 -1 1 -1 1 -1 -1 -1 1 1 1 -1 1 1 1 1 1 -1 -1 1 -1 -1 1];%a(n)=a(3)+a(2)+a(1)+a(-1)
for i=1:N
    SS(i,:)=(S1~=[S2(i:N),S2(1:i-1)]);
end
SS=SS*2-1;
S=SS(1:K,:);%K*N

C(1,:)=SS(1,:);
C(2:N,:)=rand(N-1,N)*2-1;
for i=2:N
    rou=C(i,:)*C(1:i-1,:).'./dot(C(1:i-1,:),C(1:i-1,:),2).';
    C(i,:)=C(i,:)-rou*C(1:i-1,:);
end
for i=1:N
    C(i,:)=C(i,:)/sqrt(dot(C(i,:),C(i,:)));
end
C_null=C(2:N,:);

SINR_LMS=zeros(1,step);
P1_LMS=zeros(1,step);
P_i_n_LMS=zeros(1,step);
E_ex_LMS=zeros(1,step);

SINR_RLS=zeros(1,step);
P1_RLS=zeros(1,step);
P_i_n_RLS=zeros(1,step);
E_ex_RLS=zeros(1,step);

SINR_Kalman=zeros(1,step);
P1_Kalman=zeros(1,step);
P_i_n_Kalman=zeros(1,step);
E_ex_Kalman=zeros(1,step);

SNR_db=[20,30,30,30,30,30,40,40,40,50];  %SNR (in dB)
SNR=10.^(SNR_db./10);
A=sqrt(SNR./SNR(1));%1×K

repeat_time=500;
ch_count=1;
for repeat=1:repeat_time
    repeat
    %through the channel
    b=zeros(K,step);
    b=sign(rand(K,step)-0.5);
    noiseN_step=randn(N,step);  % Noise
    noiseN_step=sqrt(10./SNR(1))*noiseN_step;
    
%     noise_real=randn(N,step);  
%     noise_imag=randn(N,step);
%     noiseN_step=noise_real+noise_imag*j;% Complex Noise
%     noiseN_step=(1/sqrt(2))*sqrt(K/SNR(1))*noiseN_step;	

    F=1.024e6;%chip rate
    fc=2e9;%fc--carry frequence
    v=3.3;%v--velorcity of MS,doppler frequence is 22Hz
    yN_step=zeros(N,step);
    for i=1:K
        %ch=rayleigh(N*step,F,fc,v);%N*step points
%         channel=(ch200(ch_count,:));
%         ch_count=ch_count+1;
        %channel=abs(ch200(ceil(rand(1)*200),:));
%         channel_temp=channel(1:step);
%         for p=1:N
%             channel(p:N:N*step)=channel_temp;
%         end
%         rand_fading	= abs(channel);
%         rand_fading	= mean(rand_fading .* rand_fading);		% sigma_u
%         channel	= channel / sqrt(rand_fading);
        channel=hdt(N*step, 2e9, 1024000, v*3.6, 0.005);
        channel=abs(channel.');
        Ab(i,:)=A(i)*b(i,:);
        y_temp=(Ab(i,:).'*S(i,:)).';%N*step
        y_temp=reshape(y_temp,1,N*step);
        y_temp=y_temp.*channel;
        y_temp=reshape(y_temp,N,step);
        yN_step=yN_step+y_temp;
    end
    %yN_step=abs(yN_step);
%     channel=abs(ch200(ceil(rand(1)*200),:));
%     y_temp=reshape(yN_step,1,N*step);
%     y_temp=y_temp.*channel;
%     yN_step=reshape(y_temp,N,step);
    yN_step=yN_step+noiseN_step;
    %yN_step=real(yN_step);


%     [P1_temp,P_i_n_temp,y,correct,E_ex_temp]=LMS_S(b,N,K,step,yN_step,S,SNR);
%     P1_LMS=P1_LMS+P1_temp;
%     P_i_n_LMS=P_i_n_LMS+P_i_n_temp;
%     E_ex_LMS=E_ex_LMS+E_ex_temp;

    [P1_temp,P_i_n_temp,y,correct,E_ex_temp]=RLS_S(b,N,K,step,yN_step,S,SNR);
    P1_RLS=P1_RLS+P1_temp;
    P_i_n_RLS=P_i_n_RLS+P_i_n_temp;
    E_ex_RLS=E_ex_RLS+E_ex_temp;
    
    [P1_temp,P_i_n_temp,y,correct,E_ex_temp]=Kalman_S(b,N,K,step,yN_step,S,SS,C_null,SNR);
    P1_Kalman=P1_Kalman+P1_temp;
    P_i_n_Kalman=P_i_n_Kalman+P_i_n_temp;
    E_ex_Kalman=E_ex_Kalman+E_ex_temp;   
end

% SINR_LMS=P1_LMS./P_i_n_LMS;
% SINR_db_LMS=10*log10(SINR_LMS);
% figure(1);
% hold on;
% plot(SINR_db_LMS,'b');
% E_ex_LMS=E_ex_LMS/repeat_time;
% E_ex_LMS=abs(E_ex_LMS);
% figure(2);
% hold on;
% plot(E_ex_LMS,'b');

SINR_RLS=abs( P1_RLS./P_i_n_RLS );
SINR_db_RLS=10*log10(SINR_RLS);
figure(1);
hold on;
plot(SINR_db_RLS,'r');
E_ex_RLS=E_ex_RLS/repeat_time;
E_ex_RLS=abs(E_ex_RLS);
figure(2);
hold on;
plot(E_ex_RLS,'r');

SINR_Kalman=P1_Kalman./P_i_n_Kalman;
SINR_db_Kalman=10*log10(SINR_Kalman);
figure(1);
hold on;
plot(SINR_db_Kalman,'g');
E_ex_Kalman=E_ex_Kalman/repeat_time;
E_ex_Kalman=abs(E_ex_Kalman);
figure(2);
hold on;
plot(E_ex_Kalman,'g');

figure(1);
title('SIR');
legend('LMS','RLS','Kalman');
figure(2);
title('EOE');
legend('LMS','RLS','Kalman');