ps=1*10^2; %码速率为100Hz
a=0.8; %成形滤波器系数为0.8
B=(1+a)*ps; %中频信号处理带宽
Fs=8*10^2; %采样速率为800Hz
fc=2*10^2; %载波频率为200Hz
N=100; %仿真数据的长度
s=randi([0 1],N,1);
s_s=rectpulse(s,Fs/ps);
subplot(4,3,1),plot(s);
xlabel('t'); ylabel('幅度');title('原来的波形');
axis([0 100 -0.5 1.5]);grid on;
%产生随机数据作为原始数据,
channel=reshape(s,2,N/2);
channel1=channel(1,:);
channel2=channel(2,:); %串并变换
ps1=ps/2;
t=0:1/Fs:(N/2*Fs/ps1-1)/Fs;
%并将绝对码变换为相对码
dchannel1=ones(1,N/2);
dchannel2=ones(1,N/2);
for i=2:N/2
if channel1(i)==1
dchannel1(i)=-dchannel1(i-1);
else
dchannel1(i)=dchannel1(i-1);
end
end
for i=2:N/2
if channel2(i)==1
dchannel2(i)=-dchannel2(i-1);
else
dchannel2(i)=dchannel2(i-1);
end
end
%对相对码数据以Fs频率采样
Ads1=upsample(dchannel1,Fs/ps1);
Ads2=upsample(dchannel2,Fs/ps1);
%设计平方升余弦滤波器
n_T=[-2 2];
rate=Fs/ps1;
T=1;
Shape_b = rcosfir(a,n_T,rate,T);%figure(4);freqz(Shape_b)
%对采样后的数据进行升余弦滤波;
rcos_Ads1=filter(Shape_b,1,Ads1);
subplot(432); plot(rcos_Ads1);
xlabel('t'); ylabel('幅度');title('升余弦后的波形a');
grid on;
rcos_Ads2=filter(Shape_b,1,Ads2);
subplot(433); plot(rcos_Ads2);
xlabel('t'); ylabel('幅度');title('升余弦后的波形b');
grid on;
%产生载频信号
f1=cos(2*pi*fc*t);
f2=-sin(2*pi*fc*t);
%产生DPSK已调信号
qdpsk1=rcos_Ads1.*f1;
qdpsk2=rcos_Ads2.*f2;
qdpsk=qdpsk1+qdpsk2;
subplot(434); plot(qdpsk);
xlabel('t'); ylabel('幅度');title('调制后的波形');
grid on;
%已调信号加入高斯白噪声%
SNR=10;
noise=awgn(qdpsk,SNR) ; %产生噪音并加入到已调信号out中,信噪比10
subplot(435); plot(noise);
xlabel('t'); ylabel('幅度');title('有噪声的波形');
grid on;
%与相干载波相乘，实现相干解调
demod_mult1=noise.*f1;
demod_mult2=noise.*f2;
fc=[ps 3.1*10^2]; %过渡带
mag=[1 0]; %窗函数的理想滤波器幅度
dev=[0.01 0.01]; %纹波
[n,wn,beta,ftype]=kaiserord(fc,mag,dev,Fs); %获取凯塞窗参数
fpm=[0 fc(1)*2/Fs fc(2)*2/Fs 1]; %firpm函数的频段向量
magpm=[1 1 0 0]; %firpm函数的幅值向量
rec_lpf=firpm(n,fpm,magpm); %firpm函数返回的最优滤波器系数
%对乘法运算后的数据进行低通滤波，输出解调后的基带信号
demod_lpf1=filter(rec_lpf,1,demod_mult1);
demod_lpf2=filter(rec_lpf,1,demod_mult2);%设计IIR-LPF
subplot(436); plot(t,demod_lpf1);
xlabel('t'); ylabel('幅度');title('通过LPF后的波形');
axis([0 1 -1 1]);grid on;
subplot(437); plot(t,demod_lpf2); %画出信号通过该低通滤波器的波形
xlabel('t'); ylabel('幅度');title('通过LPF后的波形');
axis([0 1 -1 1]);grid on;
%抽样判决%
ssf1=ones(1,N/2);
ssf2=ones(1,N/2);
q=0; %设置判决门限
for i=3:N/2
if demod_lpf1(16*(i-1)-1)>q
ssf1(i-2)=1; %若sf>判决门限，说明此时码元为1
else
ssf1(i-2)=-1; %若sf<判决门限，说明此时码元为-1
end
end
for i=3:N/2
if demod_lpf2(16*(i-1)-1)>q
ssf2(i-2)=1; %若sf>判决门限，说明此时码元为1
else
ssf2(i-2)=-1; %若sf<判决门限，说明此时码元为-1
end
end
subplot(438); plot(ssf1);
xlabel('t'); ylabel('幅度');
title('抽样判决后波形');
axis([0 50 -1.5 1.5]);grid on;
subplot(439); plot(ssf2); %画出信号通过抽样判决器的波形
xlabel('t'); ylabel('幅度');
title('抽样判决后波形');
axis([0 50 -1.5 1.5]);grid on;
%再将相对码变换为绝对码
rchannel1=ones(1,N/2);
rchannel1(:,1)=channel1(:,1);
rchannel2=ones(1,N/2);
rchannel2(:,1)=channel2(:,1);
for i=2:N/2
if ssf1(i)==-ssf1(i-1)
rchannel1(i)=1 ;
else
rchannel1(i)=0 ;
end
end
for i=2:N/2
if ssf2(i)==-ssf2(i-1)
rchannel2(i)=1 ;
else
rchannel2(i)=0 ;
end
end
rchannel=ones(2,N/2);
rchannel(1,:)=rchannel1;
rchannel(2,:)=rchannel2;
rs=reshape(rchannel,N,1);%并串变换
subplot(4,3,10),plot(rs);
xlabel('t'); ylabel('幅度');title('解调之后的波形');
axis([0 100 -0.5 1.5]);grid on;
hError = comm.ErrorRate; %计算误码率
errorStats = step(hError,s,rs);
Errorrate= errorStats(1,:);
