软件无线电1-MATLAB实现FM调制解调
1、MATLAB读取语音文件
准备一段wav的语音文件，我用笔记本自带的灌音机录制了一段自己的语音“爱福皮的姐.wav“，MATLAB读取语音文件，并获取采样率信息。
clc;
clear all;
%% ***************read file*************************************
filename ='./data_source/爱福皮的姐.wav';
= audioread(filename);% fm = 48e3; %音频的采样率
sound(m,fm);
figure(1);
plot(m);title('原始的语音信号');
这段原始语音采样率48kHz，时域波形如图1所示，
https://i-blog.csdnimg.cn/direct/04eea6452cc94a18b57ff759630f0aa6.png
2、FM调制
根据FM调制数学表达式:
https://i-blog.csdnimg.cn/direct/0b3588782dc64a528e13fd048f4b3e04.png生成FM信号。
2.1相干参数设置和语音数据源初处理
我是用的载波频率192kHz，载波采样率1.92MHz，代码里雷同这些参数大家可以根据自己的需求修改。FM宽带调制最大频偏75kHz，FM窄带调制最大频偏5kHz。由于FM调制是用基带信号控制载波的频率，即与载波信号数据相乘，所以在这之前，我们须要将基带信号的采样率和长度举行处理，使其与载波信号采样率和长度保持同等，处理代码如下：
fs = 1.92e6; %采样率 载波的采样率
fc = 192e3; %载波中心频率
df = 75e3; %最大频偏宽带75khz 窄带时5khz
kf = (df*2^32/fs)/32767;
AC= 1024;%幅度
dt=1/fs;
m_len=length(m)*fs/fm;%按照载波采样率的长度计算
m_t=zeros(1,m_len);
for i=1:length(m)
for j=1:fix(fs/fm)
m_t((i-1)*fix(fs/fm)+j)=m(i);
end
end
2.2基带信号FM调制生成基带IQ信号
查表法实现信号发生器，MATLAB产生DDS，代码如下：
%************************* FM调制 ********************
n=0:1/1024:1023/1024;
s_rom=sin(2*pi*n);
c_rom=cos(2*pi*n);
w_r=0;%相位累加器
rrom_addr=0;%
dac_i=zeros(1,m_len);%I路信号
dac_q=zeros(1,m_len);%Q路信号
for i=1:m_len
w_r = w_r + kf*m_t(i);
if(w_r > 2^32) % 做32位累加器的溢出判断
w_r= w_r - 2^32;
elseif(w_r <0)
w_r = w_r + 2^32; % 负的溢出时
end
rrom_addr=round(w_r/2^22);%读查找表的地址
if rrom_addr == 0
rrom_addr =1;
end
dac_q(i)=AC*s_rom(rrom_addr);
dac_i(i)=AC*c_rom(rrom_addr);
end
2.3 生成载波IQ信号
查表法实现信号发生器，MATLAB产生DDS，代码如下：
%*****************************载波信号 *********************************
w=fc*2^32/fs; %频率控制字，32位的相位累加器， f_out= fs*w/2^N :
n=0:1/1024:1023/1024;
s_rom=sin(2*pi*n);
c_rom=cos(2*pi*n);
%接下来产生载波地址，
w_r=0;%相位累加器
rrom_addr=0;%
cw_sin=zeros(1,m_len);%这个相当于是我们的查找表，一个周期的正弦波
cw_cos=zeros(1,m_len);
for i=1:m_len
w_r = w_r + w;
if(w_r > 2^32) % 做32位累加器的溢出判断
w_r= w_r - 2^32;
end
rrom_addr=round(w_r/2^22);%读查找表的地址
if rrom_addr == 0%matlab计数时从1开始计数
rrom_addr =1;
end
cw_sin(i)=s_rom(rrom_addr);%载波完成
cw_cos(i)=c_rom(rrom_addr);%载波
end
figure(2);
subplot(2,1,1);plot(cw_sin(1:3000),'r');title('sin载波信号');
subplot(2,1,2);plot(cw_cos(1:3000),'b');title('cos载波信号');
载波信号时域波形如图2所示。
https://i-blog.csdnimg.cn/direct/e9fb47e61ef845879d8354b0914efc83.png
图2 载波信号时域波形
2.4 生成FM信号
基带IQ信号分别和载波IQ信号相乘，之后相加，得到FM信号。
%******************* FM信号 *******************
s_t=zeros(1,m_len);
for i=1:m_len
s_t(i) = dac_i(i)*cw_cos(i) + dac_q(i)*cw_sin(i)*(-1);
end
figure(3);
subplot(2,1,1);plot(s_t);title('FM信号');
subplot(2,1,2);psd(s_t,fs);title('FM信号频谱');
FM窄带调制生成的FM信号时域波形和频谱如图3所示。
https://i-blog.csdnimg.cn/direct/9102cbace35d48a08f4f35201dc19dca.png
图3 FM窄带调制生成的FM信号时域波形和频谱
FM宽带调制生成的FM信号时域波形和频谱如图4所示。
https://i-blog.csdnimg.cn/direct/2eba8f629b1a430092ed6491204cdb9d.png
图4 FM宽带调制生成的FM信号时域波形和频谱
3 FM干系解调
3.1 去除高频分量
干系解调首先将FM信号与载波信号IQ信号相乘，然后颠末低通滤波器，滤除高频分量。
%% 解调方法一 %*******************相干解调 ok*******************
%%%%%%%%%%%%正交解调%%%%%%%%%%%%%%%%
i_data = zeros(1,m_len);
q_data = zeros(1,m_len);
for i = 1 : m_len
i_data(i) = s_t(i) * cw_cos(i);
q_data(i) = (-1)* s_t(i) * cw_sin(i);
end
%%%%%%%%%%%%%%进行低通滤波%%%%%%%%%%%%%%%%
%滤波器系数
f20k= [-4.49340747782e-05,0.00010035338004,0.00018271359187,0.000333071475,……];
adc_i=filter(f20k,1,i_data);
adc_q=filter(f20k,1,q_data);
滤波器系数获取步调：
步调1：使用MATLAB自带的滤波器设计工具，即fdatool获取滤波器系数。在MATLAB命令窗口输入fdatool回车，弹出滤波器设计工具
步调2：滤波器参数设置。
选择低通滤波器，然后设置采样率、通带、阻带等参数后，最后点击设计滤波器，期待结束，生成滤波器，如图5所示：
https://i-blog.csdnimg.cn/direct/d477eda3636c408b870683f43836e8ee.png
图5 滤波器参数设置
步调3：导出滤波器系数。
滤波器设计工具右上角“file“下拉选择”Export“，设置导出系数的名称，覆盖原数据，导出如图6所示，在MATLAB工作区找到刚刚取的名称。
https://i-blog.csdnimg.cn/direct/bb473f69c17049c99ee80999678325ac.png
图6 导出滤波器系数
3.2 FM解调
根据如下FM解调公式，得到我们须要的数据源基带信号。
https://latex.csdn.net/eq?f%5Cleft%28n%5Cright%29%3D%5Cfrac%7BX_I%5Cleft%28n%5Cright%29%5Cast%20X_Q%5Cleft%28n-1%5Cright%29-X_I%5Cleft%28n-1%5Cright%29%5Cast%20X_Q%5Cleft%28n%5Cright%29%7D%7BX_I%5Cleft%28n%5Cright%29%5Cast%20X_I%5Cleft%28n%5Cright%29+X_Q%5Cleft%28n%5Cright%29%5Cast%20X_Q%5Cleft%28n%5Cright%29%7D
%%%%%%%%%%%%%FM解调%%%%%%%%%%%%%%%%%%%%%%
c_len = length(adc_i);
cr = zeros(1,c_len);
cj = zeros(1,c_len);
for i = 2: c_len %必须从2开始，因为有I-1
cr(i) = adc_i(i)*adc_i(i) +adc_q(i)*adc_q(i);
cj(i)= adc_i(i-1) * adc_q(i) -adc_i(i)* adc_q(i-1);%I(n-1)*Q(n) -I(n)*Q(n-1)
end
angle = zeros(1,c_len);
for i = 1:c_len
if cr(i) == 0
angle(i) = 0;
else
angle(i) = (cj(i)/cr(i));
end
end
3.3 降采样率处理
上面得到的基带信号采样率依然是1.92MHz，我们使用MATLAB自带的函数resample函数举行40倍抽取处理，将基带信号的采样率降到48kHz，为了MATLAB播放出来效果更好，将信号幅值乘10处理后播放，我们可以将原语音数据播放效果和解调后的语音数据举行比较，解调后的语音听起来效果还可以，非常清晰。
%%%%%%%%%%%%%降采样%%%%%%%%%%%%%%%%%%%%%%%%%
data_out1 = resample(angle,1,40);
sound(data_out1*10,fm);
figure(4);
plot(data_out1(100:end)*10);title('FM解调信号');
FM解调信号时域波形如图7所示。
https://i-blog.csdnimg.cn/direct/7d423570de014aecbeda0745cae8b4de.png
图7 FM解调信号1
4 FM非干系解调
非干系解调不须要提取载波信息，将接收到的FM信号颠末微分器，然后颠末希尔伯特变更和包络检波，最后去偏置，得到的基带信号采样率1.92Mhz，同样使用MATLAB自带的函数resample举行40倍抽取处理，将基带信号的采样率降到48kHz，为了MATLAB播放出来效果更好，将信号幅值除以30处理后播放，我们可以将原语音数据播放效果和解调后的语音数据举行比较，非干系解调后的语音听起来也非常清晰。
%% 解调方法二%*******************非相干解调 ok*******************
for i=1:m_len-1
diff_sfm(i)=(s_t(i+1)-s_t(i))/dt;%微分器，将调频波变换成调幅调频波
end
diff_sfm=abs(hilbert(diff_sfm)); %希尔伯特变换移相，然后包络检波
diff_sfm=(diff_sfm/AC-2*pi*fc)/kf; %去偏置
data_out2 = resample(diff_sfm,1,40);
sound(data_out2/30,fm);
figure(4);
plot(data_out1(100:end)*10);title('FM解调信号');
FM解调信号时域波形如图8所示。
https://i-blog.csdnimg.cn/direct/1e3b0fd839f741529433be14bc742dcd.png
图8 FM解调信号2
5 MATLAB里FM相干库函数实现FM调制解调
MATLAB软件自带了许多库函数，可以直接调用。FM调制我们可以使用fmmod函数，FM解调使用demod函数，代码如下，同样也举行了降采样率处理后，播放终极处理的效果数据，解调后的语音听起来也非常清晰。
%% 解调方法三
%*******************matlab 库函数调制解调 ok*******************
y=fmmod(m_t,fc,fs,10);
so = demod(y,fc,fs,'fm'); %FM解调
data_out3 = resample(so,1,40);
sound(data_out3,fm);
figure(6);
subplot(2,1,1);plot(y,'b');title('FM信号');
subplot(2,1,2);plot(so,'r');title('FM解调信号');
FM信号和解调信号时域波形如图9所示。
https://i-blog.csdnimg.cn/direct/7d3ab1509b26429685c52a1b5efeb88c.png
图9 FM调制信号和解调信号
免责声明：如果侵犯了您的权益，请联系站长，我们会及时删除侵权内容，谢谢合作！更多信息从访问主页：qidao123.com:ToB企服之家，中国第一个企服评测及商务社交产业平台。
页:
[1]