QPSK通過Rayleigh信道多徑衰落的Matlab仿真

參照《通信系統(tǒng)仿真原理與無線應(yīng)用》351頁例14-1在這個例子里，我們對有3條固定路徑的AWGN多徑信道中的QPSK系統(tǒng)進行BER性能仿真，并與在理想的AWGN信道（沒有多徑）中同樣系統(tǒng)地BER性能進行比較……書上有比較詳細的數(shù)學推導，不抄了。這個例子似乎沒有考慮多普勒頻移。待我繼續(xù)學習下一個例子，這個也沒太看懂。下面是該例子的源程序，PO、Pl、P2分別是LOS路徑和兩條延遲瑞利分量的相對功率級。當p0=0且delay！=0時為瑞利頻率選擇性衰落，delay==0時為瑞利平坦衰落。主程序scriptfile：%兩徑瑞利衰落信道仿真%設(shè)定默認參數(shù)NN=256; %傳輸符號個數(shù)tb=O.5; %一比特時間fs=lO; %每符號采樣數(shù)ebnOdb=[l:2:l5];%設(shè)定Eb/NO%建立QPSK信號x=random_binary(NN,fs)+i*random_binary(NN,fs);%x為QPSK信號%輸入功率和延遲pO=O; %視距LOS分量pl=2O; %第一路徑分量p2=l; %第二路徑分量delay=l; %按照每符號采樣數(shù)決定的延遲delayO=O;delayl=O;delay2=delay;%設(shè)定復(fù)高斯(瑞利)衰減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);%仿真結(jié)束%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.相關(guān)的一些調(diào)用程序(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/