global dt t df N
close all
N=2^13; %采样点数
L=32; %每码元的采样点数
M=N/L; %码元数
Rb=2; %码速率是2Mb/s , 码速率单位是Mb/s
Ts=1/Rb; %码元间隔
dt=Ts/L; %时域采样间隔 , 本程序中时间单位是微秒
df=1/(N*dt) ; %频域采样间隔 ,频率单位为MHz
T=N*dt; %截短时间
Bs=N*df/2 ; %系统带宽
Na=4; %示波器扫描宽度为4个码元
alpha=0.5; %滚降系数=0.5
t=linspace(-T/2,T/2,N); %时域横坐标
f=linspace(-Bs,Bs,N)+eps; %频域横坐标
for loop1=1:10
Eb_N0(loop1)=(loop1-1) ; %分贝值变为真值
eb_n0(loop1)=10^(Eb_N0(loop1)/10);
n_err=0; %误码计数
for loop2=1:3
a=sign(randn(1,M)); %产生M个取值随机码，服从标准正态分布（均值为零，方差为1），后对其取符号
imp=zeros(1,N);
imp(1:L:N)=a/dt; %产生冲击序列信号
G=Ts*sin(pi*f*Ts).*exp(-j*pi*f*Ts)./(pi*f*Ts); %发送滤波器冲激响应（矩形不归零脉冲）的傅式变换
s=f2t(t2f(imp).*G); %频域的方法产生数字信号
S=t2f(s); %数字基带信号的频域表达式
hr1=sin(pi*t/Ts)./(pi*t/Ts);
hr2=cos(alpha*pi*t/Ts)./(1-(2*alpha*t/Ts).^2);
hr=hr1.*hr2;
GR=abs(t2f(hr)); %接收滤波器采用升余弦滤波器（滚降系数=0.5）
yy=real(f2t(S.*GR)); %无加性噪声基带传输
y=yy(L/2:L:N); %以L/2为起点,步长L,取样点N
y11=sign(y); %接收码元
Y=t2f(yy); %输出信号的频域表达式
n_err=n_err+length(find(y11~=a)) ; %误码个数
end
Pe(loop1)=n_err/(M*loop2); %误码率
%end
for jj=1:100
EP=zeros(size(f));
P=S.*conj(S)/T; %输入信号功率谱
EP=(EP*(jj-1)+P)/jj; %累计平均
aa=30+10*log10(EP+eps); %加eps以避免除以零
end
for jj=1:100
EP=zeros(size(f));
PY=Y.*conj(Y)/T; %输出信号功率谱
EP=(EP*(jj-1)+PY)/jj; %累计平均
aaY=30+10*log10(EP+eps); %加eps以避免除以零
end
figure(1) %输入信号码形
set(1,'position',[10,50,400,300])%设定窗口位置及大小
figure(2) %输入信号的功率谱
set(2,'position',[350,50,400,300]);%设定窗口位置及大小
figure(3) %输入信号的眼图
set(3,'position',[690,50,400,300]);%设定窗口位置及大小
figure(4) %输出信号码形
set(4,'position',[10,400,400,300]);%设定窗口位置及大小
figure(5) %输出信号码形
set(5,'position',[350,400,400,300]);%设定窗口位置及大小
figure(6) %输出信号的眼图
set(6,'position',[690,400,400,300]);%设定窗口位置及大小
figure(7) %输出信号的眼图
set(7,'position',[800,350,400,300]);%设定窗口位置及大小
figure(1) %输入信号码形
plot(t,s,'b')
grid
axis([-4,+4,-1.5,1.5]);
xlabel(('t(ms)'));
ylabel('s(t)(V)');
title('输入信号码形');
figure (2) %输入信号的功率谱
plot(f,aa,'g');
grid
axis([-10,+10,-50,+50])
xlabel('f(MHZ)')
ylabel('Ps(f)(MHZ)')
title('输入信号的功率谱');
figure(3) %输入信号的眼图
hold on
grid
xlabel('t in us')
ylabel('s(t) in V')
tt=[0:dt:Na*L*dt];
for jj=1:Na*L:N-Na*L
plot(tt,s(jj:jj+Na*L));
end
title('输入信号的眼图');
figure(4) %输出信号码形
plot(t,yy,'b')
grid
axis([-4,+4,-1.5,1.5]);
xlabel('t(ms)');
ylabel('y(t)(V)');
title('抽样前输出信号码形');
figure (5) %输出信号的功率谱
plot(f,aaY,'g');
grid
axis([-10,+10,-50,+50])
xlabel('f(MHZ)')
ylabel('Ps(f)(MHZ)')
title('抽样前输出信号的功率谱');
figure(6) %输出信号的眼图
hold on
grid
axis([0,+2,-1.5,1.5]);
xlabel('t in us')
ylabel('y(t) in V')
tt=[0:dt:Na*L*dt];
for jj=1:Na*L:N-Na*L
plot(tt,yy(jj:jj+Na*L));
end
title('抽样前输出信号的眼图');
figure (7) %误码率
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
title('误码率');
