main_random.m
上传用户:bizsale
上传日期:2015-08-05
资源大小:3k
文件大小:9k
源码类别:

邮电通讯系统

开发平台:

Matlab

  1. %%%%%%%%%%%%%%%%%%%%%%%%%
  2. %无编码同步CDMA系统的多用户检测仿真
  3. %1 单用户
  4. %2 匹配滤波检测
  5. %3 线性解相关检测器
  6. %4 线性MMSE检测器
  7. %5 PIC检测器
  8. %5 SIC检测器
  9. %%%%%%%%%%%%%%%%%%%%%%%%%
  10. clear;
  11. nUser = 8; %Number of users
  12. L_total = 100; % infomation bits plus tail bits
  13. bias_dB = 10; % user power bias in dB
  14. bias_A = sqrt(10^(bias_dB/10)); % user amplitude bias 
  15. %A = [1 bias_A 1 bias_A 1]; %Amplitude of the nUsers:Diagonal matrix 
  16. A= eye(nUser);  %Euqal power
  17. A(nUser,nUser) = bias_A;
  18. for i = 6:nUser
  19.     A(i,i) = bias_A;
  20. end
  21. %rou = 0.3;  %cross-correlation  coefficient
  22. %R = rou*ones(nUser)+(1-rou)*diag(ones(nUser,1)); %normalized cross-correlation matrix 
  23. %Generate cross-correlation matrix  for Gold code
  24. gold_seq = gold_sequence([1 0 1 0 0],[1 1 1 0 1]);
  25. gold_seq = 2*gold_seq -1;
  26. L = length(gold_seq(1,:));
  27. gold_seq = 1/sqrt(L)*gold_seq;
  28. %S = gold_seq(1:nUser,:)';
  29. %R = S'*S;
  30. % Random Phase
  31. for i =1:nUser
  32.            m = randint(1,1,31);
  33.            gold(i,:) =[gold_seq(i,m+1:L)  gold_seq(i,1:m)]; 
  34. end
  35. S = gold(1:nUser,:)';
  36. R = S'*S;
  37. invR = inv(R);
  38. eMatrix = eye(nUser);
  39. e1 = eMatrix(:,1);
  40. e2 = eMatrix(:,nUser);
  41. decor1 = S*invR*e1;
  42. decor2 = S*invR*e2;
  43. errlim = 1000;% Number of bit errors to count as a stop criterior
  44. %ferrlim = 2;% Number of frame errors to count as a stop criterior
  45. infolim = 1E8; % Number of info bits to count as a stop criterior
  46. EbN0db = [0:2:10];
  47. % Clear bit error counter and frame error counter
  48.  nCal = 2; %Number of Calculate users
  49.  errs1 = zeros(1,length(EbN0db)); %single user
  50.  errs2 = zeros(1,length(EbN0db)); %multiuser: match filter 
  51.  errs3 = zeros(1,length(EbN0db)); %multiuser: linear decorrelating multiuser detector  
  52.  errs4 = zeros(1,length(EbN0db)); % linear MMSE multiuser detector
  53.  errs5 = zeros(1,length(EbN0db)); % Multistage PIC multiuser detector
  54.  errs6 = zeros(1,length(EbN0db)); % SIC multiuser detector
  55.  nferr = zeros(1,length(EbN0db));
  56.  
  57. for nEN = 1:length(EbN0db)
  58.     en = 10^(EbN0db(nEN)/10);      % convert Eb/N0 from unit db to normal numbers
  59.     sgma = 1/sqrt(2*en);  % standard deviation of AWGN noise
  60.     
  61.     
  62.    %Single User Bound with match filter detector
  63.    % nframe = 0;    % clear counter of transmitted frames
  64.   %  while (errs1(nEN)< errlim )    
  65.  %      nframe = nframe + 1;    
  66.   %     x = round(rand(1,L_total));    % info. bits         
  67.   %    en_output = 2*x-1;  % antipodal:+1/-1
  68.   %     en_output = S(:,1)*en_output;
  69.        % en_output = reshape(en_output,1,L*L_total);
  70.         
  71.        %single user communication
  72.    %    single = en_output+sgma*randn(L,L_total);    
  73.   %      single =  S(:,1)'*single;
  74.   %     xSingle = (sign(single)+1)/2; % output of a match filter detector(Single user)              
  75.   %     err1 = length(find(xSingle ~= x));  % Number of bit errors in current frame       
  76.   %     errs1(nEN) = errs1(nEN) + err1;      % Total number of bit errors for all frame             
  77. %    end %end while
  78. %    ber1(nEN) = errs1(nEN)/nframe/L_total; % Bit error rate       
  79.     
  80.    % Multiuser Match Filter Detector
  81.    nframe = 0;    % clear counter of transmitted frames
  82.    while (errs2(nEN)< errlim )    
  83.        nframe = nframe + 1;    
  84.        x = round(rand(nUser, L_total));    % info. bits         
  85.        en_output = 2*x-ones(size(x));  % antipodal:+1/-1       
  86.        %multiuser  communication
  87.       % spread_output  = 0;
  88.       % for i = 1:nUser    
  89.      %      spread_output =spread_output + gold_seq(i,:)'*A(i)*en_output(i,:)+sgma*sqrt(L)*randn(L,L_total);
  90.      % end %end nUser 
  91.       r = S*A*en_output + sgma*randn(L,L_total);
  92.      % y =  S(:,1)'*r;
  93.        y =  S(:,nUser)'*r; %strong power
  94.       xMatch = (sign(y)+1)/2; 
  95.       err2 = length(find(xMatch ~= x(1,:)));  % Number of bit errors in current frame       
  96.       %err2 = length(find(xMatch ~= x(nUser,:)));  % Number of bit errors in current frame     
  97.       errs2(nEN) = errs2(nEN) + err2; % Total number of bit errors for all frame              
  98.    end %end while
  99.     ber2(nEN) = errs2(nEN)/nframe/L_total; % Bit error rate        
  100.     
  101.     
  102.      % linear decorrelating multiuser detector 
  103.    nframe = 0;    % clear counter of transmitted frames
  104.    while (errs3(nEN)< errlim )    
  105.        nframe = nframe + 1;    
  106.        x = round(rand(nUser, L_total));    % info. bits         
  107.        en_output = 2*x-ones(size(x));  % antipodal:+1/-1       
  108.        %multiuser  communication
  109.        r = S*A*en_output + sgma*randn(L,L_total);
  110.        y =  decor1'*r;
  111.        %y =  decor2'*r;
  112.       xDecor = (sign(y)+1)/2; 
  113.       err3 = length(find(xDecor  ~= x(1,:)));  % Number of bit errors in current frame    
  114.       %err3 = length(find(xDecor  ~= x(nUser,:)));  % Number of bit errors in current frame  
  115.       errs3(nEN) = errs3(nEN) + err3; % Total number of bit errors for all frame                    
  116.    end %end while
  117.     ber3(nEN) = errs3(nEN)/nframe/L_total; % Bit error rate    
  118.     
  119.     % linear MMSE multiuser detector 
  120.    nframe = 0;    % clear counter of transmitted frames
  121.    while (errs4(nEN)< errlim )    
  122.        nframe = nframe + 1;    
  123.        x = round(rand(nUser, L_total));    % info. bits         
  124.        en_output = 2*x-ones(size(x));  % antipodal:+1/-1       
  125.        %multiuser  communication
  126.        r = S*A*en_output + sgma*randn(L,L_total);
  127.       mmse1 = S*inv(R+diag(sgma./diag(A))^2)*e1; 
  128.       y =  mmse1'*r;
  129.       %mmse2 = S*inv(R+diag(sgma./diag(A))^2)*e2; 
  130.       %y =  mmse2'*r;
  131.       xMMSE = (sign(y)+1)/2; % output of a linear MMSE multiuser detector 
  132.      err4 = length(find(xMMSE  ~= x(1,:)));  % Number of bit errors in current frame          
  133.       %err4 = length(find(xMMSE  ~= x(nUser,:)));  % Number of bit errors in current frame   
  134.       errs4(nEN) = errs4(nEN) + err4; % Total number of bit errors for all frame                    
  135.    end %end while
  136.     ber4(nEN) = errs4(nEN)/nframe/L_total; % Bit error rate    
  137.     
  138. %Multistage parallel interference cancellation
  139.    nframe = 0;    % clear counter of transmitted frames
  140.    while (errs5(nEN)< errlim )    
  141.        nframe = nframe + 1;    
  142.        x = round(rand(nUser, L_total));    % info. bits         
  143.        en_output = 2*x-ones(size(x));  % antipodal:+1/-1       
  144.        %multiuser  communication
  145.        r = S*A*en_output + sgma*randn(L,L_total);
  146.        y =  S'*r;
  147.        %ydecor =  decor1'*r;
  148.        rParall = sign(y);
  149.        rParall = sign(y-(R-eye(nUser))*A*rParall); %stage 1
  150.        rParall = sign(y-(R-eye(nUser))*A*rParall); %stage 2  
  151.        rParall = sign(y-(R-eye(nUser))*A*rParall); %stage 3  
  152.        xParall = (sign(rParall)+1)/2;         
  153.       err5 = length(find(xParall(1,:) ~= x(1,:)));  % Number of bit errors in current frame 
  154.       %err5 = length(find(xParall(nUser,:) ~= x(nUser,:)));  % Number of bit errors in current frame  
  155.       errs5(nEN) = errs5(nEN) + err5; % Total number of bit errors for all frame                    
  156.    end %end while
  157.     ber5(nEN) = errs5(nEN)/nframe/L_total; % Bit error rate   
  158.    
  159.     
  160.     % Serial interference cancellation
  161.    nframe = 0;    % clear counter of transmitted frames
  162.    while (errs6(nEN)< errlim )    
  163.        nframe = nframe + 1;    
  164.        x = round(rand(nUser, L_total));    % info. bits         
  165.        en_output = 2*x-ones(size(x));  % antipodal:+1/-1       
  166.        %multiuser  communication
  167.        r = S*A*en_output + sgma*randn(L,L_total);
  168.        y =  S'*r;
  169.        %ydecor =  decor1'*r;
  170.        rSIC = sign(y);
  171.        for k = nUser-1:-1:1
  172.            IC = 0;
  173.            for j = k+1:nUser
  174.                IC = IC+A(j,j)*R(j,k)*rSIC(j,:);
  175.            end
  176.            rSIC(k,:) = sign(y(k,:)-IC); 
  177.        end       
  178.        xSIC = (sign(rSIC)+1)/2;         
  179.       %err6 = length(find(xSIC(1,:) ~= x(1,:)));  % Number of bit errors in current frame            
  180.       err6 = length(find(xSIC(nUser,:) ~= x(nUser,:)));  % Number of bit errors in current frame           
  181.       errs6(nEN) = errs6(nEN) + err6; % Total number of bit errors for all frame                    
  182.    end %end while
  183.     ber6(nEN) = errs6(nEN)/nframe/L_total; % Bit error rate   
  184.     
  185.       % Display intermediate results in process  
  186.    fprintf('************** Eb/N0 = %5.2f db **************n', EbN0db(nEN));
  187.    fprintf('Frame size = %d. n', L_total);
  188.    fprintf('%d frames transmitted, %d frames in error.n', nframe, nferr(nEN));   
  189.    fprintf('%d bits transmitted, %d bits in error.n', nframe*L_total, errs3(nEN));
  190.    fprintf('Bit Error Rate: %8.4e n',  ber3(nEN));   
  191.    fprintf('***********************************************nn');
  192. end %end EN
  193. load   singleuserbound;
  194. figure(1)
  195. semilogy(EbN0db,ber1(1:length(EbN0db)),'b-o',...
  196.                 EbN0db,ber2(1:length(EbN0db)),'b-s',...
  197.                 EbN0db,ber3(1:length(EbN0db)),'b-^',...
  198.                 EbN0db,ber4(1:length(EbN0db)),'b-d',...
  199.                 EbN0db0,ber5(1:length(EbN0db)),'b-<',...
  200.                 EbN0db,ber6(1:length(EbN0db)),'b-v' );
  201. legend('Single User Bound','Match Filter','Decorrelating','MMSE','PIC','SIC');         
  202. xlabel('Eb/N0'); ylabel('Bit Error Rate');
  203. title('MUD');
  204. grid on;