si.m

数字信号的传输：1、实现滚降系数分别为0

%数字基带传输系统 %本程序中时间单位是微秒 %频率单位为MHz
%码速率单位是Mb/s clear all global dt t f df N T %全局变量 
close all;
clear all;
N=2^12; %采样点数
L=32; %每码元的采样点数 
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; %滚降系数=0.5
Na=4; %示波器扫描宽度为4个码元
t=linspace(-T/2,T/2,N); %时域横坐标
f=linspace(-Bs,Bs,N)+eps; %频域横坐标
db=input('请选择信噪比[0-15]:'); 
if isempty(db) || (db<0) || (db>15)
    db=15; 
end;
x1=sin(pi*t/Ts)./(pi*t/Ts);
x2=cos(alpha*pi*t/Ts)./(1-(2*alpha*t/Ts).^2);
x=x1.*x2;
X=abs(t2f(x,dt)); 
GT=sqrt(X); 
GR=GT;
EP=zeros(size(f))+eps;
EPr=zeros(size(f))+eps;
count=100;
sa=zeros(count,N);
sr=zeros(count,N);
Eb_N0=zeros(1,20);
eb_n0=zeros(1,20);
Pe=zeros(1,20);

for loop1=1:20
    Eb_N0(loop1)=(loop1-1);%分贝值变为线性值 
    eb_n0(loop1)=10^(Eb_N0(loop1)/10);
end

Eb=1;
for loop1=1:20
    n0=Eb/eb_n0(loop1); %信道噪声谱密度
    sita=n0*Bs; %噪声功率 
    n_err=0; %误码计数
    for ii=1:count 
        code=sign(randn(1,M)); 
        imp=zeros(1,N); %产生冲激序列 
        imp(L/2:L:N)=code/dt; 
        IMP=t2f(imp,dt); 
        Sa=IMP.*GT; %升余弦信号的傅氏变换 
        sa(ii,:)=f2t(Sa,dt); %升余弦信号的时域波形 
        sa(ii)=real(sa(ii));
        P=Sa.*conj(Sa)/T; %升余弦信号的功率谱 
        EP=(EP*(ii-1)+P)/ii;
    n_ch=sqrt(sita)*randn(size(t));%信道噪声 
    nr=real(f2t(t2f(n_ch,dt).*GR,dt));%输出噪声 
    Sr=Sa.*GR; %接收信号频谱（理想情况）
    sr(ii,:)=real(f2t(Sr,dt))+nr;%接收信号
    Sr=t2f(sr(ii,:),dt); %接收信号频谱（有噪声）
    y=sign(sr(ii,L/2:L:N)); %抽样判别
    n_err=n_err+length(find(y~=code));%误码数
    Pr=Sr.*conj(Sr)/T; %平均功率 
    EPr=(EPr*(ii-1)+Pr)/ii;
    end
    tt=0:dt:Na*L*dt;
    if loop1==db+1 
        for jj=1:Na*L:N-Na*L
            figure(1) 
            set(1,'position',[20,40,930,620]) 
            hold on; 
            subplot(2,2,3); 
            grid on;
            for kk=1:count
            plot(tt,sa(kk,jj:jj+Na*L));%发送眼图
            end
            title('发送眼图');
            xlabel('t (us)') 
            ylabel('s(t) (V)')
            axis([0,2,-2.3,2.3]); 
        end
        for jj=1:Na*L:N-Na*L
            figure(2) 
            set(2,'position',[20,40,930,620]) 
            hold on; 
            subplot(2,2,3); 
            grid on; 
            for kk=1:count
            plot(tt,sr(kk,jj:jj+Na*L));%接收眼图
            end
            title('接收眼图');
            xlabel('t (us)')
            ylabel('s(t) (V)') 
            axis([0,2,-2.3,2.3]);
        end
    end
    
 if loop1==db+1,
     s=zeros(1,N);
     s=reshape(code(ones(1,L),:),1,L*M); 
     yo=zeros(1,N); 
     yo=reshape(y(ones(1,L),:),1,L*M);
     figure(1) 
     subplot(2,2,1);
     plot(t,s,'LineWidth',1.5);%发送码型 
     grid on; 
     axis([-T/10,+T/10,1.2*min(s),1.2*max(s)])
     xlabel('t (us)') 
     ylabel('s(t) (V)') 
     title('发送码形');
     subplot(2,2,2);
%      plot(t,sa,'LineWidth',1.5);%生成波形
     plot(t,sa,'b');%生成波形 
     grid on; 
     axis([-T/20,T/20,1.2*min(sa(1,:)),1.2*max(sa(1,:))]); 
     xlabel('t (us)')
     ylabel('s(t) (V)') 
     title('发送信号波形'); 
     subplot(2,2,4)
     plot(f,30+10*log10(EP),'LineWidth',1.5); %发送信号功率谱
     grid on; 
     axis([-Bs/10,Bs/10,1.2*min(30+10*log10(EP)),1.2*max(30+10*log10(EP))]);
     xlabel('f (MHz)')
     ylabel('Ps(f) (dBm/MHz)') 
     title('发送信号功率谱密度');
     
     figure(2) 
     subplot(2,2,1);
     plot(t,yo,'LineWidth',1.5);%接收码
     grid on; 
     axis([-T/10,+T/10,1.2*min(yo),1.2*max(yo)])
     xlabel('t (us)') 
     ylabel('s(t) (V)') 
     title(['接收码形 信噪比=',num2str(db),'dB']);
     subplot(2,2,2);
%      plot(t,sr,'LineWidth',1.5);%接收波形 
     plot(t,sr,'b');%接收波形 
     grid on; 
     axis([-T/20,T/20,1.2*min(sr(1,:)),1.2*max(sr(1,:))]);
     xlabel('t (us)')
     ylabel('s(t) (V)') 
     title('接收信号波形');
     subplot(2,2,4)
     plot(f,30+10*log10(EPr),'LineWidth',1.5); %接收信号功率谱
     grid on; 
     axis([-Bs/10,Bs/10,1.2*min(30+10*log10(EPr+eps)),1.2*max(30+10*log10(EPr+eps))]); 
     xlabel('f (MHz)') 
     ylabel('Ps(f) (dBm/MHz)') 
     title('接收信号功率谱密度');
end;
     
     Pe(loop1)=n_err/(M*ii)+eps; %平均误码率
     figure(3)
     semilogy(Eb_N0,Pe,'g','LineWidth',1.5);
     eb_n0=10.^(Eb_N0/10); %还原为真值 
     hold on;
     semilogy(Eb_N0,0.5*erfc(sqrt(eb_n0)),'b','LineWidth',1.5); %理论误码率曲线 
     legend('绿色：实际误码率曲线','蓝色：理想误码率曲线'); 
     grid on;
     axis([0,9,1e-4,1]) 
     title('误码率曲线'); 
     xlabel('Eb/N0 (dB)')
     ylabel('Pe') 
end;