main_random.m

多用户检测算法的比较

%%%%%%%%%%%%%%%%%%%%%%%%%
%无编码同步CDMA系统的多用户检测仿真
%1 单用户
%2 匹配滤波检测
%3 线性解相关检测器
%4 线性MMSE检测器
%5 PIC检测器
%5 SIC检测器
%%%%%%%%%%%%%%%%%%%%%%%%%


clear;
nUser = 8; %Number of users
L_total = 100; % infomation bits plus tail bits
bias_dB = 10; % user power bias in dB
bias_A = sqrt(10^(bias_dB/10)); % user amplitude bias 
%A = [1 bias_A 1 bias_A 1]; %Amplitude of the nUsers:Diagonal matrix 
A= eye(nUser);  %Euqal power
A(nUser,nUser) = bias_A;
for i = 6:nUser
    A(i,i) = bias_A;
end
%rou = 0.3;  %cross-correlation  coefficient
%R = rou*ones(nUser)+(1-rou)*diag(ones(nUser,1)); %normalized cross-correlation matrix 

%Generate cross-correlation matrix  for Gold code
gold_seq = gold_sequence([1 0 1 0 0],[1 1 1 0 1]);
gold_seq = 2*gold_seq -1;
L = length(gold_seq(1,:));
gold_seq = 1/sqrt(L)*gold_seq;
%S = gold_seq(1:nUser,:)';
%R = S'*S;

% Random Phase
for i =1:nUser
           m = randint(1,1,31);
           gold(i,:) =[gold_seq(i,m+1:L)  gold_seq(i,1:m)]; 
end
S = gold(1:nUser,:)';
R = S'*S;

invR = inv(R);
eMatrix = eye(nUser);
e1 = eMatrix(:,1);
e2 = eMatrix(:,nUser);
decor1 = S*invR*e1;
decor2 = S*invR*e2;

errlim = 1000;% Number of bit errors to count as a stop criterior
%ferrlim = 2;% Number of frame errors to count as a stop criterior
infolim = 1E8; % Number of info bits to count as a stop criterior

EbN0db = [0:2:10];

% Clear bit error counter and frame error counter
 nCal = 2; %Number of Calculate users
 errs1 = zeros(1,length(EbN0db)); %single user
 errs2 = zeros(1,length(EbN0db)); %multiuser: match filter 
 errs3 = zeros(1,length(EbN0db)); %multiuser: linear decorrelating multiuser detector  
 errs4 = zeros(1,length(EbN0db)); % linear MMSE multiuser detector
 errs5 = zeros(1,length(EbN0db)); % Multistage PIC multiuser detector
 errs6 = zeros(1,length(EbN0db)); % SIC multiuser detector
 nferr = zeros(1,length(EbN0db));

 
for nEN = 1:length(EbN0db)
    en = 10^(EbN0db(nEN)/10);      % convert Eb/N0 from unit db to normal numbers
    sgma = 1/sqrt(2*en); 	% standard deviation of AWGN noise
    
    
   %Single User Bound with match filter detector
   % nframe = 0;    % clear counter of transmitted frames
  %  while (errs1(nEN)< errlim )    
 %      nframe = nframe + 1;    
  %     x = round(rand(1,L_total));    % info. bits         
  %    en_output = 2*x-1;  % antipodal：+1/-1
  %     en_output = S(:,1)*en_output;
       % en_output = reshape(en_output,1,L*L_total);
        
       %single user communication
   %    single = en_output+sgma*randn(L,L_total);    
  %      single =  S(:,1)'*single;
  %     xSingle = (sign(single)+1)/2; % output of a match filter detector(Single user)              
  %     err1 = length(find(xSingle ~= x));  % Number of bit errors in current frame       
  %     errs1(nEN) = errs1(nEN) + err1;      % Total number of bit errors for all frame             
%    end %end while
%    ber1(nEN) = errs1(nEN)/nframe/L_total; % Bit error rate       
    
   % Multiuser Match Filter Detector
   nframe = 0;    % clear counter of transmitted frames
   while (errs2(nEN)< errlim )    
       nframe = nframe + 1;    
       x = round(rand(nUser, L_total));    % info. bits         
       en_output = 2*x-ones(size(x));  % antipodal：+1/-1       
       %multiuser  communication
      % spread_output  = 0;
      % for i = 1:nUser    
     %      spread_output =spread_output + gold_seq(i,:)'*A(i)*en_output(i,:)+sgma*sqrt(L)*randn(L,L_total);
     % end %end nUser 
      r = S*A*en_output + sgma*randn(L,L_total);
     % y =  S(:,1)'*r;
       y =  S(:,nUser)'*r; %strong power
      xMatch = (sign(y)+1)/2; 
      err2 = length(find(xMatch ~= x(1,:)));  % Number of bit errors in current frame       
      %err2 = length(find(xMatch ~= x(nUser,:)));  % Number of bit errors in current frame     
      errs2(nEN) = errs2(nEN) + err2; % Total number of bit errors for all frame              
   end %end while
    ber2(nEN) = errs2(nEN)/nframe/L_total; % Bit error rate        
    
    
     % linear decorrelating multiuser detector 
   nframe = 0;    % clear counter of transmitted frames
   while (errs3(nEN)< errlim )    
       nframe = nframe + 1;    
       x = round(rand(nUser, L_total));    % info. bits         
       en_output = 2*x-ones(size(x));  % antipodal：+1/-1       
       %multiuser  communication
       r = S*A*en_output + sgma*randn(L,L_total);
       y =  decor1'*r;
       %y =  decor2'*r;
      xDecor = (sign(y)+1)/2; 
      err3 = length(find(xDecor  ~= x(1,:)));  % Number of bit errors in current frame    
      %err3 = length(find(xDecor  ~= x(nUser,:)));  % Number of bit errors in current frame  
      errs3(nEN) = errs3(nEN) + err3; % Total number of bit errors for all frame                    
   end %end while
    ber3(nEN) = errs3(nEN)/nframe/L_total; % Bit error rate    
    
    % linear MMSE multiuser detector 
   nframe = 0;    % clear counter of transmitted frames
   while (errs4(nEN)< errlim )    
       nframe = nframe + 1;    
       x = round(rand(nUser, L_total));    % info. bits         
       en_output = 2*x-ones(size(x));  % antipodal：+1/-1       
       %multiuser  communication
       r = S*A*en_output + sgma*randn(L,L_total);
      mmse1 = S*inv(R+diag(sgma./diag(A))^2)*e1; 
      y =  mmse1'*r;
      %mmse2 = S*inv(R+diag(sgma./diag(A))^2)*e2; 
      %y =  mmse2'*r;
      xMMSE = (sign(y)+1)/2; % output of a linear MMSE multiuser detector 
     err4 = length(find(xMMSE  ~= x(1,:)));  % Number of bit errors in current frame          
      %err4 = length(find(xMMSE  ~= x(nUser,:)));  % Number of bit errors in current frame   
      errs4(nEN) = errs4(nEN) + err4; % Total number of bit errors for all frame                    
   end %end while
    ber4(nEN) = errs4(nEN)/nframe/L_total; % Bit error rate    
    
%Multistage parallel interference cancellation
   nframe = 0;    % clear counter of transmitted frames
   while (errs5(nEN)< errlim )    
       nframe = nframe + 1;    
       x = round(rand(nUser, L_total));    % info. bits         
       en_output = 2*x-ones(size(x));  % antipodal：+1/-1       
       %multiuser  communication
       r = S*A*en_output + sgma*randn(L,L_total);
       y =  S'*r;
       %ydecor =  decor1'*r;
       rParall = sign(y);
       rParall = sign(y-(R-eye(nUser))*A*rParall); %stage 1
       rParall = sign(y-(R-eye(nUser))*A*rParall); %stage 2  
       rParall = sign(y-(R-eye(nUser))*A*rParall); %stage 3  
       xParall = (sign(rParall)+1)/2;         
      err5 = length(find(xParall(1,:) ~= x(1,:)));  % Number of bit errors in current frame 
      %err5 = length(find(xParall(nUser,:) ~= x(nUser,:)));  % Number of bit errors in current frame  
      errs5(nEN) = errs5(nEN) + err5; % Total number of bit errors for all frame                    
   end %end while
    ber5(nEN) = errs5(nEN)/nframe/L_total; % Bit error rate   
   
    
    % Serial interference cancellation
   nframe = 0;    % clear counter of transmitted frames
   while (errs6(nEN)< errlim )    
       nframe = nframe + 1;    
       x = round(rand(nUser, L_total));    % info. bits         
       en_output = 2*x-ones(size(x));  % antipodal：+1/-1       
       %multiuser  communication
       r = S*A*en_output + sgma*randn(L,L_total);
       y =  S'*r;
       %ydecor =  decor1'*r;
       rSIC = sign(y);
       for k = nUser-1:-1:1
           IC = 0;
           for j = k+1:nUser
               IC = IC+A(j,j)*R(j,k)*rSIC(j,:);
           end
           rSIC(k,:) = sign(y(k,:)-IC); 
       end       
       xSIC = (sign(rSIC)+1)/2;         
      %err6 = length(find(xSIC(1,:) ~= x(1,:)));  % Number of bit errors in current frame            
      err6 = length(find(xSIC(nUser,:) ~= x(nUser,:)));  % Number of bit errors in current frame           
      errs6(nEN) = errs6(nEN) + err6; % Total number of bit errors for all frame                    
   end %end while
    ber6(nEN) = errs6(nEN)/nframe/L_total; % Bit error rate   
    
      % Display intermediate results in process  
   fprintf('************** Eb/N0 = %5.2f db **************\n', EbN0db(nEN));
   fprintf('Frame size = %d. \n', L_total);
   fprintf('%d frames transmitted, %d frames in error.\n', nframe, nferr(nEN));   
   fprintf('%d bits transmitted, %d bits in error.\n', nframe*L_total, errs3(nEN));
   fprintf('Bit Error Rate: %8.4e \n',  ber3(nEN));   
   fprintf('***********************************************\n\n');
end %end EN

load   singleuserbound;
figure(1)
semilogy(EbN0db,ber1(1:length(EbN0db)),'b-o',...
                EbN0db,ber2(1:length(EbN0db)),'b-s',...
                EbN0db,ber3(1:length(EbN0db)),'b-^',...
                EbN0db,ber4(1:length(EbN0db)),'b-d',...
                EbN0db0,ber5(1:length(EbN0db)),'b-<',...
                EbN0db,ber6(1:length(EbN0db)),'b-v' );
legend('Single User Bound','Match Filter','Decorrelating','MMSE','PIC','SIC');         
xlabel('Eb/N0'); ylabel('Bit Error Rate');
title('MUD');
grid on;