1、参照通信系统仿真原理与无线应用351 页例 14-1在这个例子里,我们对有 3 条固定路径的 AWGN 多径信道中的 QPSK 系统进行 BER 性能仿真,并与在理想的 AWGN 信道(没有多径)中同样系统地 BER 性能进行比较书上有比较详细的数学推导,不抄了。这个例子似乎没有考虑多普勒频移。待我继续学习下一个例子,这个也没太看懂。下面是该例子的源程序,P0、P1、P2 分别是 LOS 路径和两条延迟瑞利分量的相对功率级。当 p0=0 且 delay!=0 时为瑞利频率选择性衰落,delay=0 时为瑞利平坦衰落。主程序 scriptfile:% 两径瑞利衰落信道仿真% 设定默认参数NN=2
2、56; % 传输符号个数tb=0.5; % 一比特时间fs=10; % 每符号采样数ebn0db=1:2:15; % 设定 Eb/N0% 建立 QPSK 信号x=random_binary(NN,fs)+i*random_binary(NN,fs); % x 为 QPSK 信号% 输入功率和延迟p0=0; % 视距 LOS 分量p1=20; % 第一路径分量p2=1; % 第二路径分量delay=1; % 按照每符号采样数决定的延迟delay0=0;delay1=0;delay2=delay;% 设定复高斯(瑞利)衰减gain1=sqrt(p1)*abs(randn(1,NN)+i*randn
3、(1,NN);gain2=sqrt(p2)*abs(randn(1,NN)+i*randn(1,NN);for k=1:NNfor kk=1:fsindex=(k-1)*fs+kk;ggain1(1,index)=gain1(1,k);ggain2(1,index)=gain2(1,k);endendy1=x;for k=1:delay2y2(1,k)=y1(1,k)*sqrt(p0);endfor k=(delay2+1):(NN*fs)y2(1,k)=y1(1,k)*sqrt(p0)+y1(1,k-delay1)*ggain1(1,k)+y1(1,k-delay2)*ggain2(1,k)
4、;end% 匹配滤波器b=-ones(1,fs);b=b/fs;a=1;y=filter(b,a,y2);% 仿真结束% Use the semianalytic BER estimator . The following sets up the semi% analytic estimator . Find the maximun magnitude of the cross correlation% and the corresponding lag .cor lags=vxcorr(x,y);cmax=max(max(abs(cor);nmax=find(abs(cor)=cmax);t
5、imelag=lags(nmax);corrmag=cmax;theta=angle(cor(nmax);y=y*exp(-i*theta); % derotate% Noise BW calibrationhh=impz(b,a);ts=1/16;nbw=(fs/2)*sum(hh.2);% Delay the input ,and do BER estimation on the last 128 bits . Use middle% sample .Make sure the index does not exceed number of input points .Eb% should
6、 be computed at the receiver input .index=(10*fs+8:fs:(NN-10)*fs+8);xx=x(index);yy=y(index-timelag+1);n1 n2=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(Prob
7、ability of Error);grid on axis(0 14 10(-10) 1);% End of script file.相关的一些调用程序(4 个):1 vxcorr.mfunction c,lags=vxcorr(a,b)% This function calculates the unscaled cross-correlation of 2 vectors of% the same length . The output length(c) is length(a)+length(b)-1. It is a% simplified function of xcorr fu
8、nction in matlabR12 using the definition:% c(m)=Ea(n+m)*conj(b(n)=Ea(n)*conj(b(n-m)a=a(:); % convert a to column vectorb=b(:); % convert b to column vectorM=length(a); % same as length(b)maxlag=M-1; % maximum value of laglags=-maxlag:maxlag;A=fft(a,2nextpow2(2*M-1); % fft of AB=fft(b,2nextpow2(2*M-1
9、); % fft of Bc=ifft(A.*conj(B); % corsscorrelation% Move negative lags before positive lags.c=c(end-maxlag+1:end,1);c(1:maxlag+1,1);% Return row vector if a,b are row vectors.nr nc=size(a);if(nrnc)c=c.;lags=lags.;end% End of function file.2 random_binary.mfunction x,bits=random_binary(nbits,nsamples
10、)% This function generates a random binary waveform of length nbits% sampled at a rate of nsamples/bit.x=zeros(1,nbits*nsamples);bits=round(rand(1,nbits);for m=1:nbitsfor n=1:nsamplesindex=(m-1)*nsamples+n;x(1,index)=(-1)bits(m);endend% End of function file.3 qpsk_berest.m% File: psk_berest.mfunctio
11、n peideal,pesystem=psk_berest(xx,yy,ebn0db,eb,tb,nbw)% ebn0db is an array of Eb/No values in db (specified at the receiver%input); tb is the bit duration and nbw is the noise BW % xx is the reference (ideal) input; yy is the filtered output;nx=length(xx);% For comparision purposes , set the noise BW
12、 of the ideal receiver% (integrate and dump) to be equal to rs/2.nbwideal=1/(2*tb); % noise bandwidthfor m=1:length(ebn0db) peideal(m)=0.0; pesystem(m)=0.0; %initialize% find n0 and the variance of the noise.ebn0(m)=10(ebn0db(m)/10); % dB to linearn0=eb/ebn0(m); % noise powersigma=sqrt(n0*nbw*2); %v
13、ariancesigma1=sqrt(n0*nbwideal*2);% Multiply the input constellation/signal by a scale factor so that input% constellation and the constellations/signal at the input to receive% filter have the same ave power a=sqrt(2*eb/(2*tb).b=sqrt(2*eb/tb)/sqrt(sum(abs(xx).2)/nx);d1=b*abs(xx);d3=abs(yy);peideal(m)=sum(q(d1/sigma1);pesystem(m)=sum(q(d3/sigma);endpeideal=peideal/nx;pesystem=pesystem/nx;% End of function file.4 q.m% File: q.mfunction out=q(x)out=0.5*erfc(x/sqrt(2);% End of function file