rake.m

2个用户下的Rake接收机 4径-5径多径 AWGN一级瑞利衰落信道 混合了m序列、Gold序列

%
%系统描述
%异步DS-CDMA系统
%多径慢衰落信道
%矩形码片成型滤波器
%Rake接收：最大比值合并
clear all
clc
timeofprog=cputime;%cpu的时间
K=2;%用户数
Lc=128;%扩频码长度
T=1;
Tc=T/Lc;%信号比特和扩频码周期
Q=10;%矩形码片成型滤波器
Mpmax=5;%最大多径数
Mps=Mpmax*ones(K,1);%每个用户的多径数
     N=1000;%一次发送的比特数
     %%%%%%%%%%%%%%%%%%%%%%%%% creat  las %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
c1=[1 1];c2=[-1 1];
in_num=[c1 c2];
j=1;
s=[1 1 -1 1];
for i=1:5
    %%% i=1 c=4
    %%%i=2  c=8
    %%%i=3  c=16
    %%%i=4  c=32
    %%%i=5  c=64
    %%%i=6  c=128
    %%%i=7  c=256
    %%%i=8  c=512
    %%%i=9  c=1024
    in_num=s;
    t=length(j);
    c1_dyn=zeros(1,2^(i+1));
      c2_dyn=zeros(1,2^(i+1));
        d1_dyn=zeros(1,2^(i+1));
          d2_dyn=zeros(1,2^(i+1));
    for p=1:t
    [c1 c2 d1 d2]=root_fun(in_num(p,:));
    c1_dyn(p,:)=c1;c2_dyn(p,:)=c2;
    d1_dyn(p,:)=d1;d2_dyn(p,:)=d2;
    end
   j=1:2^i;
   s=zeros(1,length(j)*4);
   time=1;
   for k=1:length(j)
 s(2*k-1,:)=[c1_dyn(time,:) c2_dyn(time,:)];
 s(2*k,:)=[d1_dyn(time,:) d2_dyn(time,:)];
 time=time+1;
 if time>t break
 end
   end
end
C_part=s;
cut_size=size(s);
    C1_part=C_part(:,1:cut_size(2)/2);%%%%%%%% 树形结构所产生的上半部分
    C2_part=C_part(:,cut_size(2)/2+1:end);%%%%%%%%%树形结构所产生的下半部分
    for i=2:2:cut_size(2)
        s(:,i)=-s(:,i);
    end
      S1_part=s(:,1:cut_size(2)/2);%%%%%%%% 树形结构所产生的上半部分
    S2_part=s(:,cut_size(2)/2+1:end);%%%%%%%%%树形结构所产生的下半部分
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
C1=C1_part(1,:);C2=C2_part(1,:);
S1=S1_part(1,:);S2=S2_part(1,:);
  

%随机产生的时延，每条路径可分离
TausSam=[11 38 60 120 132;30 39 98 111 121];





%设置信噪比SNRdBs
SNRdBs=-10:0;
nSNRdBs=length(SNRdBs);%nSNRdBs=11
for SNRID=1:nSNRdBs
    SNRdB=SNRdBs(SNRID);
    %AWGN方差的产生
    N0=1/(10^(SNRdB/10))-(K*Mps(K)-1)*1/(3*N);
%AWGN的单边功率普密度
%AWGN=加性高斯白噪声
Bandwidth=Tc/2;%匹配滤波器等效带宽
sigma2=N0*Bandwidth/2;% 经过匹配滤波器后AWGN的方差
sigma=sqrt(sigma2);
nbitsAll=0;%发送的总比特数1111111111111111111111111111111111111
 nerrorbits=zeros(Mps(1),1);%错误的比特数!!!!!!!!!!!!!!!!!!!!!!!!!!1
 nerrorbits1=zeros(Mps(1),1);
 times=0;%循环次数 !!!!!!!!!!!!!!!!!!!!!!!
 while(1)
     N=1000;%一次发送的比特数
     bits=2*rem(unidrnd(2,K,N),2)-1;%BPSK信号 
     %unidrnd(x)在(1,2,3。。。。x)中随机选择一个数
     %a=rem(x,y)如果x>y则a=x;x<=y则 a为x/y的余数
     %unidrnd(n,x,y)产生一个x*y全为n的矩阵
     Signal1=zeros(1,(N+1)*136*Q);%不包括AWGN%?????????????
      Signal2=zeros(1,(N+1)*136*Q);%%%!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     for k=1:K
         bitsChip=kron(bits(k,:),codes(k,:));%bits-->chips
        % koroneicker 内积 bitsChip-->1x5000
         ChipsSam=kron(bitsChip,ones(1,Q));%chip-->samples  
         %  ChipsSam-->1x50000
         for pathID=1:Mps(k)
          Fadings(pathID,:)=kron(raylrnd(1/sqrt(2),1,N),ones(1,Lc*Q));%产生瑞利衰落振幅
         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
          %raylrnd(R = RAYLRND(B,M,N,...) or R = RAYLRND(B,[M,N,...]) returns an
   % M-by-N-by-... array.
          TauSam=round(TausSam(k,pathID));%round(a)就是取a的四舍五入
          %TauSam(k,pathID)就是取TauSam的第k行第pathID个元素
          Signal(TauSam+1:TauSam+N*Lc*Q)=Signal(TauSam+1:TauSam+N*Lc*Q)+...
                                               Fadings(pathID,:).*ChipsSam;
     end
 end
 %Mp=[1 2 3;4 5 6;7 8 9],Mp(2)=4,Mp(a)是从第一列顺序数
 
 
 
 %第一个用户的Rake接收机
 %假设时延已经准确估计
 TausSamEst=TausSam;%TausSamEst为一个1X5阵
 Rakebranch=zeros(Mps(1),N);
 for pathID=1:Mps(1)
     TauSam=round(TausSamEst(1,pathID));
     RakeChipMF=1/2*sum(reshape(Signal(TauSam+1:TauSam+N*Lc*Q),Q,N*Lc),1)*(Tc/Q);
     %解调后的信号（振幅为解调前的1/2)再经过匹配滤波器
     %reshape(a,x,y)就是把a的一组数据转换为x*y的矩阵
     %这里转换为Q X N*Lc阵，再按列相加
     Noise=sigma*randn(1,N*Lc);
     RakeChipMF=RakeChipMF+Noise;%加上噪声
     Rakebranch(pathID,:)=codes(1,:)*reshape(RakeChipMF,Lc,N)*Tc;%解扩
 end
 
 
 
 %最大比值合并
 Rakebranchpower=Rakebranch.^2;
 for pathID=1:Mps(K)
     power=sum(Rakebranchpower(1:pathID,:),1);
     for temp=1:pathID
%     MRCPara(temp,:)=ones(1,1000);
               MRCPara(temp,:)=Rakebranchpower(temp,:)./power;
         %每条路径的最大比值合并系数
     end
     orderID=pathID;
     RakeMRC(orderID,:)=sum(MRCPara(pathID,:).*Rakebranch(1:pathID,:),1);
 end
 %EGC
 for pathID=1:Mps(K)
     Rake1(pathID,:)=sum(Rakebranch(1:pathID,:),1);
 end
     
 
 % 误码率计算
 RakeRcvSign=sign(RakeMRC);
 RakeRcvSign1=sign(Rake1);
 nbitsAll=nbitsAll+N;
 nerrorbits=nerrorbits+sum(abs(sign(RakeRcvSign-kron(bits(1,:),ones(Mps(K),1)))),2);
%   nerrorbits1=nerrorbits1+sum(abs(sign(RakeRcvSign1-kron(bits(1,:),ones(Mps(K),1)))),2);
 BER(:,SNRID)=nerrorbits/nbitsAll;
%  BER1(:,SNRID)=nerrorbits1/nbitsAll;
 [SNRdB,nbitsAll/N,nerrorbits']
%   [SNRdB,nbitsAll/N,nerrorbits1']
 
 
 
 %输出数据存储
%  fid=fopen(filename,'a+');
% %  outdata=[SNRdB,nbitsAll,nerrorbits',BER'];
%  for dataID=1:length(outdata)
%      fprintf(fid,'%f',outdata(dataID));
%  end
%  fprintf(fid,'\n');
%  fclose(fid);
 times=times+1;
 
 
 
 
 %跳出while循环条件
 if nerrorbits(Mps(1))>=100|times>=10
     break
 end
end
end
timeofprog=cputime-timeofprog;
display('Program exit!!!');


semilogy(SNRdBs,BER(5,:),'k-^')
    hold on
%     semilogy(SNRdBs,BER1(5,:),'-*')
%     hold on
%   legend('MRC','EGC')
%