通信原理与MATLAB(七):2FSK的调制解调

1. 2FSK的调制原理

2FSK调制原理如下图所示，基带码元d(t)中码元为1时，波形为频率为f1的高频载波；基带码元d(t)中码元为0时，波形为频率为f2的高频载波实现2FSK信号的调制，即基带码元和f1的高频正弦波相乘生成2ASK，基带码元的反码和f2的高频正弦波相乘生成第二个2ASK，两个2ASK相加得到2FSK。

波形图如下图所示

2. 2FSK的解调原理

2FSK的解调原理如下图所示，2FSK信号经过信道传输之后，分为上下两路经过带通滤波器变成两路2ASK信号，再和对应的载波相乘，然后经过低通滤波后抽样判决恢复出原始基带码元信号。

3. 2FSK的代码

clear all;                  % 清除所有变量
close all;                  % 关闭所有窗口
clc;                        % 清屏
%% 基本参数
M=10;                       % 产生码元数    
L=100;                      % 每码元复制L次,每个码元采样次数
Ts=0.001;                   % 每个码元的宽度,即码元的持续时间
Rb=1/Ts;                    % 码元速率1K
dt=Ts/L;                    % 采样间隔
TotalT=M*Ts;                % 总时间
t=0:dt:TotalT-dt;           % 时间
Fs=1/dt;                    % 采样间隔的倒数即采样频率
%% 产生单极性波形
wave=randi([0,1],1,M);      % 产生二进制随机码,M为码元个数
fz=ones(1,L);               % 定义复制的次数L,L为每码元的采样点数
x1=wave(fz,:);              % 将原来wave的第一行复制L次，称为L*M的矩阵
jidai=reshape(x1,1,L*M);    % 产生单极性不归零矩形脉冲波形，将刚得到的L*M矩阵，按列重新排列形成1*(L*M)的矩阵
%% 基带信号取反
jidai2=ones(1,length(jidai));% 产生和jidai相同长度的jidai2
% 基带信号取反即jidai为1的时候，jidai2为0;jidai为0的时候，jidai2为1
for n=1:length(jidai)
    if jidai(n)==1
        jidai2(n)=0;
    else
        jidai2(n)=1;
    end
end
%% 2FSK调制
fc1=2000;                   % 载波1频率2kHz       
zb1=cos(2*pi*fc1*t);        % 载波1信号
fsk1=jidai.*zb1;            % 基带信号1和载波1相乘 
fc2=10000;                  % 载波2频率10kHz       
zb2=cos(2*pi*fc2*t);        % 载波2信号
fsk2=jidai2.*zb2;           % 基带信号2和载波2相乘  
fsk=fsk1+fsk2;              % 2FSK的调制
figure(1);                  % 绘制第1幅图
subplot(211);               % 窗口分割成2*1的，当前是第1个子图 
plot(t,jidai,'LineWidth',2) % 绘制基带码元波形，线宽为2
title('基带信号波形')% 标题
xlabel('时间/s');           % x轴标签
ylabel('幅度');             % y轴标签
axis([0,TotalT,-0.1,1.1])   % 坐标范围限制

subplot(212)                % 窗口分割成2*1的，当前是第2个子图 
plot(t,fsk,'LineWidth',2);  % 绘制2FSK的波形 
title('2FSK信号波形')       % 标题
xlabel('时间/s');           % x轴标签
ylabel('幅度');             % y轴标签
axis([0,TotalT,-1.1,1.1]);  % 坐标范围限制

%% 信号经过高斯白噪声信道
fsk=awgn(fsk,20);           % 信号fsk中加入白噪声，信噪比为SNR=20dB
figure(2);                  % 绘制第2幅图
subplot(311);               % 窗口分割成3*1的，当前是第1个子图 
plot(t,fsk,'LineWidth',2);  % 绘制2FSK信号加入白噪声的波形
axis([0,TotalT,-1.5,1.5]);  % 坐标范围设置
title('通过高斯白噪声信道后的信号');% 标题
xlabel('时间/s');           % x轴标签
ylabel('幅度');             % y轴标签

%% 带通滤波器设计
%% 第一个带通滤波器设计
fp1=2*(fc1-Rb);             % 第一个带通滤波器频率1
fs1=2*(fc1+Rb);             % 第一个带通滤波器频率2
b1=fir1(30, [fp1/Fs fs1/Fs], boxcar(31));% 生成fir滤波器系统函数中分子多项式的系数
% fir1函数三个参数分别是阶数，数字截止频率，滤波器类型
% 这里是生成了30阶(31个抽头系数)的矩形窗带通滤波器
[h1,w1]=freqz(b1, 1,512);   % 生成fir滤波器的频率响应
% freqz函数的三个参数分别是滤波器系统函数的分子多项式的系数，分母多项式的系数(fir滤波器分母系数为1)和采样点数(默认)512
band1=fftfilt(b1,fsk);      % 对信号进行滤波，fsk是等待滤波的信号，b1是fir滤波器的系统函数的分子多项式系数
subplot(312)                % 窗口分割成3*1的，当前是第2个子图 
plot(t,band1,'LineWidth',2) % 绘制经过低通滤波器后的信号
axis([0,TotalT,-1.1,1.1]);  % 设置坐标范围
title("fsk经过带通滤波器1后的信号");% 标题
xlabel('时间/s');           % x轴标签
ylabel('幅度');             % y轴标签

%% 第二个带通滤波器设计
fp2=2*(fc2-Rb);             % 第二个带通滤波器频率1
fs2=2*(fc2+Rb);             % 第二个带通滤波器频率2
b2=fir1(30, [fp2/Fs fs2/Fs], boxcar(31));% 生成fir滤波器系统函数中分子多项式的系数
% fir1函数三个参数分别是阶数，数字截止频率，滤波器类型
% 这里是生成了30阶(31个抽头系数)的矩形窗带通滤波器
[h2,w2]=freqz(b2, 1,512);   % 生成fir滤波器的频率响应
% freqz函数的三个参数分别是滤波器系统函数的分子多项式的系数，分母多项式的系数(fir滤波器分母系数为1)和采样点数(默认)512
band2=fftfilt(b2,fsk);      % 对信号进行滤波，fsk是等待滤波的信号，b2是fir滤波器的系统函数的分子多项式系数
subplot(313)                % 窗口分割成3*1的，当前是第3个子图 
plot(t,band2,'LineWidth',2) % 绘制经过低通滤波器后的信号
axis([0,TotalT,-1.1,1.1]);  % 设置坐标范围
title("fsk经过带通滤波器2后的信号");% 标题
xlabel('时间/s');           % x轴标签
ylabel('幅度');             % y轴标签
%% 带通滤波器的频谱
figure(3);                  % 绘制第3幅图
subplot(211);               % 窗口分割成2*1的，当前是第1个子图
plot(w1/pi*Fs/2,20*log(abs(h1)),'LineWidth',2);% 绘制滤波器的幅频响应
xlabel('频率/Hz');          % x轴标签
ylabel('幅度/dB');          % y轴标签
title('第一个带通滤波器频谱');    % 标题

subplot(212);               % 窗口分割成2*1的，当前是第2个子图
plot(w2/pi*Fs/2,20*log(abs(h2)),'LineWidth',2);% 绘制滤波器的幅频响应
xlabel('频率/Hz');          % x轴标签
ylabel('幅度/dB');          % y轴标签
title('第二个带通滤波器频谱');    % 标题

%% 解调部分
%% fsk信号经过第一个带通滤波后和载波1相乘
figure(4);                  % 绘制第4幅图
jt1=band1.*(-zb1);          % 相干解调，乘以相干载波
subplot(211)                % 窗口分割成2*1的，当前是第2个子图 
plot(t,jt1,'LineWidth',1)   % 绘制乘以相干载波后的信号
axis([0,TotalT,-1.5,1.5]);  % 设置坐标范围
title("乘以相干载波后的信号")% 标题
xlabel('时间/s');           % x轴标签
ylabel('幅度');             % y轴标签

%% fsk信号经过第一个带通滤波后和载波1相乘
jt2=band2.*(-zb2);          % 相干解调，乘以相干载波
subplot(212)                % 窗口分割成2*1的，当前是第2个子图 
plot(t,jt2,'LineWidth',1)   % 绘制乘以相干载波后的信号
axis([0,TotalT,-1.5,1.5]);  % 设置坐标范围
title("乘以相干载波后的信号")% 标题
xlabel('时间/s');           % x轴标签
ylabel('幅度');             % y轴标签

%% 加噪信号经过滤波器
%% 低通滤波器设计
fp3=2*Rb;                   % 低通滤波器截止频率
b3=fir1(30, fp3/Fs, boxcar(31));% 生成fir滤波器系统函数中分子多项式的系数
% fir1函数三个参数分别是阶数，数字截止频率，滤波器类型
% 这里是生成了30阶(31个抽头系数)的矩形窗滤波器
[h3,w3]=freqz(b3, 1,512);   % 生成fir滤波器的频率响应
% freqz函数的三个参数分别是滤波器系统函数的分子多项式的系数，分母多项式的系数(fir滤波器分母系数为1)和采样点数(默认)512
lvbo1=fftfilt(b3,jt1);      % 对信号进行滤波，jt1是等待滤波的信号，b3是fir滤波器的系统函数的分子多项式系数
lvbo2=fftfilt(b3,jt2);      % 对信号进行滤波，jt2是等待滤波的信号，b3是fir滤波器的系统函数的分子多项式系数
figure(5);                  % 绘制第5幅图
subplot(311)                % 窗口分割成3*1的，当前是第1个子图 
plot(t,lvbo1,'LineWidth',2) % 绘制经过低通滤波器后的信号
axis([0,TotalT,-0.1,1.1]);  % 设置坐标范围
title("经过低通滤波器后的信号");% 标题
xlabel('时间/s');           % x轴标签
ylabel('幅度');             % y轴标签

subplot(312)                % 窗口分割成3*1的，当前是第2个子图 
plot(t,lvbo2,'LineWidth',2) % 绘制经过低通滤波器后的信号
axis([0,TotalT,-0.1,1.1]);  % 设置坐标范围
title("经过低通滤波器后的信号");% 标题
xlabel(<