模拟角度调制系统仿真.doc

附件1：课程设计报告书格式

理学院

通 信 原 理 课 程 设 计 报 告

设 计 题 目 模拟角度调制系统仿真

学生专业班级 通信09级
学生姓名（学号）

指 导 教 师

完 成 时 间 2012-11-5

实 习（设计）地点 理信学院机房

2012 年 11 月 5 日

一、设计目的和任务分析
1、设计目的
（1）、掌握模拟系统FM、PM调制与解调的原理。
（2）、掌握模拟角度调制系统仿真，实现FM与PM调制与解调的设计方法，并要求信道为AWGN信道，画出调制信号、已调信号的波形图与频谱图，并比较两种调制的带宽、信噪比参数；
（3）、掌握matlab分析系统时域、频域特性的方法，进一步锻炼应用matlab进行编程仿真的能力。
（4）、调制与解调模拟系统，仿真实现相关功能。实现相关功能设计fm、pm调制与解调模拟系统，仿真实现相关功能。
2、任务分析
（1）、确定仿真算法
（2）、建立仿真模型
（3）、设计仿真程序，fm、pm调制与解调，先初始化，添加高斯白噪声制出调制和已调信号，运行仿真程序
（4）、输出仿真结果并进行分析
二、设计过程
1、模拟通信系统简介
通信系统是为了有效可靠的传输信息，信息由信源发出，以语言、图像、数据为媒体,通过电(光)信号将信息传输，由信宿接收。通信系统又可分为数字通信与模拟通信。基于课程设计的要求，下面简要介绍模拟通信系统。
信源是模拟信号，信道中传输的也是模拟信号的系统为模拟通信。模拟通信系统的模型如图1所示。

图1 模拟通信系统模型
调制器: 使信号与信道相匹配, 便于频分复用等。发滤波器: 滤除调制器输出的无用信号。收滤波器: 滤除信号频带以外的噪声，一般设N(t)为高斯白噪声,则Ni(t)为窄带白噪声。
FM波的幅度恒定不变，这使得它对非线性器件不甚敏感，给FM带来了抗快衰落能力。利用自动增益控制和带通限幅还可以消除快衰落造成的幅度变化效应。这些特点使得NBFM对微波中继系统颇具吸引力。WBFM的抗干扰能力强，可以实现带宽与信噪比的互换，因而WBFM广泛应用于长距离高质量的通信系统中，如空间和卫星通信、调频立体声广播、短波电台等。WBFM的缺点是频带利用率低，存在门限效应，因此在接收信号弱、干扰大的情况下宜采用NBFM，这就是小型通信机常采用NBFM的原因。
2、FM调制原理
频率调制的一般表达式为：
（2-1）
FM和PM非常相似，如果预先不知道调制信号的具体形式，则无法判断已调信号是调频信号还是调相信号。

图 2

图 3
图（2）所示的产生调频信号的方法称为直接调频法，图（3）所示的产生调频信号的方法称为间接调频法[4]。由于实际相位调制器的调节范围不可能超出，因而间接调频的方法仅适用于相位偏移和频率偏移不大的窄带调制情形，而直接调频则适用于宽带调制情形。
根据调制后载波瞬时相位偏移的大小，可将频率调制分为宽带调频（WBFM）与窄带调频（NBFM）。宽带与窄带调制的区分并无严格的界限，但通常认为由调频所引起的最大瞬时相位偏移远小于30°时，
  （2-2）
称为窄带调频。否则，称为宽带调频。
 为方便起见，无妨假设正弦载波的振幅A＝1，则由式（2-1）调频信号的一般表达式，得

= （2-3）
通过化解，利用傅立叶变化公式可得NBFM信号的频域表达式：
（2-4）   
在NBFM中，由于下边频为负，因而合成矢量不与载波同相，而是存在相位偏移，当最大相位偏移满足式（2-2）时，合成矢量的幅度基本不变，这样就形成了FM信号。

图4 NBFM信号频谱
3、PM调制原理
在模拟调制中，一个连续波有三个参数可以用来携带信息而构成已调信号。当幅度和频率保持不变时，改变载波的相位使之随未调信号的大小而改变，这就是调相的概念。
角度调制信号的一般表示形式为：
S (t)=Acos[ωt+φ(t)]
式中，A是载波的恒定振幅；[ωt+φ(t)]是信号的瞬时相位，而φ(t)称为瞬时相位偏移；d[ωt+φ(t)]/dt为信号的瞬时频率，而dφ(t)/dt称为瞬时频率偏移，即相对于ω的瞬时频率偏移。
设高频载波为u=Ucosωt，调制信号为UΩ(t)，则调相信号的瞬时相位
φ(t)=ω+KUΩ(t)
瞬时角频率 ω(t)==ω+K
调相信号 u=Ucos［ωt+KuΩ(t)］
将信号的信息加在载波的相位上则形成调相信号，调相的表达式为：
S(t)=Acos[ωt+Kf(t)+φ]
这里K称为相移指数，这种调制方式，载波的幅度和角频率不变，而瞬时相位偏移是调制信号f(t)的线性函数，称为相位调制。
调相与调频有着相当密切的关系，我们知道相位与频率有如下关系式：
ω==ω+Kf(t)
φ(t)=ωt+K
所以在调相时可以先将调制信号进行微分后在进行频率调制，这样等效于调相，此方法称为间接调相，与此相对应，上述方法称为直接调相。调相信号的产生如图5所示：

图5 PM调相信号的产生
实现相位调制的基本原理是使角频率为ω的高频载波u(t)通过一个可控相移网络, 此网络产生的相移Δφ受调制电压uΩ(t)控制, 满足Δφ=KuΩ(t)的关系, 所以网络输出就是调相信号，可控相移网络调相原理图如图6所示：

图6 可控相移网络调相原理图
三、源程序代码
fm：
dt=0.001; %设定时间步长
t=0:dt:1.5; %产生时间向量
am=5; %设定调制信号幅度
fm=5; %设定调制信号频率
mt=am*cos(2*pi*fm*t); %生成调制信号
fc=50; %设定载波频率
ct=cos(2*pi*fc*t); %生成载波
kf=10; %设定调频指数
int_mt(1)=0;
for i=1:length(t)-1
int_mt(i+1)=int_mt(i)+mt(i)*dt; %求信号m(t)的积分
end %调制，产生已调信号
sfm=am*cos(2*pi*fc*t+2*pi*kf*int_mt); %调制信号
%*****************************************
%*************添加高斯白噪声**************
sn1=10; %设定信躁比(小信噪比)
sn2=30; %设定信躁比(大信噪比)
sn=0; %设定信躁比(无信噪比)
db=am^2/(2*(10^(sn/10))); %计算对应的高斯白躁声的方差
n=sqrt(db)*randn(size(t)); %生成高斯白躁声
nsfm=n+sfm; %生成含高斯白躁声的已调信号（信号通
%过信道传输）
%*****************************************
%****************FM解调*******************
for i=1:length(t)-1 %接受信号通过微分器处理
diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
end
diff_nsfmn = abs(hilbert(diff_nsfm)); %hilbert变换，求绝对值得到瞬时幅度（包络检波）
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
diff_nsfmn1=diff_nsfmn-zero;
%*****************************************
%**************时域到频域转换**************
ts=0.001; %抽样间隔
fs=1/ts; %抽样频率
df=0.25; %所需的频率分辨率，用在求傅里叶变换
%时，它表示FFT的最小频率间隔
%*****对调制信号m(t)求傅里叶变换*****
m=am*cos(2*pi*fm*t); %原调信号
fs=1/ts;
if nargin==2
n1=0;
else
n1=fs/df;
end
n2=length(m);
n=2^(max(nextpow2(n1),nextpow2(n2)));
M=fft(m,n);
m=[m,zeros(1,n-n2)];
df1=fs/n; %以上程序是对调制后的信号u求傅里变换
M=M/fs; %缩放，便于在频铺图上整体观察
f=[0:df1:df1*(length(m)-1)]-fs/2; %时间向量对应的频率向量

%************对已调信号u求傅里变换**********
fs=1/ts;
nargin=6；
if nargin==2
n1=0;
else
n1=fs/df;
end
n2=length(sfm);
n=2^(max(nextpow2(n1),nextpow2(n2)));
U=fft(sfm,n);
u=[sfm,zeros(1,n-n2)];
df1=fs/n; %以上是对已调信号u求傅里变换
U=U/fs; %缩放
%******************************************
%***************显示程序******************
disp('按任意键可以看到原调制信号、载波信号和已调信号的曲线')
pause
%**************figure(1)******************
figure(1)
subplot(3,1,1);plot(t,mt); %绘制调制信号的时域图
xlabel('时间t');
title('调制信号的时域图');
subplot(3,1,2);plot(t,ct); %绘制载波的时域图
xlabel('时间t');
title('载波的时域图');
subplot(3,1,3);
plot(t,sfm); %绘制已调信号的时域图
xlabel('时间t');
title('已调信号的时域图');
%******************************************
disp('按任意键可以看到原调制信号和已调信号在频域内的图形')
pause
%************figure(2)*********************
figure(2)
subplot(2,1,1)
plot(f,abs(fftshift(M))) %fftshift:将FFT中的DC分量移到频谱中心
xlabel('频率f')
title('原调制信号的频谱图')
subplot(2,1,2)
plot(f,abs(fftshift(U)))
xlabel('频率f')
title('已调信号的频谱图')
%******************************************
disp('按任意键可以看到原调制信号、无噪声条件下已调信号和解调信号的曲线')
pause
%**************figure(3)******************
figure(3)
subplot(3,1,1);plot(t,mt); %绘制调制信号的时域图
xlabel('时间t');
title('调制信号的时域图');
subplot(3,1,2);plot(t,sfm); %绘制已调信号的时域图
xlabel('时间t');
title('无噪声条件下已调信号的时域图');
nsfm=sfm;
for i=1:length(t)-1 %接受信号通过微分器处理
diff_nsfm(i)=(nsfm(i+1)-nsfm(i))./dt;
end
diff_nsfmn = abs(hilbert(diff_nsfm)); %hilbert变换，求绝对值得到瞬时幅度（包络检波）
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
diff_nsfmn1=diff_nsfmn-zero;
subplot(3,1,3); %绘制无噪声条件下解调信号的时域图
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
xlabel('时间t');
title('无噪声条件下解调信号的时域图');
%*****************************************
disp('按任意键可以看到原调制信号、小信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')
pause
%**************figure(4)******************
figure(4)
subplot(3,1,1);plot(t,mt); %绘制调制信号的时域图
xlabel('时间t');
title('调制信号的时域图');
db1=am^2/(2*(10^(sn1/10))); %计算对应的小信噪比高斯白躁声的方差
n1=sqrt(db1)*randn(size(t)); %生成高斯白躁声
nsfm1=n1+sfm; %生成含高斯白躁声的已调信号（信号通
%过信道传输）
for i=1:length(t)-1 %接受信号通过微分器处理
diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;
end
diff_nsfmn1 = abs(hilbert(diff_nsfm1)); %hilbert变换，求绝对值得到瞬时幅度（包络检波）
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
diff_nsfmn1=diff_nsfmn1-zero;
subplot(3,1,2);
plot(1:length(diff_nsfm),diff_nsfm); %绘制含小信噪比高斯白噪声已调信号的时域图
xlabel('时间t');
title('含小信噪比高斯白噪声已调信号的时域图');
subplot(3,1,3); %绘制含小信噪比高斯白噪声解调信号的时域图
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
xlabel('时间t');
title('含小信噪比高斯白噪声解调信号的时域图');
%*****************************************
disp('按任意键可以看到原调制信号、大信噪比高斯白噪声条件下已调信号和解调信号已调信号的曲线')
pause
%**************figure(5)******************
figure(5)
subplot(3,1,1);plot(t,mt); %绘制调制信号的时域图
xlabel('时间t');
title('调制信号的时域图');
db1=am^2/(2*(10^(sn2/10))); %计算对应的大信噪比高斯白躁声的方差
n1=sqrt(db1)*randn(size(t)); %生成高斯白躁声
nsfm1=n1+sfm; %生成含高斯白躁声的已调信号（信号通过信道传输）
for i=1:length(t)-1 %接受信号通过微分器处理
diff_nsfm1(i)=(nsfm1(i+1)-nsfm1(i))./dt;
end
diff_nsfmn1 = abs(hilbert(diff_nsfm1)); %hilbert变换，求绝对值得到瞬时幅度（包
%络检波）
zero=(max(diff_nsfmn)-min(diff_nsfmn))/2;
diff_nsfmn1=diff_nsfmn1-zero;
subplot(3,1,2);
plot(1:length(diff_nsfm1),diff_nsfm1); %绘制含大信噪比高斯白噪声已调信号
%的时域图
xlabel('时间t');
title('含大信噪比高斯白噪声已调信号的时域图');
subplot(3,1,3); %绘制含大信噪比高斯白噪声解调信号
%的时域图
plot((1:length(diff_nsfmn1))./1000,diff_nsfmn1./400,'r');
xlabel('时间t');
title('含大信噪比高斯白噪声解调信号的时域图');
%*****************************************
%******************结 束*******************
pm：
function [v,phi]=env_phas(x,ts,f0)
if nargout==2 %nargout为输出变数的个数
z=loweq(x,ts,f0); %产生调制信号的正交分量
phi=angle(z); %angle是对一个复数求相角的函数
end
v=abs(hilbert(x)); %abs用来求复数hilbert(x)的模

function [M,m,df]=fftseq(m,ts,df)
fs=1/ts;
if nargin==2 n1=0; %nargin为输入参量的个数
else n1=fs/df;
end
n2=length(m);
n=2^(max(nextpow2(n1),nextpow2(n2))); %nextpow2(n)取n最接近的较大2次幂
M=fft(m,n); %M为信号m的傅里叶变换，n为快速傅里叶变换的点数，及基n-FFT变换
m=[m,zeros(1,n-n2)]; %构建新的m信号
df=fs/n;
%重新定义频率分辨率
function x1=loweq(x,ts,f0)
t=[0:ts:ts*(length(x)-1)];
z=hilbert(x); %希尔伯特变换对的利用---通过实部来求虚部
x1=z.*exp(-j*2*pi*f0*t); %产生信号z的正交分量，

t0=0.2; %信号的持续时间,用来定义时间向量
ts=0.001; %抽样间隔
fs=1/ts; %抽样频率
fc=300; %载波频率,fc可以任意改变
t=[-t0/2:ts:t0/2]; %时间向量
kf=100; %偏差常数
df=0.25; %所需的频率分辨率，用在求傅里叶变换时，它表示FFT的最小频率间隔
m=sin(100*t); %调制信号，m(t)可以任意更改
int_m(1)=0; %求信号m(t)的积分
for i=1:length(t)-1
int_m(i+1)=int_m(i)+m(i)*ts;
end
[M,m,df1]=fftseq(m,ts,df); %对调制信号m(t)求傅里叶变换
M=M/fs; %缩放，便于在频谱图上整体观察
f=[0:df1:df1*(length(m)-1)]-fs/2; %时间向量对应的频率向量
u=cos(2*pi*fc*t+2*pi*kf*int_m); %调制后的信号
[U,u,df1]=fftseq(u,ts,df); %对调制后的信号u求傅里叶变换
U=U/fs; %缩放
%通过调用子程序env_phas和loweq来实现解调功能
[v,phase]=env_phas(u,ts,fc); %解调，求出u的相位
phi=unwrap(phase); %校正相位角，使相位在整体上连续，便于后面对该相位角求导
dem=(1/(2*pi*kf))*(diff(phi)*fs); %对校正后的相位求导
%再经一些线性变换来恢复原调制信号
%乘以fs是为了恢复原信号，因为前面使用了缩放
subplot(3,2,1) %子图形式显示结果
plot(t,m(1:length(t))) %现在的m信号是重新构建的信号，
%因为在对m求傅里叶变换时m=[m,zeros(1,n-n2)]
axis([-0.1 0.1 -1 1]) %定义两轴的刻度
xlabel('时间t')
title('原调制信号的时域图')
subplot(3,2,2)
plot(t,u(1:length(t)))
axis([-0.1 0.1 -1 1])
xlabel('时间t')
title('已调信号的时域图')
subplot(3,2,3)
plot(f,abs(fftshift(M))) %fftshift：将FFT中的DC分量移到频谱中心
axis([-600 600 0 0.04])
xlabel('频率f')
title('原调制信号的频谱图')
subplot(3,2,4)
plot(f,abs(fftshift(U)))
axis([-600 600 0 0.04])
xlabel('频率f')
title('已调信号的频谱图')
subplot(3,2,5)
plot(t,m(1:length(t)))
axis([-0.1 0.1 -1 1])
xlabel('时间t')
title('原调制信号的时域图')
subplot(3,2,6)
plot(t,dem(1:length(t)))
axis([-0.1 0.1 -1 1])
xlabel('时间t')
title('解调后信号的时域波形')
四、运行结果
FM

PM

五、结论与体会

通过这一次课程设计，增长了通信原理的基本知识，使自己充分体会到了在设计过程中的成功喜悦，也让我明白在以后的工作和生活中都应该不断学习，努力提高自己的知识和综合素质。
六、参看书目
[1] 樊昌信. 通信原理（第6版）.国防工业出版社，2006，09
[2] 任嘉伟. 数字频带通信系统计算机仿真[J].电脑知识与技术，2008，07
[4] 吕跃广 通信系统仿真.电子工业出版社,2010.03

课程设计成绩评定表
学生姓名

专业班级
通信工程09级
设计题目
模拟角度调制系统仿真
指导教师评语及意见：

指导教师评阅成绩： 指导教师签字：

年 月 日