通信系统仿真上机实验指导书
指导教师：刘健
2013年10月
1
最终报告内容：
根据提供的例子，学习通信系统仿真基本知识，搭建包含：随机数生成、调制、编码、AWGN信道、解码、解调、差错统计等模块的通信系统仿真M文件或是Simulink仿真模型，最终得出评估通信系统性能的SNR/BER曲线、频谱利用率等曲线。
2
目录
目录................................................................................................................................. 1 设计一 随机信号分析................................................................................................... 4 设计二 模拟信号的数字化........................................................................................... 8 设计三 数字基带传输系统......................................................................................... 14 设计四 模拟线性调制解调系统................................................................................. 18 设计五 2FSK调制解调系统 ...................................................................................... 26 设计六 2PSK和2DPSK调制解调系统 .................................................................... 32 设计七 数字通信系统的抗噪性能分析..................................................................... 37 设计八 载波同步......................................................................................................... 43 设计九 信道编码和译码............................................................................................. 46 参考文献....................................................................................................................... 50
3
设计一 随机信号分析
一、设计目的
1． 利用MATLAB产生各种随机信号。
2． 利用MATLAB计算随机信号的自相关函数和功率谱密度。 3． 掌握随机信号的自相关函数与功率谱密度的关系及其推导。
二、设计原理
1．随机信号的产生
在MATLAB中，提供了大量的随机数发生函数，其中rand和randn是两类核心函数，rand函数产生在[0，1]区间服从均匀分布的随机信号；randn函数产生服从均值为0，方差为1的高斯分布的随机信号。其它类型的随机信号可以通过这两类随机信号变换而得到。
(1) rand函数 格式：X=rand(N)
X=rand(M,N) (2) randn函数 格式：X=randn(N)
X=randn(M,N) 2． 随机信号的自相关函数
在MATLAB中，xcorr函数用于估计随机序列的自相关函数，其调用格式如下： (1) c=xcorr(x)可计算出矢量x的自相关函数，矢量x表示序列x(n)。 (2) c=xcorr(x,'option')中的option可用来指定相关选项：
·当option=biased时，xcorr函数可计算自相关函数的有偏估计，即
1Rx(m)?NN?m?1?xn?0n?mxn
·当option=unbiased时，xcorr函数可计算自相关函数的无偏估计，即
1Rx(m)?N?m即m?0的自相关函数值归一化为1。
N?m?1?xn?0n?mxn
·当option=coeff时，xcorr函数对序列进行归一化处理，使零滞后的自相关函数为1.0，·当option=none时，即为缺省情况，不做归一化处理。
4
3． 随机信号的功率谱密度
随机信号?(t)的自相关函数与其功率谱密度之间互为傅立叶变换关系。即
?P?(?)??R(?)e?j??d?
?? R(?)?12?????P?(?)ej??d?
【例1-1】在（-1/2，1/2）内产生一均匀分布的随机数N=1000的离散时间序列，计算该序列的自相关，定义为
1N?mRx(m)?XnXn?m,m?0,1,...,M ?N?mn?1N1 ??XnXn?m,m??1,?2,...,?M
N?mn?m同时，用计算Rx(m)的离散傅里叶变换（DFT）求序列{xn}的功率谱。DFT定义为
DFT?Rx(m)??m??M?RMx(m)e?j2?fm?2?M?1??
解 产生序列{xn}，计算自相关和功率谱DFT?Rx(m)?的MATLAB源程序在下面给出。应该值得注意的是，自相关函数和功率谱都呈现出显著的波动，因此有必要在几次实现上对样本的自相关求平均。本程序在随机过程的10次实现上采用平均自相关得到自相关函数和功率谱。源程序如下：
echo on N=1000; M=50; t=0:1:50; s=-0.5:0.02:0.5; Rx_av=zeros(1,M+1); Sx_av=zeros(1,M+1);
for j=1:10; % 取10次实现的整体平均 X=rand(1,N)-1/2; % X在 -1/2 和 1/2之间取值 Rx=Rx_est(X,M); % 本次实现的自相关函数 Sx=fftshift(abs(fft(Rx))); % 本次实现的功率谱密度 Rx_av=Rx_av+Rx; % 自相关函数之和
5
Sx_av=Sx_av+Sx; % 功率谱密度之和 end;
Rx_av=Rx_av/10; % 自相关函数的整体平均 Sx_av=Sx_av/10; % 功率谱密度的整体平均 subplot(2,1,1); plot(t,Rx_av)
xlabel( '(a) 自相关函数' ); subplot(2,1,2); plot(s,Sx_av)
xlabel( '(b) 功率谱密度' );
%-------------------------------------------------------------------------------------------------------------- function [Rx]=Rx_est(X,M) % 估计随机序列X的自相关函数 % 计算出Rx(0),Rx(1),...,Rx(M) % Rx(m)代表实际中的Rx(m-1). N=length(X); Rx=zeros(1,M+1); for m=1:M+1 for n=1:N-m+1
Rx(m)=Rx(m)+X(n)*X(n+m-1); end;
Rx(m)=Rx(m)/(N-m+1); end;
程序运行结果如图1.1所示。
6
图1.1 例1-1的运行结果
三、设计内容
1．已知瑞利分布随机信号的概率密度函数为
?xx2?exp(?2) x?0 f(x)???22?? 0 x?0?用randn函数产生??2的瑞利分布随机变量。
提示：两个独立分布、均值为0、方差为?的高斯随机变量的平方和开根号所得的随机变量服从功率为2?的瑞利分布。
2．产生均值为0，方差为1的高斯随机信号，长度为0.5s，计算该信号的自相关函数，并绘制高斯信号及自相关函数的图形。
提示：计算机只能产生各种离散信号序列，对于模拟连续信号，本质上都是产生信号的采样序列。只要采样频率足够高，模拟就可以非常准确。在这里，假定信号的采样频率
22
2fs?1000Hz。
3．在（-1，1）内产生一均匀分布的随机数N=500的离散时间序列，计算该序列的自相关函数和功率谱密度，并绘制图形。
4．一个带限随机过程X(t)的功率谱为
?f?B?N02 P(f)??
0 f?B?? 计算它的自相关函数，并绘制自相关函数与功率谱密度的图形。
提示：X(t)的自相关函数Rx(?)?
?B?BN0j2?ftedf?N0BSa(2?B?)。 2
7
设计二 模拟信号的数字化
一、设计目的
1．掌握脉冲编码调制的基本原理。 2．理解均匀PCM与非均匀PCM。
3．利用MATLAB对抽样值进行A律13折线编码。
二、设计原理
在模拟信号数字化方式中，出现最早且应用最广泛的是脉冲编码调制，即PCM编码，它经过抽样（时间离散化）、量化（幅值离散化）和编码（多电平转换为二电平）3个步骤，将一个时间和幅值都连续变换的模拟信号变成二进制数字信号。
1．抽样
抽样是将模拟信号在时间上离散化。对于低通型信号，当抽样频率fs?2fH时，属于正常抽样，不会发生频谱混叠；当fs?2fH时，属于欠抽样，已抽样信号频谱发生混叠。下面通过举例来看如何利用MATLAB分析已抽样信号的频谱。
【例2-1】有一信号m(t)?50Sa2(200t)，当抽样频率fs1=100Hz时，绘出原始信号和抽样后信号的波形及其频谱。
解 由题意可知，原始信号的最高频率fH?64Hz，因此用fs1=100Hz抽样时，会发生频谱混叠，MATLAB实现的源程序如下：
clear
t0=10; % 定义时间长度 ts=0.001; % 抽样周期 ts1=0.01; % 欠抽样周期 fs=1/ts;
fs1=1/ts1; % 抽样频率 df=0.5; % 定义频率分辨力 t=[-t0/2:ts:t0/2]; % 定义时间序列 x=sin(200*t); % 定义信号序列 m=x./(200*t);
w=t0/(2*ts)+1; % 确定t=0的点
m(w)=1; % 确定t=0的点对应的信号值为1 m=m.*m; m=50.*m;
[M,mn,dfy]=fftseq(m,ts,df); % 傅里叶变换 M=M/fs;
8
f=[0:dfy:dfy*length(mn)-dfy]-fs/2; % 定义频率序列 % 欠抽样过程
t1=[-t0/2:ts1:t0/2]; % 定义抽样时间序列
x1=sin(200*t1); % 计算对应抽样序列的信号序列 m1=x1./(200*t1); % 计算Sa函数序列
w1=t0/(2*ts1)+1; % 由于除0产生了错误值，计算该值的标号 m1(w1)=1; % 将错误值修正 m1=m1.*m1; m1=50.*m1;
[M1,mn1,df1]=fftseq(m1,ts1,df); % 对已抽样信号进行傅里叶变换 M1=M1/fs1;
% 将一个周期的图像扩展到多个周期
N1=[M1,M1,M1,M1,M1,M1,M1,M1,M1,M1,M1,M1,M1]; f1=[-7*df1*length(mn1):df1:6*df1*length(mn1)-df1]-fs1/2; % 绘图 subplot(2,2,1) plot(t,m)
axis([-0.15,0.15,-1,50]);
xlabel('(a) 原始信号的时域波形'); subplot(2,2,2) plot(f,abs(fftshift(M))); axis([-500,500,0,1]);
xlabel('(b) 原始信号的频谱'); subplot(2,2,3) stem(t1,m1)
axis([-0.15,0.15,-1,50]);
xlabel('(c) 欠抽样信号的时域波形'); subplot(2,2,4)
plot(f1,abs(fftshift(N1))); axis([-500,500,0,1]);
xlabel('(d) 欠抽样信号的频谱');
%----------------------------------------------------------------------------------------------------------- function [M,m,df]=fftseq(m,ts,df)
9
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)]; df=fs/n;
程序运行结果如图2.1所示。
50403020100-0.1-0.0500.050.1(a) 原始信号的时域波形50403020100-0.1-0.0500.050.1(c) 欠抽样信号的时域波形0-5000(b) 原始信号的频谱5000.5110.80.60.40.20-5000(d) 欠抽样信号的频谱500
图2.1 例2-1运行结果
2．量化
在MATLAB中，quantiz函数可以实现信号的均匀量化，compand函数可以实现信号的压缩与扩张，二者结合起来，可以实现信号的非均匀量化。它们的调用格式如下。
(1) quantiz函数
① index=quantiz(sig,partition)，参数sig表示输入信号，参数partition表示分区向量，分区向量给出了量化间隔端点的向量，长度为M?1，将信号分成M个区域。根据向量partition对输入信号sig产生量化索引index，index的长度与矢量sig的长度相同。向量partition是由
10
若干个边界判断点且各边界点的大小严格按升序排列组成的实矢量。若partition的矢量长度小于M?1，则索引向量index中的每个元素的大小为[0，M?1]范围内的一个整数。量化方法如下：
·当sig≤partition(1)时，输出0；
·当partition(i)
② [index,quants]=quantiz(sig,partition,codebook)，参数codebook表示码本向量，存放每个partition的量化值。根据码本codebook产生量化索引index和信号的量化值quants。如果partition的矢量长度为M?1，那么codebook的矢量长度为M。
③ [index,quants,distor]=quantiz(sig,partition,codebook)，产生量化索引index、信号的量化值quants和量化误差distor。
(2) compand函数
① out=compand(in,Mu,v,'mu/compressor')，对输入向量in进行?律压缩，Mu为压缩参数，v为输入信号的峰值。
② out=compand(in,Mu,v,'mu/expander')，对输入向量in进行?律扩张。
③ out=compand(in,A,v,'A/compressor')，对输入向量in进行A律压缩，A为压缩参数，v为输入信号的峰值。
④ out=compand(in,A,v,'A/expander')，对输入向量in进行A律扩张。 【例2-2】比较均匀量化与非均匀量化的量化误差。
解 在计算均匀量化的量化误差时，直接调用quantiz函数求解；在计算非均匀量化的量化误差时，先对原始信号进行压缩，再进行均匀量化，最后把经过扩张得到的恢复信号与与原始信号比较，计算其量化误差。源程序如下：
Mu=255; % 定义压缩参数 t=-4:0.1:4; % 定义时间序列 sig=exp(t); % 定义原始信号 v=max(sig); % 原始信号的最大幅值
[index1,quants1,distor1]=quantiz(sig,0:floor(v),0:ceil(v)); % 对原始信号进行均匀量化 compsig=compand(sig,Mu,v,'mu/compressor'); % 对原始信号进行压缩 [index,quants]=quantiz(compsig,0:floor(v),0:ceil(v)); % 对压缩信号进行均匀量化 newsig=compand(quants,Mu,max(quants),'mu/expander'); % 对非均匀量化信号扩张 distor2=sum((newsig-sig).^2)/length(sig); % 非均匀量化的量化误差 distor1 distor2
程序运行结果为
11
distor1 = 0.5348 distor2 = 0.0397
3．A律13折线编码
下面将通过MATLAB程序实现A律13折线编码，并输出8位码。
【例2-3】设输入信号抽样值为+1270个量化单位，按照A律13折线特性编成8位码。 解 按照A律13折线的编码规则，用MATLAB实现的源程序如下： clear all x=+1270; if x>0 out(1)=1; else
out(1)=0; end
if abs(x)>=0 & abs(x)<16
out(2)=0;out(3)=0;out(4)=0;step=1;st=0; elseif 16<=abs(x) & abs(x)<32
out(2)=0;out(3)=0;out(4)=1;step=1;st=16; elseif 32<=abs(x) & abs(x)<64
out(2)=0;out(3)=1;out(4)=0;step=2;st=32; elseif 64<=abs(x) & abs(x)<128
out(2)=0;out(3)=1;out(4)=1;step=4;st=64; elseif 128<=abs(x) & abs(x)<256
out(2)=1;out(3)=0;out(4)=0;step=8;st=128; elseif 256<=abs(x) & abs(x)<512
out(2)=1;out(3)=0;out(4)=1;step=16;st=256; elseif 512<=abs(x) & abs(x)<1024
out(2)=1;out(3)=1;out(4)=0;step=32;st=512; elseif 1024<=abs(x) & abs(x)<2048
out(2)=1;out(3)=1;out(4)=1;step=64;st=1024; else
out(2)=1;out(3)=1;out(4)=1;step=64;st=1024; end
12
if abs(x)>=2048
out(2:8)=[1 1 1 1 1 1 1]; else
tmp=floor((abs(x)-st)/step);
t=dec2bin(tmp,4)-48; %函数dec2bin输出的是ASICC字符串，48对应0 out(5:8)=t(1:4); end
out=reshape(out,1,8) 程序运行结果为 out =
1 1 1 1 0 0 1 1
三、设计内容
1．原始信号为m(t)?sin200?t?2cos100?t，现分别用3种频率对其进行抽样，
fs1=100Hz，fs2=200Hz，fs3=300Hz，绘出原始信号和抽样后信号的波形及其频谱。
2．利用MATLAB中的Simulink模块进行抽样与恢复过程的仿真，当输入信号为一频率为10Hz的正弦波时，观察对于同一输入信号有不同输入频率时，已抽样信号及恢复信号的不同形态。
提示：仿真中用到的主要模块有“Sine Wave”模块、“Pulse Generator”模块、“Analog Filter Design”模块和“Gain”模块。
3．对一个正弦信号m(t)?sint进行均匀量化，在图上同时显示出原始信号和量化后信号，当增大量化级时，观察量化级与量化误差之间的关系。
4．量化间隔端点的向量partition取值为0、1、3，每个区间的取值为-1、0.5、2、3，输入量化的离散信号为[-2.5，-1，-0.4，0，0.3，1，1.3，1.9，2，2.6，3，3.5，4.5]，先按照规则进行计算，再用MATLAB中的quantiz函数进行验证，最后再用“Scalar Quantizer”模块构造一个量化模型，观察量化前后取值变化。
5．输入向量in为1～10，输入信号的峰值为10，进行?=255的?律压缩与扩张，输出压缩后的向量和扩张后的向量。
6．画出A律13折线近似的压缩特性曲线和A=87.6时的压缩特性曲线，并将二者进行比较。
7．设输入信号抽样值为-350个量化单位，按照A律13折线特性编成8位码。
13
设计三 数字基带传输系统
一、设计目的
1．利用MATLAB画出数字基带信号的波形图。 2．利用MATLAB画眼图。
3．利用MATLAB分析无码间干扰基带系统的抗噪性能。 4．掌握无码间干扰系统的频谱特性。
二、设计原理
1．数字基带信号
数字基带信号就是消息代码的电波形，它是用不同的电平或脉冲来表示相应的消息代码。数字基带信号的波形和码型很多，最常用的由矩形脉冲组成的基带信号有：单极性归零及不归零波形，双极性归零及不归零波形，差分波形和多电平波形等。
【例3-1】用单极性归零码来表示信息序列100110000101，画出波形示意图。 解 MATLAB源程序如下： t=[1 0 0 1 1 0 0 0 0 1 0 1]; srz(t);
%--------------------------------------------------------------------------------------------------------------- function y=srz(x)
% 本函数实现将输入的一段二进制代码变为相应的单极性归零码输出 % 输入x为二进制码，输出y为编好的码 grid=300;
t=0:1/grid:length(x); % 定义对应的时间序列 for i=1:length(x) % 进行码型变换 if(x(i)==1) % 如果输入信息为1 for j=1:grid/2
y(grid/2*(2*i-2)+j)=1; % 定义前半时间值为1 y(grid/2*(2*i-1)+j)=0; % 定义后半时间值为0 end else
for j=1:grid/2 % 反之，输入信息为0 y((i-1)*grid+j)=0; % 定义所有时间值为0 end
14
end end
y=[y,x(i)]; % 给序列y加上最后一位 M=max(y); m=min(y); plot(t,y);
axis([0,i,m-0.1,M+0.1]); 程序运行结果如图3.1所示。
10.80.60.40.20024681012
图3.1 单极性归零信号波形
2．无码间干扰系统
满足无码间干扰时、频域条件的基带系统有3类，分别是理想LPF系统、滚降系统和部分响应系统。前两个系统的理论基础是奈奎斯特第一准则，第三个系统的理论基础是奈奎斯特第二准则。理想LPF系统可达到理论最大频带利用率，但是难以实现，且对位定时精度要求高。滚降系统可实现，且对位定时精度要求降低，但频带利用率降低。部分响应系统兼具了前两个系统的优点，即频带利用率高和低位定时精度，但其可靠性降低。
3．基带信号的眼图
眼图是指利用实验的方法估计和改善基带系统性能时，在示波器上所观察到的像人的眼睛一样的图形。
① MATLAB函数
在MATLAB中，eyediagram函数用来绘制眼图，其调用格式如下： eyediagram(x,n,period,offset,plotstring)
其中x是信号；n是每个轨迹包括的采样点数；period是指水平轴的坐标范围，即[-period/2,period/2]；offset是偏置因子，信号的第（offset+1）个采样点之后每n个值为一周期，且该周期为period的整数倍，offset必须是非负整数，其范围是[0,n-1]；plotstring是绘制眼图时采用的符号、线形和颜色，其格式与plot函数相同，如不设置，采用系统缺省值。
【例3-2】产生一个二进制随机方波信号，绘出通过升余弦滤波器后，方波的高频分量被滤掉后的眼图。
15
解 源程序如下：
x=randint(3000,1,2); % 产生长度为3000的二进制随机序列 y=[0;rcosflt(x,1,10)]; % x通过一个升余弦滤波器得到y t=1:length(y); % 定义时间序列 figure(1)
plot(t,y); % 绘出y的时域波形 axis([1,300,-0.5,1.5]) grid on
eyediagram(y,20,4); % 绘出y的眼图 程序运行结果如图3.2所示。
Eye Diagram 1.5
1.511Amplitude501001502002503000.50.5
00-0.5-0.5
-2-1.5-1-0.50Time0.511.52(a) 时域波形图 (b) 眼图
图3.2 通过升余弦滤波后的二进制数据流图形
② Simulink模块
在Simulink模块库中，显示眼图的模块为“Discrete-Time Eye Diagram Scope”，图形及参数设置界面如图3.3所示。
图3.3 眼图模块及其参数设置
16
4．无码间干扰基带系统的抗噪性能
在0、1等概的条件下，双极性信号和单极性信号的误码率计算公式分别如下：
Pe双?1Aerfc() 22?n1Aerfc() 222?nP?e单可以看出，在相同条件下，双极性误码率比单极性低，抗噪性能好。
三、设计内容
1．当输入二元信息序列为101001110001时，画出单极性不归零码和双极性不归零码的波形示意图，当占空比分别为1/2、1/3时，画出单极性归零码和双极性归零码的波形示意图。
2．当输入的二元信息序列为100101011时，画出差分码的波形示意图。 3．当输入的二元信息序列为100110000101时，画出数字双相码的波形示意图。 4．画出a=0，0.5，1的升余弦滚降系统的频谱图，并画出其各自对应的时域波形。 5．在例3-2中，利用Simulink建模，画出眼图。
6．根据单极性信号和双极性信号的误码率计算公式，作图比较两种信号的抗噪性能。 7．设有双极性不归零数字基带信号an，码元周期为Ts，g(t)??性高斯白噪声的双边功率谱密度为N0/2?0，画出眼图。
(1) 经过理想低通H(f)???1 0?t?Ts，设加
?0 其他?1 f?5/(2Ts)?后的眼图。
??0 其他??1 f?1/(Ts)(2) 经过理想低通H(f)??后的眼图。
0 其他??
17
设计四 模拟线性调制解调系统
一、设计目的
1．掌握模拟线性调制解调原理及MATLAB编程实现方法。
2．利用Simulink设计模拟线性调制解调系统，观察模拟线性调制信号的时域波形和频谱图。
二、设计原理
1．模拟线性调制原理
模拟线性调制是指调制后的信号频谱为调制信号频谱的平移或线性变换，而非线性调制则没有这个性质。下面介绍的模拟信号的连续波线性调制，主要指幅度调制，它包括常规幅度调制、抑制双边带调幅、单边带调幅和残留边带调幅。幅度调制器的一般模型如图4.1所示。其输出已调信号的时域和频域一般表达式为
sm(t)??m(t)cos?ct??h(t) （4.1） Sm(?)?
m(t)1?M(???c)?M(???c)?H(?) （4.2） 2?cos?ctSm(t)H(?)Sm'(t)
图4.1 幅度调制器的一般模型
式中，M(?)为基带调制信号m(t)的频谱；H(?)?h(t)；?c为载波频率。
在该模型中，根据选择的滤波器特性H(?)及m(t)的频谱成分不同，便可以产生各种幅度调制信号。
(1) 调幅信号（AM）
在图4.1中，若调制信号m(t)外加直流A0，滤波器H(?)=1为全通网络，则输出为调幅（AM）信号：
sAM(t)??A0?m(t)?cos?ct?A0cos?ct?m(t)cos?ct （4.3）
SAM(?)??A0??(???c)??(???c)??1?M(???c)?M(???c)? （4.4） 2从波形上看：当满足条件m(t)max?A0时，AM信号的包络与输入基带信号m(t)成正
18
比，故AM信号的解调可采用包络检波。
从频谱上看：AM的频谱由频频分量和上、下对称的两个边带组成，因此 BAM?2f其中，fHH （4.5）
是调制信号的最高频率。
(2) 双边带信号（DSB-SC）
若m(t)中无直流分量，H(?)为全通网络，则输出为抑制载波的双边带调幅信号： sDSB(t)?m(t)cos?ct （4.6） SDSB(?)?1?M(???c)?M(???c)? （4.7） 2从波形上看：DSB信号的包络不与m(t)成正比，故不能进行包络检波，需采用相干解
从频谱上看：DSB信号虽节省了载波发射功率，但仍具有上、下两个边带，频带宽度与
调。
AM信号相同。由于这两个边带是完全对称的，它们携带了调制信号的全部信息，故传输其中一个边带即可，这种方式就是单边带调制。
(3) 单边带信号（SSB）
单边带调制中只传输双边带信号中的一个边带，因此产生SSB信号的最直观方法是将图4.1中的滤波器设计成截止频率为fc的理想低通或理想高通特性。
SSB信号的时域表示式为
11?(t)sin?ct （4.8） m(t)cos?ct?m22?(t)是m(t)的希尔伯特变换。若M(?)为m(t)的傅式中，“+”为下边带，“-”为上边带。m?(?)为 ?(t)的傅氏变换M氏变换，则m sSSB(t)??(?)?M(?)??jsgn?? （4.9） M式中，sgn?为符号函数。
综上所述：SSB调制方式可节省载波发射功率，占用的频带宽度为BSSB?fH，只有AM、DSB的一半，即频带利用率提高一倍；SSB信号的解调仍需采用相干解调。
(4) 残留边带信号（VSB）
残留边带调制从频域上来看是介于SSB与DSB之间的一种调制方式，它保留了一个边带和另一个边带的一部分。它既克服了DSB信号占用频带宽，又解决了单边带滤波器不易实现的难题。
重要结论：为了保证相干解调时无失真恢复基带信号。必须要求残留边带滤波器的传输函数在载频处具有互补对称性。即
19
H(???c)?H(???c)?常数，???H， （4.10）
【例4-1】一有限长度信号m(t)，其表达式为：
?t 0?t?t０4?m(t)???t?t04 t04?t?3t04
?t?t 3t4?t?t 00?0将其调制在载波c(t)?cos2?fct上，假设t0?0.5，fc?50Hz，调制系数a?0.8，求出已调制信号的时域波形，未调信号和已调信号的频谱关系图，计算出已调信号和未调信号的功率，并且考虑有噪声的情况下，假设信噪比为10dB，求出噪声功率。
解 已调信号的时域表达式为：
u(t)??1?a?m(t)/mm?cos2?fct
因为调制时要进行归一化，所以用m(t)除以最大值mm。
未调信号的归一化功率为：Pmn?Pm/mm，Pm为未调信号的功率。 调制效率为：??2a2Pmn1?a2Pmn
已调信号的功率为：P1?amn(t)?/2 u?E?在给定信噪比条件下，可计算噪声功率为：Pn??Pu/10。源程序如下： t0=0.5; % 定义信号的持续时间 ts=0.001; % 定义采样时间 fs=1/ts; % 定义采样频率 df=0.2; % 频率分辨力
snr=10; % 定义信噪比，用dB表示 snr_lin=10^(snr/10); % 信噪比的数值 fc=50; % 定义载波频率
a=0.8; % 定义调制系数（Am/Ao） t=[0:ts:t0]; % 定义出采样点数据 % 定义信号m m=zeros(1,501);
for i=1:1:125 m(i)=i; end for i=126:1:375 m(i)=m(125)-i+125; end for i=376:1:501 m(i)=m(375)+i-375; end
20
m=m/1000;
c=cos(2*pi*fc.*t); % 载波信号 m_n=m/max(abs(m)); % 归一化 [M,m,df1]=fftseq(m,ts,df); % 傅里叶变换 M=M/fs; % 频率缩放，便于作图 f=[0:df1:df1*(length(m)-1)]-fs/2; % 定义频率向量
u=(1+a*m_n).*c; % 将调制信号调制到载波上 [U,u,df1]=fftseq(u,ts,df); % 对已调信号进行傅里叶变化 U=U/fs; % 频率缩放 signal_power=spower(u(1:length(t))); % 计算信号功率 pmn=spower(m(1:length(t)))/(max(abs(m)))^2; % 计算调制信号功率 eta=(a^2*pmn)/(1+a^2*pmn); % 计算调制效率 noise_power=eta*signal_power/snr_lin; % 计算噪声功率 noise_std=sqrt(noise_power); % 噪声标准差
noise=noise_std*randn(1,length(u)); % 产生高斯分布的噪声 r=u+noise; % 总接收信号，即加噪信号 [R,r,df1]=fftseq(r,ts,df); % 对总信号进行傅里叶变换 R=R/fs; % 频率缩放 % 以下为结果显示 signal_power pmn noise_power eta figure(1)
subplot(2,2,1); plot(t,m(1:length(t))); xlabel('时间'); title('调制信号'); subplot(2,2,3); plot(t,u(1:length(t))); xlabel('时间'); title('已调信号'); subplot(2,2,2); plot(f,abs(fftshift(M))); xlabel('频率'); title('调制信号的频谱'); subplot(2,2,4); plot(f,abs(fftshift(U))); xlabel('频率'); title('已调信号的频谱'); figure(2)
subplot(2,2,1); plot(t,noise(1:length(t))); xlabel('时间'); title('噪声'); subplot(2,2,3); plot(t,r(1:length(t))); xlabel('时间'); title('加噪信号'); subplot(2,2,2); plot(f,abs(fftshift(U))); xlabel('频率'); title('信号频谱'); subplot(2,2,4); plot(f,abs(fftshift(R))); xlabel('频率'); title('加噪信号频谱'); %---------------------------------------------------------------------------------------------------------------
21
function [M,m,df]=fftseq(m,ts,df) 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)]; df=fs/n;
%--------------------------------------------------------------------------------------------------------------- function p=spower(x) % 计算信号的功率 p=(norm(x)^2)/length(x); 程序运行结果如下 signal_power = 0.6074 pmn = 0.3327 noise_power = 0.0107 eta =
0.1755
调制信号0.20.030.0250.10.0200.0150.01-0.10.005-0.200.10.20.3时间已调信号20.40.40.50-5000频率已调信号的频谱500调制信号的频谱10.300.2-10.1-200.10.20.3时间0.40.50-5000频率500 图4.2 调幅信号及其频谱
22
噪声0.40.4信号频谱0.20.300.2-0.20.1-0.400.10.20.3时间加噪信号0.40.50-5000频率加噪信号频谱50020.410.300.2-10.1-200.10.20.3时间0.40.50-5000频率500 图4.3 加噪调幅信号及其频谱
2．模拟线性解调原理 (1) 包络检波法
对AM信号，当满足m0?m(t)'max时，不会发生过调制现象，此时用包络检波的方法
很容易恢复原始基带信号m(t)。如图4.4所示。
SAM(t)包络检波器Km(t)
图4.4 包络检波法
(2) 同步解调法
双边带信号不能用包络检波来解调，只能采用同步解调法，即相干解调，其原理如图4.5所示。
SDSB(t)?cos?ct图4.5 同步解调法
低通Km(t)
相干解调的关键是必须产生一个同频同相的本地载波，如果同频同相的条件得不到保证，则会破坏基带信号的恢复。所以，所需的设备比包络检波法复杂。
同理，单边带(SSB)信号和残留边带(VSB)信号也可采用同步解调的方法恢复出基带信号
m(t)。
3．MATLAB函数
23
在MATLAB中，modulate函数实现对信号的调制，demod函数是modulate函数的反过程，即实现对信号的解调。
(1) modulate函数
格式：y=modulate(x,Fc,Fs,'method') y=modulate(x,Fc,Fs,'method',opt) [y,t]=modulate(x,Fc,Fs,'method',opt)
其中，x是要调制的信号序列，Fc是载波频率，Fs是采样频率（采样频率Fs至少是载波频率Fc的2倍），method是采用的调制方法，opt是method指定的调制方法所需的参数（并不是每种方法都需要附加参数），t可以输出计算时间。表4.1列出method对应的调制方式。
表4.1 modulate中method所对应的调制方式 method amdsb-sc amdsb-tc或am
amssb fm pm pwm ppm qam
调制方式 双边带抑制载波调制 双边带带载波调制 单边带调制 频率调制 相位调制 脉冲宽度调制 脉冲相位调制 正交幅度调制
(2) demod函数
格式：x=demod(y,Fc,Fs,'method') x=demod(x,Fc,Fs,'method',opt)
其中，y是已调制信号，其余参数与modulate函数中的一样。
三、设计内容
1．已知未调制信号为
??sinc(200t) t?t0m(t)??
0 其余t?? 若t0取2s，载波为c(t)?cos2?fct，fc?100Hz，用抑制载波调幅来调制信号，画出调制信号和已调信号的时域波形及频谱图。
2．已知未调制信号为
2??[sinc(200t)] t?t0m(t)??
0 其余t?? 若t0取2s，载波为c(t)?cos2?fct，fc?100Hz，画出单边带调幅波的时域波形及频谱图。
24
3．若消息信号m(t)定义为
0?t?t03?1 ?m(t)??-2 t03?t?2t03
?0 其余t?将其调制在载波c(t)?cos2?fct上，假设t0?0.15s，fc?250Hz，调制系数a?0.85，分别采用AM、DSB和SSB三种调制方式，画出调制信号和已调信号的时域波形及频谱图。
4．对于第3题中的DSB信号采用同步解调法解调，画出解调信号的时域波形及频谱图，并与原始消息信号m(t)进行比较。
5．对于第3题中的AM信号采用包络检波法解调，画出解调信号的时域波形，并与原始消息信号m(t)进行比较。
6．对抑制载波的双边带调制解调系统进行仿真。设传输的基带信号为正弦波，幅度为1，频率范围为1～10Hz，载波频率为100Hz，传输信道为高斯白噪声信道，其信噪比为10dB，系统仿真采样频率设置为1000Hz。
提示：仿真中用到的主要模块有“Signal Generator”模块、“DSBSC AM Modulator Passband”模块、“DSBSC AM Demodulator Passband”模块和“AWGN Channel”模块。
25
设计五 2FSK调制解调系统
一、设计目的
1．掌握2FSK信号的调制解调原理及MATLAB编程实现方法。 2．利用Simulink设计2FSK信号的调制解调系统。 3．画出2FSK信号的时域波形和频谱图。
二、设计原理
1．2FSK信号的调制解调原理
数字移频键控是用载波的频率来传送数字消息，或者说用所传送的数字消息控制载波的频率。2FSK信号便是符号“1”对应于载频?1，而符号“0”对应于载频?2（与?1不同的另一载频）的已调波形，而且?1与?2之间的改变是瞬间完成的。从原理上讲，数字调频可用模拟调频法来实现，也可用键控法来实现，如图5.1所示。
载波
?f1开关Se0(t)载波s(t)e0(t)模拟调频器 (a) 模拟法?f2s(t)(b) 键控法
图5.1 2FSK信号的调制
根据以上对2FSK信号的调制原理的分析，已调信号的数字表达式可以表示为 e0(t)? ?ag(t?nT)cos(?t??)??ag(t?nT)cos(?t??) （5.1）
ns1nns2nnn其中，g(t)是单个矩形脉冲，脉宽为Ts，且
???0 概率为P?0 概率为?1?P?an?? an??
???1 概率为?1?P??1 概率为P2FSK信号的常用解调方法采用如图5.2所示的非相干检测法和相干检测法。这里的取样判决器是判定哪一个输入样值大，此时可以不专门设置门限电平。
26
?1带通滤波器y1(t)包络检波器v1(t)输入yi(t)?取样脉冲取样判决器输出s'(t)?2带通滤波器y2(t)包络检波器v2(t)
(a) 非相干方式
?1带通滤波器y1(t)乘法器z1(t)低通滤波器v1(t)输入yi(t)?cos?1t取样脉冲取样判决器?2带通滤波器y2(t)乘法器z2(t)低通滤波器v2(t)
cos?2t
(b) 相干方式 图5.2 2FSK信号的解调
2FSK调制属于非线性调制，其频谱特性的研究常用把2FSK信号看成是两个2ASK信号相叠加的方法。2FSK信号的功率谱密度为：
22Tsin?(f?f1)TsPE(f)?s[16?(f?f1)Tssin?(f?f2)Ts?(f?f2)Ts2sin?(f?f1)Ts??(f?f1)Ts2?
sin?(f?f2)Ts?]?
?(f?f2)Ts1[?(f?f1)??(f?f1)??(f?f2)??(f?f2)] （5.2） 16传输2FSK信号所需频带约为
B2FSK?f2?f1?2fs （5.3） 【例5-1】用MATLAB产生独立等概的二进制信源，画出2FSK信号的波形及其功率谱。 解 首先产生随机的二进制数字基带信号，然后根据2FSK信号的表达式产生二进制数字调制信号，最后通过FFT变换求解调制信号的功率谱。源程序如下：
A=1; % 调制信号幅度 fc=2; % 载波频率
27
N_sample=8; % 每码元的采样点数 N=500; % 码元数目 Ts=1; % 采样间隔 dt=Ts/fc/N_sample; % 波形采样间隔 t=0:dt:N*Ts-dt; % 定义时间序列 L=length(t); % 计算时间序列长度 % 产生二进制信源 d=sign(randn(1,N));
dd=sigexpand((d+1)/2,fc*N_sample); gt=ones(1,fc*N_sample); d_NRZ=conv(dd,gt);
[f,d_NRZf]=T2F(t,d_NRZ(1:L)); % 数字基带信号的傅里叶变换 sd_2fsk=2*d_NRZ-1; % 双极性基带信号 s_2fsk=A*cos(2*pi*fc*t+2*pi*sd_2fsk(1:L).*t); % 产生2FSK信号 [f,s_2fskf]=T2F(t,s_2fsk); % 2FSK信号的傅里叶变换 figure(1)
subplot(2,2,1); plot(t,d_NRZ(1:L)); axis([0,10,0,1.2]); xlabel('(a) 输入信号');
subplot(2,2,2); plot(f,10*log10(abs(d_NRZf).^2/Ts)); axis([-2,2,-40,40]); xlabel('(b) 输入信号功率谱密度'); subplot(2,2,3); plot(t,s_2fsk); axis([0,10,-1.2,1.2]); xlabel('(c) 2FSK信号');
subplot(2,2,4); plot(f,10*log10(abs(s_2fskf).^2/Ts)); axis([-fc-4,fc+4,-40,40]); xlabel('(d) 2FSK信号功率谱密度');
%--------------------------------------------------------------------------------------------------------------- function [out]=sigexpand(d,M)
% 将输入的序列扩展成间隔为N-1个0的序列 N=length(d); out=zeros(M,N); out(1,:)=d;
out=reshape(out,1,M*N);
%--------------------------------------------------------------------------------------------------------------- function [f,sf]=T2F(t,st)
% 利用FFT对信号进行傅里叶变换
28
% t,st分别为输入时间和信号 % f,sf分别为输出频率和信号频谱 dt=t(2)-t(1); T=t(end); df=1/T; N=length(st);
f=-N/2*df:df:N/2*df-df; sf=fft(st); sf=T/N*fftshift(sf);
程序运行结果如图5.3所示。
401200-200-40-20.505(a) 输入信号10-101(b) 输入信号功率谱密度210.50-0.5-105(c) 2FSK信号1040200-20-40-505(d) 2FSK信号功率谱密度
图5.3 2FSK信号的波形及其频谱
2．MATLAB函数
在MATLAB中，fskmod函数和fskdemod函数分别实现FSK信号的调制和FSK信号的非相干解调，awgn函数实现在信号中加入高斯白噪声，symerr函数实现错误符号数的统计和误符号率的计算。
(1) fskmod函数
格式：y=fskmod(x,M,freq_sep,nsamp) y=fskmod(x,M,freq_sep,nsamp,Fs)
其中，x是消息信号；M是消息的符号数，必须是2的整数次幂，消息信号是0～M-1之间的整数；freq_sep是两载波频率之间的频率间隔，单位为Hz；nsamp是输出信号y中每符号的采样数，必须是大于1的正整数；Fs是采样频率，freq_sep和M必须满足（M-1）*freq_sep<=Fs。
29
(2) fskdemod函数
格式：z=fskdemod(y,M,freq_sep,nsamp) z=fskdemod(y,M,freq_sep,nsamp,Fs)
其中，y是已调信号；其余参数与fskmod函数中的一样。 (3) awgn函数 格式：y=awgn(x,snr) y=awgn(x,snr,sigpower)
y=awgn(…,powertype)
其中，x是信号；snr是信噪比，单位为dB；sigpower是信号功率；powertype指定信噪比snr的单位。
(4) symerr函数
格式：[number,ratio]=synerr(x,y)
其中，x、y是待比较信号，二元序列；number是x、y对应元素比较后不相同的符号个数；ratio是误符号率。
3．Simulink模块
在Simulink模块库中，用于实现FSK信号调制解调的模块分别是“M-FSK Modulator Baseband”和“M-FSK Demodulator Baseband”，图形及参数设置如图5.4所示。
图5.4 FSK信号调制解调模块及参数设置
30
三、设计内容
1．对二元序列01101010，画出2FSK信号的波形，其中发“0”码时，载频与码元速率相等；发“1”码时，载频是码元速率的2倍。
2．利用MATLAB提供的fskmod函数实现2FSK调制，消息信号是1000个二进制的随机数序列，两载波频率之间的频率间隔为8Hz，已调制信号中每符号的采样数为8，采样频率为32Hz，画出2FSK信号的时域波形和频谱图。
3．利用MATLAB提供的FSK调制、解调、误码率分析函数实现2FSK信号的解调与检测，消息信号是5000个二进制的随机数序列，两载波间的频率间隔为8Hz，已调制信号中每符号的采样数为17，采样频率为16Hz，信道中的噪声为加性高斯白噪声，当系统信噪比为5dB时，求误符号率。
4．利用Simulink中的“M-FSK Modulator Baseband”模块和“M-FSK Demodulator Baseband”模块建立2FSK调制解调系统模型，并观察调制、解调前后的波形，计算差错率。
5．根据图5.1分别建立模拟调频法和键控法产生2FSK信号的系统模型，观察调制前后的信号波形。
6．根据图5.2分别建立2FSK信号的相干解调和非相干解调解调系统模型，观察解调前后的信号波形。
31
设计六 2PSK和2DPSK调制解调系统
一、设计目的
1．掌握2PSK和2DPSK调制解调原理及MATLAB编程实现方法。 2．利用Simulink设计2PSK和2DPSK信号的调制解调系统。 3．画出2PSK和2DPSK信号的时域波形和频谱图。
二、设计原理
1．2PSK和2DPSK调制解调原理
2PSK（二进制绝对移相键控）信号是利用载波相位直接表示数字信号，2DPSK（相对移相键控）信号是利用前后相邻载波的相位差来表示数字信号，2DPSK信号可以看作是将绝对码变成相对码后再进行绝对移相而得到的。2PSK和2DPSK信号的产生方法如图6.1所示。
双极性不归零极性变换器带通滤波器?an?差分编码器?bn?2DPSK信号cos?ct(a) 模拟调相法产生2DPSK信号
0相振荡器门电路1
e(t)加法器 ?相倒相器门电路2倒相器s(t)?(b) 相位选择法产生2PSK信号 图6.1 2PSK和2DPSK调制
2PSK和2DPSK信号的解调方法如图6.2所示。
e0(t)带通e0(t)z(t)低通抽样判决器x(t)输出?a'n?cos?ct(a) 2PSK信号的相干解调定时脉冲
32
2DPSK带通滤波器信号相乘器低通滤波器?bn?抽样判决器 码（反）变换器数据输出?an?本地载波(b) 2DPSK信号的解调（极性比较法）2DPSK带通滤波器信号
抽样判决器数据输出y1(t)相乘器z(t)低通滤波器x(t)?an?延迟Tsy2(t)(c) 2DPSK信号的解调（差分检测法）
图6.2 2PSK和2DPSK信号的解调
当0、1等概出现时，2PSK信号和2DPSK信号的功率谱密度为
?sin?(f?f)TTcSS? Pe(f)?4??(f?fc)TS?2PSK信号的带宽为
22sin?(f?fc)TS?? （6.1） ??(f?fc)TS?? B2PSK?2fs （6.2）
【例6-1】用MATLAB产生独立等概的二进制信源，画出2PSK信号的波形及其功率谱。 解 首先产生随机的二进制数字基带信号，然后根据2PSK信号的表达式产生二进制数字调制信号，最后通过FFT变换求解各种调制信号的功率谱。源程序如下，用到的子函数sigexpand.m和T2F.m见例5-1。
A=1; fc=2; N_sample=8; N=500; Ts=1; dt=Ts/fc/N_sample; t=0:dt:N*Ts-dt; L=length(t); % 产生二进制信源 d=sign(randn(1,N));
dd=sigexpand((d+1)/2,fc*N_sample); gt=ones(1,fc*N_sample); d_NRZ=conv(dd,gt);
[f,d_NRZf]=T2F(t,d_NRZ(1:L)); % 数字基带信号的傅里叶变换 ht=A*cos(2*pi*fc*t); % 调制信号 d_2psk=2*d_NRZ-1; % 双极性基带信号 s_2psk=d_2psk(1:L).*ht; % 产生2PSK信号
33
[f,s_2pskf]=T2F(t,s_2psk); % 2PSK信号的傅里叶变换 figure(1)
subplot(2,2,1); plot(t,d_NRZ(1:L)); axis([0,10,0,1.2]); xlabel('(a) 输入信号');
subplot(2,2,2); plot(f,10*log10(abs(d_NRZf).^2/Ts)); axis([-2,2,-40,40]); xlabel('(b) 输入信号功率谱密度'); subplot(2,2,3); plot(t,s_2psk); axis([0,10,-1.2,1.2]); xlabel('(c) 2PSK信号');
subplot(2,2,4); plot(f,10*log10(abs(s_2pskf).^2/Ts)); axis([-fc-4,fc+4,-40,40]); xlabel('(d) 2PSK信号功率谱密度'); 程序运行结果如图6.3所示。
401200-200-40-20.505(a) 输入信号10-101(b) 输入信号功率谱密度210.50-0.5-105(c) 2PSK信号1040200-20-40-505(d) 2PSK信号功率谱密度
图6.3 2PSK调制波形及其频谱
2．MATLAB函数
在MATLAB中，pskmod函数和pskdemod函数分别实现PSK信号的调制和解调，scatterplot函数用于画数字调制信号的星座图。
(1) pskmod函数 格式：y=pskmod(x,M)
y=pskmod(x,M,ini_phase)
其中，x是消息信号；M是消息的符号数，必须是2的整数次幂，消息信号是0～M-1之间的整数；ini_phase是初始相位，单位为rad。
(2) pskdemod函数
34
格式：z=pskdemod(y,M)
z=pskdemod(y,M,ini_phase)
其中，y是已调信号；z是还原信号，其余参数与pskmod函数中的一样。 (3) scatterplot函数
scatterplot函数是画信号散点图函数，常用于画数字调制信号星座图，其调用格式如下： scatterplot(x)
x是信号，如果x是一个两列的实矩阵，该矩阵的第一列作为同相分量（习惯用二维空
Scatter plot间的X轴表示），第二列作为正交分量（习惯用二维空间的Y轴表示）；如果x是一个复矢量，Scatter plot1该复矢量的实部作为同相分量，虚部作为正交分量；信号x在X-Y平面（同相-正交平面）1的位置就是星座图。 0.8星座图可以在信号空间展示信号所处的位置，为系统的传输特性分析提供直观、具体的0.80.60.60.40.4显示结果。利用scatterplot函数画已调信号星座图需执行以下几步：
① 如果调制中使用的符号数是M，那么创建信号[0:M-1]。该信号表示对调制器所有可tuaaaturedrudrreQuaQ0.2 能的输入。0.2② 0使用适当的调制函数完成对信号的调制。 ③ 0对上步产生的输出信号使用scatterplot函数画图。 【例6-2】用scatterplot函数画2PSK和8PSK信号的星座图。 解 源程序如下： -0.6M1=2; x1=[0:M1-1]; -0.6M2=8; x2=[0:M2-1]; -0.8-0.4-0.2-0.4-0.2scatterplot(pskmod(x1,M1)) -0.8scatterplot(pskmod(x2,M2)) -1-1-0.5-10程序运行结果如图6.4所示。 In-Phase-1-0.502PSKIn-Phase10.80.60.50.51118PSK
Quadrature0.80.60.40.20-0.2-0.4-0.6-0.8-1
Quadrature0.40.20-0.2-0.4-0.6-0.8-1-1-0.50In-Phase0.51-1-0.5
图6.4 2PSK和8PSK信号的星座图
0In-Phase0.51
35
3．Simulink模块
在Simulink模块库中，用于实现2PSK信号调制解调的模块分别是“BPSK Modulator Baseband”和“BPSK Demodulator Baseband”，用于实现2DPSK信号调制解调的模块分别是“DBPSK Modulator Baseband”和“DBPSK Demodulator Baseband”，如图6.5所示。
图6.5 2PSK和2DPSK信号调制解调模块
三、设计内容
1．对二元序列11010011，画出2PSK和2DPSK信号的波形，假设0相表示数字信息“0”，，且载波频率是码元速率的2倍。 ?相表示数字信息“1”
2．利用MATLAB提供的pskmod函数和pskdemod函数实现2PSK调制解调，分别画出消息信号、已调信号和还原信号的时域波形，并利用symerr函数计算误符号率。
3．利用Simulink中的“BPSK Modulator Baseband”模块和“BPSK Demodulator Baseband”模块建立2PSK调制解调系统模型，并观察调制、解调前后的波形，计算差错率。
4．利用Simulink中“DBPSK Modulator Baseband”模块和“DBPSK Demodulator Baseband”模块，建立2DPSK调制解调系统模型，并观察调制、解调前后的波形，计算差错率。
5．根据图6.1分别建立模拟调相法产生2DPSK信号和相位选择法产生2PSK信号的系统模型，观察调制前后的信号波形。
6．根据图6.2分别建立2PSK和2DPSK信号的解调系统模型，观察解调前后的信号波形。
36
设计七 数字通信系统的抗噪性能分析
一、设计目的
1．了解蒙特卡罗仿真算法的思想。
2．利用蒙特卡罗仿真分析数字通信系统的抗噪性能。
二、设计原理
蒙特卡罗仿真算法又称随机性模拟算法，是通过计算机仿真来解决问题的算法，同时可以通过模拟来检验自己模型的正确性。它的基本思想是，为了求解数学、物理、工程技术以及管理等方面的问题，首先建立一个概率模型或随机过程，使它们的参数，如概率分布或数学期望等是所求问题的解；然后通过对模型或过程的观察或抽样试验来计算所求参数的统计特征，并用算术平均值作为所求解的近似值。对于随机性问题，有时还可以根据实际物理背景的概率法则，用电子计算机直接进行抽样试验，从而对问题解答。在通信系统误码率计算中，由于计算公式复杂，甚至在很多情况下无法得到解析解，因此通过蒙特卡罗方法模拟实际的通信过程，得到仿真的通信系统误码率就成为一种方便的手段。
为了用蒙特卡罗方法估计通信系统的误比特率，让N个符号通过系统（实际上是系统的计算机仿真模型），并计算发送差错的个数Ne。如果在N次的符号发送中有Ne次差错，则比特差错率为
Pe?Ne （7.1） N下面给出通信系统的蒙特卡罗仿真模型，如图7.1所示。首先作以下假设：
噪声源数据源d?n?xd?n?发射机nd?n?yd?n?接收机??n?dxq?n?nq?n?噪声源yq?n?时延比较符号Pe图7.1 通信系统的蒙特卡罗仿真模型
37
① 在发射机中没有进行脉冲成形；
② 假设信道是加性高斯白噪声信道AWGN； ③ 信源输出端的数据符号是相互独立和等概率的； ④ 在系统中没有滤波处理，因此不存在码间串扰。
有了这些假设，该系统模型的分析和方针都及其简单。由于没有滤波，系统的时延为零，但在图7.1中仍然用虚线表示出时延模块，因为几乎所有的仿真都需要这个重要部分。下面通过举例说明蒙特卡罗仿真估计数字通信系统差错率的过程。
【例7-1】对一个使用相关器或匹配滤波器的二进制通信系统，用蒙特卡罗仿真估计Pe，并画出误码率Pe对信噪比SNR的图。该系统模型如图7.2所示。
均匀随机数发生器高斯随机数发生器n00/?二进制数据源r0检测器输出数据1/?n1r1高斯随机数发生器比 较差错计数器
图7.2 例7-1的仿真模型
解 误码率计算公式如下：
Pe?12??x???e2?x22?xdx?12????N0e?x22????dx?Q??N? （7.2） ?0?比值?/N0称为信号噪声比（SNR）。
先仿真产生随机变量r0和r1，它们构成了检测器的输入。首先产生一个等概率出现且互为统计独立的二进制0和1的序列。为了实现这一点，使用一个产生范围在（0，1）内的均匀随机数的随机数发生器，如果产生的随机数在（0，0.5）范围以内，二进制信源的输出就是0；否则输出是1。如果产生一个0，那么r0???n0，r1?n1；如果产生1，r0?n0，
r1???n1。
利用两个高斯噪声发生器产生加性噪声分量n0和n1，它们的均值是零，方差是
?2??N0/2。为了方便，可以将信号能量?归一化到1而改变?2。应该注意，这样SNR（定
义为?N0）就等于12?2。检测器的输出与二进制发送序列进行比较，差错计数器用来对
38
比特差错数计数。源程序如下：
SNRindB1=0:1:12; SNRindB2=0:0.1:12; for i=1:length(SNRindB1)
smld_err_prb(i)=smldPe(SNRindB1(i)); % 仿真误码率 end
for i=1:length(SNRindB2)
SNR=exp(SNRindB2(i)*log(10)/10);
theo_err_prb(i)=Qfunct(sqrt(SNR)); % 理论误码率 end
% 以下为绘图命令
semilogy(SNRindB1,smld_err_prb,'*'); hold
semilogy(SNRindB2,theo_err_prb); xlabel( '信噪比/dB'); ylabel( '误码率'); legend('仿真值','理论值');
%--------------------------------------------------------------------------------------------------------------- function [p]=smldPe(snr_in_dB)
% 计算以dB为单位的给定信噪比的误码率 E=1;
SNR=exp(snr_in_dB*log(10)/10); % 计算信噪比的数值 sgma=E/sqrt(2*SNR); % 噪声的方差 N=10000; % 符号数 % 产生二进制数据源 for i=1:N
temp=rand; % 产生一个区间在(0,1)的均匀随机变量 if (temp<0.5),
dsource(i)=0; % 如果产生的随机数在(0,0.5)以内，二进制源输出为0 else
dsource(i)=1; % 反之，二进制源输出为1 end end
% 判决，计算误码率
39
numoferr=0; for i=1:N
% 在判决器的接收端的信号，对于第i个符号为 if (dsource(i)==0) r0=E+gngauss(sgma);
r1=gngauss(sgma); % 如果信源输出为“0” else
r0=gngauss(sgma);
r1=E+gngauss(sgma); % 如果信源输出为“1” end % 判决如下 if (r0>r1)
decis=0; % 判决为“0” else
decis=1; % 判决为“1” end
if (decis~=dsource(i))
numoferr=numoferr+1; % 误比特数加1 end end
p=numoferr/N; % 误码率
%--------------------------------------------------------------------------------------------------------------- function grsv=gngauss(m,sgma)
% 函数生成两个统计独立的高斯分布随机数，以m为均值，sgma为方差 if(nargin==0) m=0; sgma=1; else if nargin==1 sgma=m;m=0; end end
u=rand; % 产生一个（0，1）间均匀分布的随机数u z=sgma*(sqrt(2*log(1/(1-u)))); % 利用上面的u产生一个瑞利分布的随机数 u=rand; % 重新产生（0，1）间均匀分布的随机数u
40
grsv(1)=m+z*cos(2*pi*u); grsv(2)=m+z*sin(2*pi*u);
%--------------------------------------------------------------------------------------------------------------- function [y]=Qfunct(x) % Q函数的计算 y=(1/2)*erfc(x/sqrt(2));
在几个不同的SNR值下，传输N?10000个比特时的仿真结果如图7.3所示。可以看出，仿真结果与式（7.2）给出的理论值Pe之间的一致性。也应该注意到，N?10000个数据比特
?3的仿真可以可靠地对差错概率估计在Pe?10以下，也就是说，用N?10000的数据比特，
在对Pe可靠的估计下应该至少有10个差错。
100仿真值理论值10-110误码率-210-310-410-50246信噪比/dB81012
图7.3 例7-1的运行结果
三、设计内容
1．观察仿真图7.3，说明仿真结果和理论计算结果在低信噪比下为什么情况完全一致，而在高信噪比下一致性稍差？思考应该如何改变仿真过程以使在高信噪比下也能得到更好的一致性？
2．实现一个2DPSK通信系统的蒙特卡罗仿真。设传输信号的相位是??0和??180，
????0?的相位变化对应传输0，??180?的相位变化对应传输1。要求完成在不同SNR参数?N0值下，传输N?10000个比特的仿真，画出测得的2DPSK的误码率及理论差错概率曲
线。
3．完成一个2FSK通信系统的蒙特卡罗仿真，其中信号波形由式（7.3）和（7.4）给出
41
u(t)?1 u2(t)?2?b cos(2?f1t)， ０?t?Tb （7.3）
Tb2?b（7.4） cos(2?f2t)， ０?t?Tb,
Tb式中，?b是每比特信号能量，Tb是码元周期，?f?f2?f1?1Tb。
4．完成M=4的PSK系统的蒙特卡罗仿真，信号波形由下式给出
2?m)， m?0，1，?，M-1 （7.5） um(t)?AgT(t)cos(2?fct?M其中gT(t)是发送滤波器的脉冲，它决定了传输信号的频谱特性，A是信号的幅度。5．完成M=4的DPSK系统的蒙特卡罗仿真。
42
设计八 载波同步
一、设计目的
1．掌握载波同步的方法，即直接法和插入导频法。 2．理解在DSB信号中插入正交载波作为导频信号。 3．利用Simulink设计载波提取系统，实现信号的解调。
二、设计原理
载波同步的方法有直接法（自同步法）和插入导频法（外同步法）两种。直接法不需要专门发送导频，而在接收信号中设法提取载波。插入导频法是在发送有用信号的同时，在适当的频率位置上插入正弦波作为导频，在接收端由导频提取出载波。
1．直接法
直接法适用于接收信号中具有载波分量，或接收信号（如DSB-SC信号、PSK信号）经过某种非线性变换后，具有载波的谐波分量的情况。
(1) 非线性变换——滤波法
这种方法是对不含有载波分量的信号先进行非线性处理，使其携带相应的载波分量，然后再用窄带滤波器或锁相环进行滤波，提取出载波分量。对于抑制载波的双边带调制信号常采用平方变换法和平方环法提取载波，如图8.1所示。
s(t)输入已调信号平方律部件e(t)2fc窄带滤波器 2fc二分频fc载波输出
(a) 平方变换法提取载波
s(t)输入已调信号平方律部件e(t)鉴相器（PD）环路滤波器（LF）压控振荡器（VCO）2fc二分频fc载波输出锁相环（PLL）
(b) 平方环法提取载波 图8.1 非线性变换——滤波法
(2) 同相正交环法
同相正交环又叫科斯塔斯（Costas）环，其原理如图8.2所示。由于加于两个乘法器相乘的本地载波分别是压控振荡器的输出信号和它的正交信号，因此，称这种环路为同相正交环。
同相正交环有两个主要的优点：一是当环路锁定后，同相之路的输出v5就是解调所要得到的信号。二是环路的工作频率是载波频率本身，而平方环的工作频率是载波频率的两倍，
43
因此当载波频率很高时，工作频率较低的同相正交环易于实现。因此，同相之路的乘法器兼有提取载波和相干解调的两种功能。同相正交环的缺点是电路较复杂，同时存在相位模糊问题。
同相支路v3v1低通滤波器v5m?(t)解调输出m(t)cos?ct输入已调信号载波输出90o相移压控振荡器环路滤波器v7v2正交支路v4低通滤波器v6
图8.2 同相正交环法提取载波
2．插入导频法
抑制载波的双边带信号本身不含有载波分量；残留边带信号虽含有载波分量，但很难从已调信号的频谱中把它分离出来；单边带信号更是不存在载波分量。对这些信号的载波提取，可以在发送端插入导频载波。插入导频的方法有两种，频域插入法和时域插入法。
(1) 频域插入导频法
插入导频的位置应该在信号频谱为零的位置，否则导频与信号频谱成分重叠在一起，接收时不易取出。对于抑制载波的双边带调制而言，在载频处，已调信号的频谱分量为零，同时对调制信号m(t)进行适当处理，就可以使已调信号在载频附近的频谱分量很小，这样便于插入导频以及解调时易于取出。图8.3所示为插入的导频和已调信号的频谱示意图。插入导频的频率就是fc，但它的相位与被调制载波正交，称为“正交载波”。在接收端提取这一导频，移相后作为相干载波。这种在频域插入导频的原理如图8.4所示。
S(f)0fc?fm导频fcfc?fmf
图8.3 抑制载波双边带信号的导频插入
m(t)调制信号带通滤波器uo(t)信道uo(t)带通滤波器v(t)低通滤波器m(t)解调输出acsin?ct90o相移acsin?ct发送端?accos?ctfc窄带 滤波器?accos?ct90o相移接收端
图8.4 插入导频法原理图
44
(2) 时域插入导频法
时域插入法是对被传输数据信号和导频信号在时间上加以区别，例如按图8.5那样分配。把一定数目的数字信号分作一组，称为一帧。在每一帧中，除有一定数目的数字信号外，在
t0～t1的时隙中传送位同步信号，在t1～t2的时隙内传送帧同步信号，在t2～t3的时隙内传送载波同步信号，而在t3～t4时间内才传送数字信息；以后各帧都如此。这种插入的结果只
是在每帧的一小段时间内才出现载波标准，在接收端用相应的控制信号将载波标准取出。
位同步帧同步载波同步信息位同步帧同步载波同步第二帧信息t0t1t2t3t4第一帧
图8.5 时域插入导频法
三、设计内容
1．构建一个抑制载波的双边带调制解调系统。载波频为10kHz，被调制信号为1kHz的正弦波，试用平方变换法恢复载波并进行解调。
2．用平方环法实现第1题中的双边带解调。 3．用科斯塔斯环法实现第1题中的双边带解调。
4．在图8.4所示的插入导频法中，假设输入调制信号为sin100?t，载波信号为
sin2000?t，画出接收端的解调输出信号的波形。
5．在图8.4所示的插入导频法中，假设输入调制信号为sin100?t，载波信号为
sin2000?t，如果将载波信号作为导频信号，即不经过90°相移，直接与已调信号相加后输
出，画出接收端的解调输出信号的波形。
45
设计九 信道编码和译码
一、设计目的
1．掌握信道编译码的基本原理。
2．利用MATLAB实现汉明码、线性分组码和循环码的编译码。
二、设计原理
信道编码的原理是在传输信息的同时加入信息冗余，通过信息冗余来达到信道差错控制的目的。信道编码一般可以分成两大类，即分组码和卷积码。分组码编码时将输入信息分成不同的组，对各组信息分别进行独立编码，加入冗余信息，组与组之间是独立的，其译码也是分组独立译码。卷积码编码时将输入信息与一固定结构的编码器进行卷积，卷积的输出作为传输信息。由于卷积的关系，卷积码的输出信息是前后关联的，因此译码时，卷积码一般采用序列译码的方式。
在MATLAB中，encode函数和decode函数实现信道编码和译码，hammgen函数和cyclgen函数分别产生汉明码和循环码的校验矩阵和生成矩阵，cyclpoly函数生成循环码的生成多项式，gen2par函数实现生成矩阵和校验矩阵的转换。信道编译码函数可完成6种主要的信道编译码：汉明码、线性分组码、循环码、BCH码、R-S码和卷积码，下面仅对这些函数在汉明码、线性分组码和循环码中的使用方法进行介绍。
1．encode函数
格式：code=encode(msg,N,K)
code=encode(msg,N,K,'method',opt) [code,added]=encode(...)
其中，msg是信息，信息位必须以二进制矢量或列矩阵的形式表示；method注明编码方式；N是码字长度；K是信息位的长度；opt是有些编码方式需要的参数，具体含义参见表9.1；added为使输入信息位数达到K时所需添加的列数。
表9.1 encode函数的参数用法
method hamming linear cyclic
含义 汉明编码 线性分组码 循环码
opt
可用来指定一个原始多项式，如省略，则使用默认多项式 必须指定一个校验矩阵 必须指定一个生成多项式
2．decode函数
格式：msg=decode(code,N,K)
46
msg=decode(code,N,K,'method',opt1,opt2,opt3,opt4) [msg,err]=decode(...) [msg,err,ccode]=decode(...) [msg,err,ccode,cerr]=decode(...)
decode函数对接收到的码字code，按method指定的方式进行译码，opt1，?，opt4是可选项参数，其用法如表9.2所示。
表9.2 decode函数的参数用法
method
译码方法
hamming 汉明译码，opt1可用来指定一个原始多项式，也可省略不用；opt2到opt4不用。 线性分组码译码，opt1用来指定一个校验矩阵；opt2用来指定一个检错逻辑电路，如
linear
省略，则默认单个错纠正逻辑；opt3和opt4不用。
循环码译码，opt1是必须指定的生成多项式，可使用cyclpoly函数选择合适的循环多
cyclic
项式；opt2指定检错逻辑，如省略，则采用默认的单个纠错逻辑；opt3和opt4不用。
译码器输出msg与码字code的格式匹配，当码字是一个N列矩阵时，输出msg信息以K列矩阵表示。当decode函数输入的码字与encode函数输出的格式不一样时，这个函数停止工作。
用[msg,err]=decode(...)可以输出译码过程检测出的错误数。当err为负数时，它表示纠错逻辑对差错程度无能为力。
用[msg,err,ccode]=decode(...)输出纠正的码字。
用[msg,err,ccode,cerr]=decode(...)输出code每行的错误数。 3．hammgen函数
hammgen函数产生汉明码的校验矩阵和生成矩阵，调用格式如下：
h=hammgen(M) h=hammgen(M,p) [h,g]=hammgen(...) [h,g,N,K]=hammgen(...)
其中，M是校验位的长度；p是原始多项式；h是校验矩阵；g是生成矩阵；N是码字长度；K是信息位长度。
例如，对M=3的校验阵和生成阵可以用[h,g]=hammgen(3)很容易得到： h =
1 0 0 1 0 1 1 0 1 0 1 1 1 0
47
0 0 1 0 1 1 1 g =
1 1 0 1 0 0 0 0 1 1 0 1 0 0 1 1 1 0 0 1 0 1 0 1 0 0 0 1 4．cyclgen函数
cyclgen函数产生循环码的校验阵和生成阵，调用格式如下：
h=cyclgen(N,p) h=cyclgen(N,p,'opt') [h,g]=cyclgen(...) [h,g,K]=cyclgen(...)
其中，N是码字长度；p是生成多项式；h是校验阵；g是生成阵；K是信息位长度；opt用于产生给定选项的校验阵。
·opt=nonsys时，生成一个非系统的循环校验阵。 ·opt=system时，生成一个系统的循环校验阵。 5．cyclpoly函数
cyclpoly函数可产生一个给定码长N和信息位长度K的生成多项式p，调用格式如下：
p=cyclpoly(N,K)
p=cyclpoly(N,K,'fd_flag')，可根据fd_flag寻找生成多项式。
·fd_flag=min时，寻找最小阶数循环生成多项式。 ·fd_flag=max时，寻找最大阶数循环生成多项式。 ·fd_flag=all时，对所有给定阶数的生成多项式进行搜索。 ·fd_flag=L时，搜索有L项的生成多项式。 如果没找到满足条件的多项式，p为空输出。 6．gen2par函数
gen2par函数是生成矩阵和校验矩阵的转换函数，调用格式如下： 格式：h=gen2par(g)
g=gen2par(p)
其中，g是生成矩阵，h是校验矩阵。
三、设计内容
1．（7，4）汉明码的校验矩阵为
48
1 1 0 1 0 0??1 ?
H??0 1 1 1 0 1 0???1 0 1 0 0 1??1 ?用MATLAB求（7，4）汉明码的码字。
2．假设（7，4）线性分组码的生成矩阵为
1 0 1 0 0 0??1 ?0 ?1 1 0 1 0 0? G???1 1 1 0 0 1 0???1 0 1 0 0 0 1??用encode函数和decode函数实现编译码。
3．给出（15，9）循环码的所有生成多项式。
4．求x?1的所有因子，构造（15，4）循环码的所有可能的生成多项式，挑出一个作为（15，4）的生成多项式，得到所有的码字。
5．设（7，4）循环码的生成多项式为g(x)?x3?x2?1，求出所有许用码组，并用encode函数和decode函数实现编译码，计算错误数。
15
49
参考文献
[1] 普罗克斯，萨勒斯．现代通信系统——使用MATLAB．西安：西安交通大学出版社．2001 [2] 刘敏，魏玲．MATLAB通信仿真与应用．北京：国防工业出版社．2001
[3] 李建新，刘乃安，刘继平．现代通信系统分析与仿真——MATLAB通信工具箱．西安：
西安电子科技大学出版社．2000
[4] 赵静，张瑾，高薪科．基于MATLAB的通信系统仿真．北京：北京航空航天大学出版
社，2007
[5] 郭文彬，桑林．通信原理——基于Matlab的计算机仿真．北京：北京邮电大学，2006 [6] 邵玉斌．Matlab/Simulink通信系统建模与仿真实例分析．北京：清华大学出版社．2008 [7] 罗卫兵，孙桦，张捷．SystemView动态系统分析及通信系统仿真设计．西安：西安电
子科技大学出版社，2001
[8] 樊昌信，曹丽娜．通信原理（第六版）．北京：国防工业出版社．2006 [9] 李白萍，吴冬梅．通信原理与技术．北京：人民邮电出版社．2003 [10] 李白萍，张鸣．通信原理上机实验指导．西安：西安科技大学，2006
50