final.m

本程序设计了一个最佳升余弦滚降系统

%例五 最佳升余弦滚降系统

%本程序中时间单位是微秒
%频率单位为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 	%滚降系数
Na=4

t=linspace(-T/2,T/2,N); %时域横坐标
f=linspace(-Bs,Bs,N)+eps;   %频域横坐标

figure(1)
 set(1,'Position',[10,50,500,400])                 
figure(2)       %设定窗口位置及大小
 set(2,'Position',[200,300,600,400])
figure(3)       %设定窗口位置及大小
 set(3,'Position',[10,200,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;
G=t2f(hr);
HR=abs(t2f(hr));%取模是为了忽略时延
GT=sqrt(HR);    %根号升余弦
GR=GT;

EPi1=zeros(size(f))+eps;
EPi2=zeros(size(f))+eps;
EPr1=zeros(size(f))+eps;
EPr2=zeros(size(f))+eps;

for loop1=1:50
 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:10 
  %************************输入信号部分*******************
  a=sign(randn(1,M));
  imp=zeros(1,N);    %产生冲激序列
  imp(L/2:L:N)=a/dt;
  IMP=t2f(imp);
  Pi1=IMP.*conj(IMP)/T;      %升余弦信号的功率谱
  EPi1=(EPi1*(loop2-1)+Pi1+eps)/loop2;    %取平均，使功率谱曲线更平滑
  
  si=real(f2t(IMP.*GR));          %输入信号经过根升余弦
  SI=t2f(si);
  Pi2=SI.*conj(SI)/T;      %升余弦信号的功率谱
  EPi2=(EPi2*(loop2-1)+Pi2+eps)/loop2;    %取平均，使功率谱曲线更平滑
  
  %************************信道及噪声*********************
  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);       %判决，大于零判1，小于零判-1   
  SR=t2f(sr);
  Pr1=SR.*conj(SR)/T;        %输出信号的功率谱；conj是取复共轭
  EPr1=(EPr1*(loop2-1)+Pr1)/loop2;      %取平均
  
  %************************构造最终输出，方波信号及功率*****
  bb=zeros(L,length(aa));
  for loop=1:L,bb(loop,:)=aa; end  %每一行都赋为aa的值，共M列，L行即每码元采样L点
  rr=reshape(bb,1,N);              %转化成1行N列，以便作为横坐标
  RR=t2f(rr);       %输出信号的傅氏变换
  Pr2=RR.*conj(RR)/T;        %输出信号的功率谱；conj是取复共轭
  EPr2=(EPr2*(loop2-1)+Pr2)/loop2;      %取平均
  
  %************************误码计数***********************
  n_err=n_err+length(find(aa~=a));  %共累加loop2次
 end
 
 Pe(loop1)=n_err/(M*loop2);  %平均每个码元的错误概率
 
  %*************************图像部分**********************
  %误码率曲线
 figure(1) 
  subplot(2,1,1)
  semilogy(Eb_N0,Pe,'r');   %纵坐标为log，横坐标为线性；共取loop1个点
  grid
  eb_n0=10.^(Eb_N0/10);
  hold on           %叠加，共loop1次
  semilogy(Eb_N0,0.5*erfc(sqrt(eb_n0)));    %计算误码率值
  axis([0,10,1e-4,1])
  xlabel('Eb/N0')
  ylabel('Pe')
  title('误码率曲线')
end

 subplot(2,1,2)      %眼图
 hold on
 grid 
 tt=[0:dt:Na*L*dt];
 for jj=1:Na*L:N-Na*L
  plot(tt,sr(jj:jj+Na*L));
 end
 title('经信道后信号的眼图')

 %输入、输出波形，取样前眼图
 figure(2)
 subplot(2,2,1)    %输入信号波形
 plot(t,imp);
 grid
 axis([0,5,-5,5])
 xlabel('t（us）')
 ylabel('Input（V）')
 title('输入（冲激序列）信号波形')
 
 subplot(2,2,2)    %输入信号经过发送滤波器波形
 plot(t,si);
 grid
 axis([0,5,-3,3])
 xlabel('t（us）')
 ylabel('Input（V）')
 title('输入经过根升余弦滤波器后的信号波形')
 
 subplot(2,2,3)     %抽样前接受信号波形
 plot(t,sr);
 grid
 axis([0,5,-3,3])
 xlabel('t（us）')
 ylabel('Output（V）')
 title('抽样前接收到的信号波形')
 
 subplot(2,2,4)     %输出信号波形
 plot(t,rr);
 grid
 axis([0,5,-3,3])
 xlabel('t（us）')
 ylabel('Output（V）')
 title('输出（方波，双极性不归零）信号波形')
 
%输入、输出信号功率谱
figure(3)
 subplot(2,2,1)     %输入窄脉冲信号功率谱
 plot(f,30+10*log10(EPi1),'g');
 grid
 axis([-3,+3,-50,50])
 xlabel('f (MHz)')
 ylabel('Ps(f)  (dBm/MHz)')
 title('输入（根号升余弦）信号功率谱')

 subplot(2,2,2)     %输入经过根号升余弦滤波器后的信号功率谱
 plot(f,30+10*log10(EPi2),'g');
 grid
 axis([-3,+3,-50,50])
 xlabel('f (MHz)')
 ylabel('Ps(f)  (dBm/MHz)')
 title('输入（根号升余弦）信号功率谱')
 
 subplot(2,2,3)      %抽样前接收到的信号功率谱
 plot(f,30+10*log10(EPr1),'g');
 grid
 axis([-5,+5,-50,50])
 xlabel('f (MHz)')
 ylabel('Ps(f)  (dBm/MHz)')
 title('抽样前接收到的信号功率谱')
 
 subplot(2,2,4)    %输出（方波，双极性不归零码）信号功率谱
 plot(f,30+10*log10(EPr2),'g');
 grid
 axis([-5,+5,-50,50])
 xlabel('f (MHz)')
 ylabel('Ps(f)  (dBm/MHz)')
 title('输出（方波，双极性不归零）信号功率谱')