%實驗4 多徑信道
%本程序中時間單位是微秒
%頻率單位為MHz
%碼速率單位是Mb/s
global dt t f df N T
close all
clear Eb_N0 Pe
N=2^13;         %采樣點數
L=8;            %每碼元的采樣點數
M=N/L           %碼元數
Rb=2;		    %碼速率是2Mb/s
Ts=1/Rb;        %碼元間隔
dt=Ts/L;        %時域采樣間隔
df=1/(N*dt)     %頻域采樣間隔
T=N*dt          %截短時間
Bs=N*df/2       %系統帶寬
alpha=0.5 	    %滾降系數

t=linspace(-T/2,T/2,N);     %時域橫坐標
f=linspace(-Bs,Bs,N)+eps;   %頻域橫坐標

figure(1)
 set(1,'Position',[10,50,500,400]) %設定窗口位置及大小
                 
%升余弦
hr1=sin(pi*t/Ts)./(pi*t/Ts);
hr2=cos(alpha*pi*t/Ts)./(1-(2*alpha*t/Ts).^2);
hr=hr1.*hr2;        %升余弦脈沖波形
HR1=abs(t2f(hr));   %取模是為了忽略時延
GT=sqrt(HR1);
GR=GT;
tao=Ts/2;
C=abs(1-0.5*exp(-j*2*pi*f*tao));%多徑信道
HR=HR1.*C

for loop1=1:20
 Eb_N0(loop1)=(loop1-1) %Eb/N0 in dB
 eb_n0(loop1)=10^(Eb_N0(loop1)/10);
 Eb=1;	
 n0=Eb/eb_n0(loop1); %信道的噪聲譜密度   
 sita=n0*Bs;         %信道中噪聲功率
 n_err=0;            %誤碼計數

for loop2=1:5  
  a=sign(randn(1,M));
  imp=zeros(1,N);                 %產生沖激序列
  imp(L/2:L:N)=a/dt;
  IMP=t2f(imp);
  n_ch=sqrt(sita)*randn(size(t)); %信道噪聲
  nr=real(f2t(t2f(n_ch).*GR));    %輸出噪聲
  sr=real(f2t(IMP.*HR))+nr;   %接收信號
  y=sr(L/2:L:N);		      %取樣時間無偏差
  aa=sign(y);                 %判決                 
  n_err=n_err+length(find(aa~=a))  
 end
 Pe(loop1)=n_err/(M*loop2)
 figure(1) 
  semilogy(Eb_N0,Pe,'b');
  eb_n0=10.^(Eb_N0/10);
  hold on
  semilogy(Eb_N0,0.5*erfc(sqrt(eb_n0)),'r');
  axis([0,9,1e-4,1])
  xlabel('Eb/N0')
  ylabel('Pe')
 end
end