模拟调制与解调
前言
主要内容来自参考资料[2],在学习记录的同时勘正了书中代码的部分错误。
一、拟合
常用函数:polyfit(),polyval(),功能如下图所示:
二、调制与解调仿真
代码如下(示例):
SSB(单边带调制)
希尔伯特变换
m ^ ( t ) = m ( t ) / π t M ^ ( f ) = − j × s g n ( f ) M ( f ) \hat{m}(t)=m(t)/\pi t\,\,\,\,\hat{M}(f)=-j\times sgn(f)M(f) m^(t)=m(t)/πtM^(f)=−j×sgn(f)M(f)其中 m ^ ( t ) \hat{m}(t) m^(t)为 m ( t ) m(t) m(t)的希尔伯特变换,而 M ^ ( f ) \hat{M}(f) M^(f), M ( f ) M(f) M(f)为两者的傅里叶变换
假设
m
(
t
)
m(t)
m(t)为消息信号,则对其进行单边带调制,得到调制信号
u
(
t
)
u(t)
u(t)为
u
(
t
)
=
A
c
m
(
t
)
c
o
s
(
2
π
f
c
t
)
/
2
∓
A
c
m
^
(
t
)
c
o
s
(
2
π
f
c
t
)
/
2
u(t)=A_cm(t)cos(2\pi f_ct)/2\mp A_c\hat{m}(t)cos(2\pi f_ct)/2
u(t)=Acm(t)cos(2πfct)/2∓Acm^(t)cos(2πfct)/2使用MATLAB进行仿真。对正弦信号 y = sin(t) 进行希尔伯特变换,并绘图。
clear all;close all
t = 0:0.1:30;
y = sin(t);
s_y = hilbert(y);
plot(t,y,'linewidth',2);hold on
plot(t,real(s_y),'ro',t,imag(s_y),'bs');
legend('原信号','希尔伯特变换器实部(仍然是原信号)','希尔伯特变换器虚部')
绘图结果如下
SSB调制
基带信号是在150Hz~400Hz内,幅度随频率逐渐递减的音频信号,对其进行SSB调制,仿真时间为1s。载波频率为1000Hz,采样频率为10000Hz。
示例MATLAB代码:
% 基带信号是在150Hz~400Hz内,幅度随频率逐渐递减的音频信号,对其进行SSB调制,仿真时间为1s
clear all;close all;
fs = 10000;
fc = 1000;
t = 1/fs:1/fs:1;%仿真时间
f_v = 0:length(t)-1;
m(fs)=0;%创建长度为fs的消息信号
for f = 150:400
m = m + 0.01*sin(2*pi*f*t)*(400-f);%m(t)包含许多频率分量,且频率分量的幅度随频率递减
end
m_shift = imag(hilbert(m));%基带信号的希尔伯特变换
carriercos = cos(2*pi*fc*t);%载波
carriersin = sin(2*pi*fc*t);%载波的正交变换
ssb_u = m.*carriercos - m_shift.*carriersin; %上边带SSB
ssb_l = m.*carriercos + m_shift.*carriersin; %上边带SSB
figure;
subplot(421);plot(t,carriercos,t,carriersin,':m');axis([0 0.01 -1 1]);title('载波及其正交变换')%载波
subplot(422);plot(f_v,abs(fft(carriercos)));axis([0 2000 -500 12000]);title('载波频谱图') %载波频谱
subplot(423);plot(t,m);axis([0 0.01 -500 500]);title('基带信号波形图') %基带信号
subplot(424);plot(f_v,abs(fft(m)));axis([0 2000 -500 12000]);title('基带信号频谱图') %基带信号频谱图
subplot(425);plot(t,ssb_u);axis([0 0.01 -500 500]);title('%SSB上边带波形图') %SSB上边带波形图
subplot(426);plot(f_v,abs(fft(ssb_u)));axis([0 2000 -500 12000]);title('SSB上边带频谱图') %SSB上边带频谱图
subplot(427);plot(t,ssb_l);axis([0 0.01 -500 500]);title('SSB下边带波形图') %SSB下边带波形图
subplot(428);plot(f_v,abs(fft(ssb_l)));axis([0 2000 -500 12000]);title('SSB下边带频谱图') %SSB下边带频谱图
suptitle('波形图与频谱图')
仿真结果:
SSB解调
在实际应用中,一般采用高稳定度的晶体振荡器或频率合成器来产生本地解调载波,无需像DSB解调那样用PLL恢复载波,大大降低了单边带接收机的技术复杂度与成本。
在本例中,将采用相干解调法,即产生本地振荡波形与接收信号相乘,再经过低通滤波,恢复原始信号,但由于本地振荡与载波不严格同频同相,此法存在一定的频率偏移和相位偏移。
示例代码:
out = ssb_u.*carriercos;%相干解调
[a,b] = butter(4,500/(fs/2)); %设计4阶低通滤波器,且截止频率为500Hz
out1 = filter(a,b,out); %滤波输出
figure;
subplot(323);plot(t,out);axis([0 0.01 -500 500]);title('解调信号波形图')
subplot(324);plot(f_v,abs(fft(out)));axis([0 2000 -500 12000]);title('解调信号频谱图')
subplot(325);plot(t,out1);axis([0 0.01 -500 500]);title('滤波信号波形图')
subplot(326);plot(f_v,abs(fft(out1)));axis([0 2000 -500 12000]);title('滤波信号频谱图')
subplot(321);plot(t,ssb_u);axis([0 0.01 -500 500]);title('%SSB上边带波形图') %SSB上边带波形图
subplot(322);plot(f_v,abs(fft(ssb_u)));axis([0 2000 -500 12000]);title('SSB上边带频谱图') %SSB上边带频谱图
suptitle('波形图与频谱图')
仿真结果图:
角度调制
频率调制
示例代码
%% 频率调制 FM
clear all;close all;
t0 = 0.15;%信号保持时间
ts = 0.0005;%采样时间
fc = 200;%载波频率
kf = 50;%调制系数
fs = 1/ts;%采样频率
df = 0.25;%频率分辨率
t_v = [0:ts:t0];%时间序列
m = [ones(1,t0/(3*ts)),-2*ones(1,t0/(3*ts)),zeros(1,t0/(3*ts)+1)];%消息信号m(t),t/ts是序列的数目。注意端点,要+1
m_int(length(m)) = 0;
%以下两种方法都可实现对m(t)的积分。法二更直观,一定要乘以ts!
for i = 1:length(m)-1
m_int(i+1) = m_int(i) + m(i)*ts;
end
% for i = 1:length(m)
% m_int(i) = ts*sum(m(1:i));
% end
[M,m,df1] = fftseq(m,ts,df);%傅里叶变换
M = M/fs;
f = [0:df1:df1*(length(m)-1)]-fs/2;%画[-fs/2,fs/2]的图像
u = cos(2*pi*fc*t_v + 2*pi*kf*m_int);%调制信号
figure
plot(t_v,u)
[U,u,df1] = fftseq(u,ts,df);%傅里叶变换
U = U/fs;
%画图
figure;
subplot(221);plot(t_v,m(1:length(t_v)));%经过fftseq函数后,m的长度与t_v的长度不一定相等,但有效长度相等
axis([0 0.15 -2.1 2.1]);
xlabel('时间');title('消息信号')
subplot(222);plot(t_v,u(1:length(t_v)))
axis([0 0.15 -2.1 2.1]);
xlabel('时间');title('调频信号')
subplot(223);plot(f,abs(fftshift(M)));
xlabel('频率');title('消息信号频谱')
subplot(224);plot(f,abs(fftshift(U)));
xlabel('频率');title('调制信号频谱')
仿真结果图:
相位调制
示例代码
clear all;close all;
t0 = 0.25;%信号保持时间
ts = 0.0005;%采样时间
fc = 200;%载波频率
fs = 1/ts;%采样频率
df = 0.25;%频率分辨率
t_v = [0:ts:t0];%时间序列
m(length(t_v)) = 0;
lent = length(t_v)-1;
for i = 1:(t0/ts/4)
m(i) = i;
end
for i = (lent/4+1):3*(lent/4)
m(i) = m(lent/4) - i+lent/4;
end
for i = (3*lent/4+1):lent+1
m(i) = m(3*lent/4) + i-lent*3/4
end
m = m/50;
[M,m,df1] = fftseq(m,ts,df);%傅里叶变换
M = M/fs;
f = [0:df1:df1*(length(m)-1)]-fs/2;%画[-fs/2,fs/2]的图像
u = cos(2*pi*fc*t_v + m(1:length(t_v)));%调制信号
[U,u,df1] = fftseq(u,ts,df);%傅里叶变换
U = U/fs;
%画图
figure;
subplot(221);plot(t_v,m(1:length(t_v)));%经过fftseq函数后,m的长度与t_v的长度不一定相等,但有效长度相等
axis([0 0.25 -3 3]);
xlabel('时间');title('消息信号')
subplot(222);plot(t_v,u(1:length(t_v)))
axis([0 0.15 -2.1 2.1]);
xlabel('时间');title('调频信号')
subplot(223);plot(f,abs(fftshift(M)));
xlabel('频率');title('消息信号频谱')
subplot(224);plot(f,abs(fftshift(U)));
xlabel('频率');title('调制信号频谱')
仿真结果图:
相关函数
fft()和fftshift()的区别
M = fft(m)是对消息信号 m(t) 进行快速FFT变换,得到取值范围为 (0,fs) 的函数 M(f) ,如果想画出取值范围为 (-fs/2,fs/2) 的函数时,需要用到 fftshift() 函数。
需要注意的是,使用时直接 fftshift(M) 即可,该函数只是实现对函数 M(f) 的移位,并未对其进行FFT变换。
详情见博客: https://www.cnblogs.com/limanjihe/p/10014142.html
参考资料
[1] ttps://blog.csdn.net/u014630987/article/details/70156489
[2] 《MATLAB通信系统建模与仿真(第2版)》 邓奋发 编著