%在這個例子中，我們對有三條固定路徑的AWGN多徑信道中的QPSK系統進行BER性能仿真，并與在理想的AWGN信道（沒有多徑）中同樣系統的BER性能進
%行比較。為了簡化模型，我們做如下假設：
%1，信道中有三條路徑，包括一條是衰落的LOS路徑和兩條具有瑞利分量的路徑。與各路徑相應的接收功率級以及路徑間的延遲差是仿真的參數。
%2，信道的瑞利衰落僅影響發送信號的幅度，不影響瞬時相位
%在符號間隔內各多徑分量的衰減幅度是常數，并與相鄰間隔無關（不需要多普勒頻譜成型）
%4，沒有使用發射機濾波，接收機模型也是理想的積分—清除接收機。
%主程序：
%file c14_threeray.m
%
%Default parameters
%
NN=256;      %number of symbols
tb=0.5;        %bit time
fs=16;         %samples /symbol
ebn0db=[1:2:14];      %Eb/N0 vector 
%
%Establish Qpsk signals
%
x=random_binary(NN,fs)+i*random_binary(NN,fs); %QPSK signal
%
%Input powers and delays
%
p0=input('Enter p0>');
p1=input('Enter p1>');
p2=input('Enter p2>');
delay=input('Enter tau>');
delay0=0;delay1=0;delay2=delay;
%
%Set up the Complex Gaussian (Rayleigh) gains
%
gain1=sqrt(p1)*abs(randn(1,NN)+i*randn(1,NN));
gain2=sqrt(p2)*abs(randn(1,NN)+i*randn(1,NN));
for k=1:NN
    for kk=1:fs
        index=(k-1)*fs+kk;
        ggain1(1,index)=gain1(1,k);
        ggain2(1,index)=gain2(1,k);
    end
end
y1=x;
for k=1:delay2
    y2(1,k)=y1(1,k)*sqrt(p0);
end
for 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);
end
%
%Match filter
%
b=-ones(1,fs); b=b/fs; a=1;
y=filter(b,a,y2);
%
%ENd of simulation
%
%Use the semianalytic BER estimator .The following sets up the semi
%analytic estimator. Find the maximum magnitude of the cross correlation
%and the corresponding lag.
%
[cor lags]=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
%
%Noise BW calibration
%
hh=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 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);
figure
semilogy(ebn0db,peideal,'b*-',ebn0db,pesystem,'r+-')
xlabel('E_b/N_0(dB)');ylabel('Probability of Error');grid
axis([0 14 10^(-10) 1])
%END of script file