QPSK通过Rayleigh信道多径衰落的Matlab仿真

参照《通信系统仿真原理与无线应用》351页例14-1在这个例子里，我们对有3条固定路径的AWGN多径信道中的QPSK系统进行BER性能仿真，并与在理想的AWGN信道（没有多径）中同样系统地BER性能进行比较……书上有比较详细的数学推导，不抄了。这个例子似乎没有考虑多普勒频移。待我继续学习下一个例子，这个也没太看懂。下面是该例子的源程序，PO、Pl、P2分别是LOS路径和两条延迟瑞利分量的相对功率级。当p0=0且delay！=0时为瑞利频率选择性衰落，delay==0时为瑞利平坦衰落。主程序scriptfile：%两径瑞利衰落信道仿真%设定默认参数NN=256; %传输符号个数tb=O.5; %一比特时间fs=lO; %每符号采样数ebnOdb=[l:2:l5];%设定Eb/NO%建立QPSK信号x=random_binary(NN,fs)+i*random_binary(NN,fs);%x为QPSK信号%输入功率和延迟pO=O; %视距LOS分量pl=2O; %第一路径分量p2=l; %第二路径分量delay=l; %按照每符号采样数决定的延迟delayO=O;delayl=O;delay2=delay;%设定复高斯(瑞利)衰减gainl=sqrt(pl)*abs(randn(l,NN)+i*randn(l,NN));gain2=sqrt(p2)*abs(randn(l,NN)+i*randn(l,NN));fork=l:NNforkk=l:fsindex=(k-l)*fs+kk;ggainl(l,index)=gainl(l,k);ggain2(l,index)=gain2(l,k);endendyl=x;fork=l:delay2y2(l,k)=yl(l,k)*sqrt(pO);endfork=(delay2+l):(NN*fs)y2(l,k)=yl(l,k)*sqrt(pO)+yl(l,k-delayl)*ggainl(l,k)+yl(l,k-delay2)*ggain2(l,k);end%匹配滤波器b=-ones(l,fs);b=b/fs;a=l;y=filter(b,a,y2);%仿真结束%UsethesemianalyticBERestimator.Thefollowingsetsupthesemi%analyticestimator.Findthemaximunmagnitudeofthecrosscorrelation%andthecorrespondinglag.[corlags]=vxcorr(x,y);cmax=max(max(abs(cor)));nmax=find(abs(cor)==cmax);timelag=lags(nmax);corrmag=cmax;theta=angle(cor(nmax));y=y*exp(-i*theta);%derotate%NoiseBWcalibrationhh=impz(b,a);ts=l/16;nbw=(fs/2)*sum(hh42);%Delaytheinput,anddoBERestimationonthelast128bits.Usemiddle%sample.Makesuretheindexdoesnotexceednumberofinputpoints.Eb%shouldbecomputedatthereceiverinput.index=(10*fs+8:fs:(NN-10)*fs+8);xx=x(index);yy=y(index-timelag+1);[n1n2]=size(y2);ny2=n1*n2;eb=tb*sum(sum(abs(y2).人2))/ny2;eb=eb/2;[peideal,pesystem]=qpsk_berest(xx,yy,ebn0db,eb,tb,nbw);figuresemilogy(ebn0db,peideal,'b*-',ebn0db,pesystem,'r+-')xlabel('Eb/N0(db)');ylabel('ProbabilityofError');gridonaxis([01410人(-10)1]);%Endofscriptfile.相关的一些调用程序(4个)：vxcorr.mfunction[c,lags]=vxcorr(a,b)%Thisfunctioncalculatestheunscaledcross-correlationof2vectorsof%thesamelength.Theoutputlength(c)islength(a)+length(b)-1.Itisa%simplifiedfunctionofxcorrfunctioninmatlabR12usingthedefinition:%c(m)=E[a(n+m)*conj(b(n))]=E[a(n)*conj(b(n-m))]a=a(:);%convertatocolumnvectorb=b(:);%convertbtocolumnvectorM=length(a);%sameaslength(b)maxlag=M-1;%maximumvalueoflaglags=[-maxlag:maxlag]';A=fft(a,2Anextpow2(2*M-1));%fftofAB=fft(b,2mextpow2(2*M-l));%fftofBc=ifft(A.*conj(B)); %corsscorrelation%Movenegativelagsbeforepositivelags.c=[c(end-maxlag+l:end,l);c(l:maxlag+l,l)];%Returnrowvectorifa,barerowvectors.[nrnc]=size(a);if(nr>nc)c=c.';lags=lags.';end%Endoffunctionfile.random_binary.mfunction[x,bits]=random_binary(nbits,nsamples)%Thisfunctiongeneratesarandombinarywaveformoflengthnbits%sampledatarateofnsamples/bit.x=zeros(l,nbits*nsamples);bits=round(rand(l,nbits));form=l:nbitsforn=l:nsamplesindex=(m-l)*nsamples+n;x(1,index)=(-1)Abits(m);endend%Endoffunctionfile.qpsk_berest.m%File:psk_berest.mfunction[peideal,pesystem]=psk_berest(xx,yy,ebn0db,eb,tb,nbw)%ebn0dbisanarrayofEb/Novaluesindb(specifiedatthereceiver%input);tbisthebitdurationandnbwisthenoiseBW%xxisthereference(ideal)input;yyisthefilteredoutput;nx=length(xx);%Forcomparisionpurposes,setthenoiseBWoftheidealreceiver%(integrateanddump)tobeequaltors/2.nbwideal=1/(2*tb);%noisebandwidthform=1:length(ebn0db)peideal(m)=0.0;pesystem(m)=0.0;%initialize%findn0andthevarianceofthenoise.ebn0(m)=10A(ebn0db(m)/10);%dBtolinearn0=eb/ebn0(m);%noisepowersigma=sqrt(n0*nbw*2);%variancesigma1=sqrt(n0*nbwideal*2);%%Multiplytheinputconstellation/signalbyascalefactorsothatinput%constellationandtheconstellations/signalattheinputtoreceive%filterhavethesameavepowera=sqrt(2*eb/(2*tb)).b=sqrt(2*eb/tb)/sqrt(sum(abs(xx).人2)/nx);d1=b*abs(xx);/r/