基于音频信号的AM、SSB、DSB、FM调制MATLAB仿真

一、读取音频

MATLAB自带一系列音频，例如：chirp（鸟叫声）、train（火车鸣笛)……调用方法如下：

load chirp;
% load train;

load上述信号后，可以看到工作区出现两个变量y（音频信号），Fs（音频信号采样频率）。

注：MATLAB自带的音频信号采样率Fs = 8192Hz。

二、采样率的那些坑

---

---

错误代码1：
如果你使用了下面的代码，程序结果一定会让你十分疑惑。

load chirp;		% y, Fs
N = length(y - 1);
fc = 20e3;		% 载波频率
t = (0:N - 1) / Fs;
st = (y + 0.2) .* cos(2 * pi * fc * t);		% AM信号为例

当你对上述代码生成的信号st进行频谱分析时，会发现频率与理论值不符。问题在于，音频信号的采样频率 Fs = 8192Hz，但采用20KHz的载波，不符合奈奎斯特采样定律。

---

---

错误代码2：
如果你使用了下面的代码，程序结果依然会让你十分疑惑。

load chirp;		% y, Fs
N = length(y - 1);
fc = 20e3;		% 载波频率
fs = 4 * fc;	% AM信号采样频率
t = (0:N - 1) / fs;
st = (y + 0.2) .* cos(2 * pi * fc * t);		% AM信号为例

虽然进行了一定修改，但上述代码仍然存在一定问题。问题在于，上述代码仅仅提高了载波 cos(2 * pi * fc * t) 的采样率，音频信号 y 的采样率为 Fs ，但你却将 y 的采样率硬生生看做 fs = 4 * fc。

---

---

正确的打开方式：
最终，我们明白了问题所在，并解决了问题。

fc = 20e3;
fs = 4 * fc;               	% 对载波的采样频率

y = resample(y,fs,Fs);     	% 关键，对音频信号进行再采样
N = length(y);
t = (0:N - 1) / fs;
st = (y + 0.2) .* cos(2 * pi * fc * t);		% AM信号为例

---

---

三、AM调制与解调

调制：

k = 0.2;		% 控制调幅指数
st = (y + k) .* cos(2 * pi * fc * t);		% AM信号为例

解调：

sr = st .* cos(2 * pi * fc * t);		% 相干解调
point = 6000 / (Fsr / 2);      % 6000Hz略大于chirp音频信号的最高频率
bhi = fir2(34,[0 point point 1],[1 1 0 0]);
freqz(bhi,1,Nr);
yo = filter(bhi,1,yo);                  % 低通滤波
yo = yo - mean(yo);                     % 滤除直流分量

四、DSB调制与解调

调制：

st = y .* cos(2 * pi * fc * t);		% AM信号为例

解调：
同上

五、SSB调制与解调

调制：

y1 = imag(hilbert(y));	% MATLAB函数问题，虚部才是希尔伯特变换

st = 0.5 * y .* cos(2 * pi * fc * t) + 0.5 * y1 .* sin(2 * pi * fc * t);	% +,-均可，对应不同边频

证明如下：

解调：
同上

六、FM调制与解调

看了网上很多文章，大多是对单频的信号进行FM调制解调仿真的。采用单频信号仿真具有一定的局限性，其仿真结果可能会使得我们对FM信号的特点理解产生误差。下面，我们采用一段复杂频率的音频信号，实现FM调制、解调仿真，并阐述其中一系列注意事项。

调制：
FM调制的原理、公式，本文不再赘述，大家可以轻易地在网上查找到FM调制的公式及原理。在此，我们重点介绍，面对一段实际的音频信号，如何合理地确定FM调制的“调频灵敏度Kf”。

我们基于以下两个公式进行讨论：

B = 2 * (dfm  /  Fmax + 1) * Fmax = 2 * (dfm + Fmax)
dfm = kf * max(abs(mt)) / (2 * pi)

其中，
B为调频信号的近似带宽（事实上，调频信号带宽无限大），
mt为基带信号，
Fmax为基带信号mt的上限频率，
dfm为使用mt进行调频后，最大的调频频率。

我们的目的是求解 Kf 。通常，我们可以通过频谱分析，获得基带信号的上限频率 Fmax ，因此，Fmax 是已知的。而调频信号带宽B，应满足如下公式：

B / 2 < fc			---(1)
B / 2 + fc < fs / 2	---(2)

公式(1)保证调频信号正、负频谱不会重叠
公式(2)保证调频信号满足采样率要求

由此，我们可以求解出 dfm ，进而得到 Kf ，公式如下：

B = min(2 * fc, fs - 2 * fc);
B = B * 0.9;				% 留一定的带宽余量，系数<=1
dfm = B / 2 - Fmax;			% Fmax通过对mt频谱分析得到

kf = 2 * pi * dfm / max(abs(mt));

现在，我们就求解出了合理的 Kf 。下面展示FM调制MATLAB仿真代码：

Fmax = 4800;           	% chirp信号的最高频率大概是4800Hz,Fmax对应音频信号最高频率
Bfm = min(2 * fc,Fst - 2 * fc);
Bfm = Bfm * 0.9;       	% 留下一定的带宽余量
dfm = Bfm / 2 - Fmax;  	% dfm表示由音频信号调频产生的最大频偏
kf = 2 * pi * dfm / max(abs(y));
phi = cumsum(y) / Fst;	% 积分
yi = cos(2 * pi * fc * tt + kf * phi);

解调：
本文只列出非相干解调（包络检波）的方法（相关理论推导以及相干解调的方法可以参考其它网络博客），代码如下：

% FM信号解调,利用希尔伯特变换实现包络检波（非相干解调）
ydt = diff(out);
yo = [0 abs(hilbert(ydt))];
yo = yo - mean(yo);			% 滤除直流分量

注意事项：
1、FM 信号的仿真， Kf 的选择尤其重要，直接影响了仿真的结果。
2、FM 信号包络不一定是恒定的，因为FM信号的带宽无限大，很可能出现各个分量叠加、抵消的情况。（网上大多用单频信号进行 FM 调制，他们的实验结果给我造成很大的困惑。）

七、仿真结果

将原始的音频信号、解调后的音频信号，利用sound函数对比，基本没有质量损失，下面展示时域信号及其频谱。

AM:

DSB：

SSB：

FM:

八、结语

仿真文件代码：【radioStation.m】
文章观点或有偏颇之处，如有问题，请多多交流。