bpsk_sw.m

BPSK调制程序

%BPSK的频谱分析，一个码元对应一个周期正弦波（Fs＝2000）
clear
df=0.25;
fs=20;%20MHz
fc=2;
rb=0.2;
ts=1/fs;
m=1024;%m是码元数，自己随便定义
n=40;%n是每个码元取点数，即对应的正弦波4个周期取多少个点
a=sj(m,n);
%sj(m,n)的解释参见sj.m,sj(m,20)的作用相当于上面的for循环
t=0:1:m*n-1;
Fc1=sin(8*pi*t/n);%码元0发射的载波
Fc2=sin(8*pi*t/n+pi);%码元1发射的载波
b=1-a;
u0=a.*Fc2+b.*Fc1;

%产生白噪声
snr=13;%S/N
snr_lin=10^(snr/10);%线性信噪比
signal_power=(norm(u0)^2)/length(u0);
noise_power=signal_power/snr_lin;
noise_std=sqrt(noise_power);
noise=noise_std*randn(1,length(t));
u1=u0+noise;%加白噪声
u=u1;

%求瞬时平方包络
v=blpf(u);
%均值
miu=mean(v);
%方差
delt=var(v);
%求R参数
r=delt/(miu^2);


%信号平方
sqrtu=u.*u;

%求功率谱
nfft=1024;
u2=Rx_set(u,nfft);  %这个100是取了100个自相关系数
u3=fftshift(abs(fft(u2)));%求功率谱
ff=-1*fs/2:fs/nfft:fs/2;

%求平方功率谱及共振峰个数
sqrtu2=Rx_set(sqrtu,nfft);  %这个100是取了100个自相关系数a
sqrtu3=fftshift(abs(fft(sqrtu2)));%求功率谱
ls=length(sqrtu3)
lsl=round((fs/2+2*fc-1)/fs*ls);
lsr=round((fs/2+2*fc+1)/fs*ls);
c=20;
gzf=0;
for iii=lsl:lsr
    if sqrtu3(iii)>20
        gzf=gzf+1;
    end
end
gzf
%频谱分析
[M,m,df1]=fftseq(a,ts,df);
M=M*ts;
f=[0:df1:df1*(length(m)-1)]-fs/2;
[U,u,df1]=fftseq(u,ts,df);
U=U*ts;
[SQRTU,sqrtu1,df1]=fftseq(sqrtu,ts,df);
SQRTU=SQRTU*ts;
% phase=angle(U);
%解调接口，因为解调部分是FM的，所以不能正确解调。
% [v,phase]=env_phas(u,ts,fc);
% phi=unwrap(phase);
% dem=(1/(2*pi*kf))*(diff(phi)/ts);


figure (4)
subplot(5,2,1)
plot(t,m(1:length(t)))
xlabel('time')
title('the message signal')
subplot(5,2,2)
plot(t,u(1:length(t)))
xlabel('time')
title('the modulated signal')



subplot(5,1,2)
plot(f,abs(fftshift(U))) %(f,abs(fftshift(M)))
xlabel('frequency(MHz)')
ylabel('频谱')

subplot(5,1,3)
plot(ff,u3) 
xlabel('frequency(MHz)')
ylabel('功率谱')

subplot(5,1,4)
plot(f,abs(fftshift(SQRTU)))
xlabel('frequency(MHz)')
ylabel('平方频谱')

subplot(5,1,5)
plot(ff,sqrtu3) 
xlabel('frequency(MHz)')
ylabel('平方功率谱')