Matlab通信仿真系列——角度调制之调频FM

微信公众号上线，搜索公众号小灰灰的FPGA,关注可获取相关源码，定期更新有关FPGA的项目以及开源项目源码，包括但不限于各类检测芯片驱动、低速接口驱动、高速接口驱动、数据信号处理、图像处理以及AXI总线等

本节目录

一、调频FM
1、公式推导
2、调频FM的matlab仿真验证
二、FM信号的解调
三、Matlab源码
1、调频FM的matlab仿真验证源码
2、FM信号的解调源码

本节内容
一、调频FM
频率调制FM，使高频载波的频率按调制信号的规律变化而振幅保持恒定的调制方式，可以看作是载波角度的变化。
已调信号频谱不再是原调制信号频谱的线性搬移，而是频谱的非线性变换，会产生与频谱搬移不同的新的频率成分，故又称为非线性调制。
1、公式推导
瞬时频率直接随调制信号线性变化，瞬时角频率和瞬时相位

其中Kf表示调频器的灵敏度，单位是Hz/V
若m(t)=Amcos2πfmt,用其对载波作频率调制，已调信号表达式为

调频波的频谱包含个无穷多个分量。理论上，调频波的频带宽度为无限宽。当取适当的n值使分频变量小到可以忽略不计的程度，调频信号可近似认为具有有限频谱。
当βf≥1，取边频数n=βf+1。
调频波的带宽为Bfm=2(βf+1)W。
信号功率Pfm=A^2/2。
2、调频FM的matlab仿真验证

matlab代码示例：
消息信号是[-3 3]均匀分布的随机整数，
产生的时间间隔为1/10s,
用FM方法调制载波cos2πfct。
Kf=50,fc=250,0≤t≤10,消息信号的带宽W=50Hz；
求：
①消息信号和已调信号的频谱
②已调信号的功率、消息信号的功率、
调制指数及调制信号的带宽

上图是A=1，功率计算Pc=0.5=A^2/2。
下图是将A的值设置为4，功率计算Pc=8。

调制指数βf=kf×max|m(t)|/B=3。
带宽为Bfm=2(βf+1)W=2×(3+1)×50=400。
二、FM信号的解调
调频信号的瞬时频率正比于调制信号的幅度，调频信号的解调器必须能产生正比于输入频率的输出电压。解调器可以采用具有频率-电压转化特性的鉴频器。

理想鉴频器可以看作一个带微分器的包络检波器，微分器输出为

幅度、频率均含调制信息的调幅调频信号，包络检波器将其幅度变化取出，并滤去直流后得到mo(t),其中Kd称为鉴频器灵敏度。

matlab代码示例：
消息信号是[-3 3]均匀分布的随机整数，
产生的时间间隔为1/2s,用FM调制载波cos2πfct。
fc=300,0≤t≤5，kf=50；求：
①用鉴频法解调该信号，绘制原始信号和解调信号
②调制信号通过awgn信道，信噪比为20dB，
绘制解调后的信号与原始信号

三、Matlab源码
1、调频FM的matlab仿真验证源码

clear all;
clc;
%信号抽样时间间隔
ts=0.001
%时间矢量
t=0:ts:10-ts;
%抽样频率
fs=1/ts;
%fft的频率分辨率
df=fs/length(t);
%生成随机数种子的消息序列
msg=randi([-3 3],100,1);
%扩展成取样信号形式
msg1=msg*ones(1,fs/10);
msg2=reshape(msg1.',1,length(t));
%求消息信号的频谱
Pm=fft(msg2)/fs;
f=-fs/2:df:fs/2-df;
subplot(2,1,1);
plot(f,fftshift(abs(Pm)));
title('消息信号频谱');
axis([-200 200 0 5]);
%消息信号积分
int_msg(1)=0;
for ii=1:length(t)-1
    int_msg(ii+1)=int_msg(ii)+msg2(ii)*ts;
end
A=4;
%调频器灵敏度
kf=50;
%载波频率
fc=250;
%已调信号
Sfm=A*cos(2*pi*fc*t+2*pi*kf*int_msg);
%已调信号的频谱
Pfm=fft(Sfm)/fs;
subplot(2,1,2);
plot(f,fftshift(abs(Pfm)));
title('FM信号频谱');
axis([-1000 1000 0 5]);
%Sussb已调信号功率
Pc=sum(abs(Sfm).^2)/length(Sfm);
%消息信号功率
Ps=sum(abs(msg2).^2)/length(msg2);
%消息信号带宽    
B=50;
%调制指数
betaf=kf*max(msg)/B;
%调制信号带宽
W=2*(betaf+1)*B;

2、FM信号的解调源码

clear all;
clc;
%信号抽样时间间隔
ts=0.001
%时间矢量
t=0:ts:5-ts;
%抽样频率
fs=1/ts;
%fft的频率分辨率
df=fs/length(t);
f=-fs/2:df:fs/2-df;
%生成随机数种子的消息序列
msg=randi([-3 3],10,1);
%扩展成取样信号形式
msg1=msg*ones(1,fs/2);
msg2=reshape(msg1.',1,length(t));
subplot(3,1,1);
plot(t,msg2);
title('消息信号');
axis([0 5 -5 5]);
%消息信号积分
int_msg(1)=0;
for ii=1:length(t)-1
    int_msg(ii+1)=int_msg(ii)+msg2(ii)*ts;
end
A=1;
%调频器灵敏度
kf=50;
%载波频率
fc=250;
%已调信号
Sfm=A*cos(2*pi*fc*t+2*pi*kf*int_msg);
%FM调制信号相位
phase=angle(hilbert(Sfm).*exp(-j*2*pi*fc*t));
%相位角
phi=unwrap(phase);
%求相位微分，得到消息信号
dem=(1/(2*pi*kf)*diff(phi)/ts);
dem(length(t))=0;
subplot(3,1,2);
plot(t,dem);
title('无噪声的解调信号FM');
%FM调制信号经过awgn信道
y1=awgn(Sfm,20,'measured');
axis([0 5 -5 5]);
%调制信号限幅
y1(find(y1>1))=1;
y1(find(y1<-1))=-1;
%FM调制信号相位
phase1=angle(hilbert(y1).*exp(-j*2*pi*fc*t));
%相位角
phi1=unwrap(phase1);
%求相位微分，得到消息信号
dem1=(1/(2*pi*kf)*diff(phi1)/ts);
dem1(length(t))=0;
subplot(3,1,3);
plot(t,dem1);
title('信噪比为20dB的解调信号FM');
axis([0 5 -5 5]);