基于MATLAB的DPSK通信系统仿真

目录

1 设计题目························································3 2 设计要求························································3 3 设计原理························································3

3.1 DPSK调制····················································3 3.2 DPSK解调····················································4

4 设计流程························································5 5 代码实现························································5

5.1 基带信号的获取···············································5 5.2 差分编码·····················································6 5.3 DPSK调制···················································7 5.4 AWGN信道···················································8 5.5 限噪处理·····················································8 5.6 差分相干解调·················································9 5.4 AWGN信道··················································11 5.7 功率谱密度图的绘制

6 测试结果·······················································12

6.1 调制解调正确性的测试········································12 6.2 系统抗噪性能的测试··········································14

7 总结····························································17 8 心得体会·······················································17 附：完整代码·····················································18

1

1 设计题目

用MATLAB仿真DPSK通信系统

2 设计要求

基本参数如下： 1) Fc=1800Hz 2) Rb=1200bps

3) 考虑加入AWGN信道，Eb/N0=20dB 4) 自定义一小段信息序列 结果要求：

1) 绘制各阶段信号波形 2) 绘制信号的PSD

3 设计原理

二进制差分相移键控常简称为二相相对调相，记为2DPSK。它不是利用载波相位的绝对数值传送数字信息，而是用前后码元的相对载波相位值传送数字信息。所谓相对载波相位是只本码元初相与前一码元初相之差。 3.1 DPSK调制

DPSK方式即是利用前后相邻码元的相对载波相位值去表示数字信息的一种方式。例如，假设相位值用相位偏移△φ表示（△φ定义为本码元初相与前一码元初相只差），并设

△φ=π→数字信息1 △φ=0→数字信息0

则数字信息序列与2DPSK信号的码元相位关系可举例表示如如下： 数字信息： 0 0 1 1 1 0 0 1 0 1 2DPSK信号相位： 0 0 0 π 0 π π π 0 0 π 或 π π π 0 π 0 0 0 π π 0 画出的2PSK及DPSK信号的波形如图3-1所示。

2

图3-1 2PSK及2DPSK信号的波形

2DPSK的产生基本类似于2PSK，只是调制信号需要经过码型变换，将绝对码变为相对码。2DPSK产生的原理框图如图2-2所示[1]，图（a）为模拟调制法，图（b）为键控法。

(a) 模拟调制法 (b) 键控法

图3-2 2DPSK信号的调制原理图

从上面分析可见，无论接收信号是2DPSK还是2PSK信号，单从接收端看是区分不开的。因此2DPSK信号的功率谱密度和2PSK信号的功率谱密度是完全一样的。 3.2 DPSK解调

DPSK信号可以采用相干解调法（极性比较法）和差分相干解调法（相位比较法）。本次设计采用差分相干解调法，图3-3为差分相干解调法原理图。

差分相干解调法的思想是直接比较前后码元的相位差。由于此时的解调已同时完成了码变换作业，故无需另加码变换器。由于这种解调方法又无需专门的相干载波，故事一种比较实用的方法。

3

图3-3 2DPSK差分相干解调解调原理图

4 设计流程

设计流程如下框图所示，箭头上的信号表示实际传递的信号名称。

bn dn mod modz

差分编码 DPSK调制 AWGN信道 BPF 对比 r Bn mt demod 抽样判决 LPF 延时相乘

差分相干解调

5 代码实现

5.1 基带信号的获取

利用MATLAB自带函数randint 产生定长度的二进制随机序列sign，为了将sign转换为单极性NRZ信号，要对sign的每个随机码进行T点增样，最后得到基带数字信号bn

代码实现如下：

%%%%%%%%%%%% 初始参量定义 %%%%%%%%%% fc=1800; %载波频率

Rb=1200; %基带数据传输速率 Tb=1/Rb;

T=50; %每个Tb内进行的增样点数 M=8; %原始码字序列的长度 N=9; %差分码字序列的长度

sign=randint(1,M); %产生长度为M的随机二进制序列 SNR=20; %信噪比

%%%%%%%%%%%% 对sign进行T点增样后的基带数字信号 %%%%%%%%%% bn=[]; for n=1:M

4

if sign(n)==0; F=zeros(1,T); else

F=ones(1,T); end

bn=[bn,F]; end

绘图

5.2 差分编码

设差分码cfm（二进制序列），差分编码原理为cfm[n] = b[n] xor cfm[n-1] 例如，如果 sign= [1 0 0 1 0 0 1 1 ],那么cfm=[1 0 0 0 1 1 1 0 1 ];

同时，还要对得到的差分码cfm进行T点增样，得到差分信号dn（单极性NRZ信号）。

代码实现如下：

%%%%%%%%%%%%% 差分编码得到差分码cfm %%%%%%%%%% cfm=[];

cfm(1)=1; %设cfm的初值为1 for i=1:M

cfm(i+1)=xor(sign(i),cfm(i)); end

%%%%%%%%%%%%% 对cfm进行T点增样 %%%%%%%%%%%%%% dn=[]; %差分信号 for n=1:N;

if cfm(n)==0;

A=zeros(1,T); %每个symbol进行T点增样 else

A=ones(1,T); end

dn=[dn A]; end

绘图

5

5.3 DPSK调制

DPSK调制由于PSK调制的主要原因是DPSK调制解决了“不确定性反相”问题，即DPSK信号通过相邻时隙载波相位的变化与否来“携带”信息。 在代码实现这一步的过程中，我假设了一个相邻时隙相位差变量delta，易知，当sign的一个码字为‘1’时，delta=pi，否则delta=0.通过delta就能反映前后时隙的载波相位变化，进而实现DPSK调制。

代码实现如下：

%%%%%%%%%%%%%%%% DPSK调制 %%%%%%%%%%%%%%%

t=Tb/T:Tb/T:Tb; %对Tb=1/Rb间隔进行T点划分 mod=cos(2*fc*pi*t); %调制后的信号 delta=0; %相位差 for n=1:M

if sign(n)==1;

delta=delta+pi; else

delta=delta; end

c2=cos(2*fc*pi*t+delta); mod=[mod c2]; end

绘图

5.4 AWGN信道

6

利用MATLAB自带函数awgn实现对已调信号mod进行加噪处理。

代码实现如下：

%%%%%%%%%%% 调制信号过高斯白噪声的信道 %%%%%%% modz=awgn(mod,SNR); %过AWGN信道

绘图

5.5 限噪处理

利用MATLAB自带函数fir1构造FIR I型 带通滤波器，使其对准信号频带，让信号几乎无失真通过，同时抑制带外噪声，此处取BPF的通带宽度为2*Rb，中心频率为fc。（由于设计的是数字滤波器，所以不要忘记对频率进行归一化处理）

值得注意的是，由于设计的滤波器是100阶，所以通过BPF后的信号会有出现接近于一个symbol宽度的时延。所以在信号modz接入BPF前，要对其进行一个拓宽处理，对modz最后一个symbol进行复制，得到modz_broad，再接入到BPF上。

同时，不要忘记对从BPF出来的信号r0进行截取，因为r0的第一个symbol是滤波器时延带来的。

代码实现如下：

%%%%%%%%%%% BPF %%%%%%%%%%%%%

modz_broad=[modz modz((N-1)*T+1:N*T)]; fs=Rb*T;

hBPF=fir1(100,2*pi*[fc-Rb fc+Rb]/fs); r0=filter(hBPF,1,modz_broad);

r=r0(T+1:(N+1)*T); %对r0进行截取

绘图

7

5.6 差分相干解调

差分相干解调的核心思想是比较两个相邻时隙上信号的相位，从而直接还原出信号Bn，分为三步： (1) 延时相乘

代码实现如下：

%%%%%%%%%%% 延时相乘 %%%%%%% r_shift=r(T+1:N*T); demod=[]; for i=1:M*T;

a=r_shift(i) .* r(i); demod=[demod a]; end

绘图

(2) LPF

利用MATLAB自带函数fir1构造FIR I型 低通滤波器，此处取LPF的通带宽度为Rb。

同BPF一样，由于设计的滤波器是100阶，所以通过LPF后的信号会有出现接近于一个symbol宽度的时延。所以要对信号demod接入LPF前，进行拓宽处理，得到demod_broad，再接入到LPF上。

同样，不要忘记对从LPF出来的信号mt0进行截取，因为mt0的第一个symbol

8

是滤波器时延带来的。

代码实现如下：

%%%%%%%%%%%%% LPF %%%%%%%%%%

demod_broad=[demod demod((M-1)*T+1:M*T)] hLPF=fir1(100,2*pi*[Rb]/fs); mt0=filter(hLPF,1,demod_broad);

mt=mt0(T+1:(M+1)*T); %对mt0进行截取

绘图

(3) 抽样判决

实现思想：对mt一个symbol的T个采样点进行求和，如果sum大于0，则Bn对应这个symbol内的所有点的值为1，反之为0.

代码实现如下：

%%%%%%%%%%%%% 抽样判决 %%%%%%%%%% Bn=zeros(1,M*T); for i=1:M;

if sum( mt( (i-1)*T+1 :i*T ) ) > 0; Bn((i-1)*T+1:i*T)=0; else

Bn((i-1)*T+1:i*T)=1; end end

绘图

9

从上图可知，最终经DPSK通信系统调制解调后得到的信号Bn与一开始的基带

数字信号bn完全一致。

5.7 功率谱密度图的绘制

采用Welch法绘制DPSK调制信号mod的功率谱密度图（Power spectrum magnitude)。选择适当的窗函数w(n)，并再周期图计算前直接加进去，加窗的优点是无论什么样的窗函数均可使谱估计非负。二是在分段时，可使各段之间有重叠，这样会使方差减小。本设计中尝试了三种窗，分别是矩形窗、海明窗和Blacman窗。

代码实现如下：

%%%%%%%%%%%%%%% PSD %%%%%%%%%%%%%%% window=boxcar(100); %矩形窗 window1=hamming(100); %海明窗

window2=blackman(100); %Blackman窗 Nfft=N*T*Rb;

noverlap=0; %数据无重叠

range='onesided'; %计算单边带功率谱密度

[Pxx,f]=pwelch(mod,window,noverlap,Nfft,fs,range); [Pxx1,f]=pwelch(mod,window1,noverlap,Nfft,fs,range); [Pxx2,f]=pwelch(mod,window2,noverlap,Nfft,fs,range);

figure(3) subplot(311)

plot(f,Pxx);axis([0 5000 0 0.0008]);

title('rectangular window PSD');xlabel('f/Hz');

subplot(312)

plot(f,Pxx1);axis([0 5000 0 0.0006]); title('Hamming PSD');xlabel('f/Hz');

subplot(313)

plot(f,Pxx2);axis([0 5000 0 0.0006]); title('Blackman PSD');xlabel('f/Hz');

10

绘图 x 10500x 105-4-4rectangular window PSD50010001500200025003000f/HzHamming PSD350040004500500000x 10-450010001500200025003000f/HzBlackman PSD35004000450050005005001000150020002500f/Hz30003500400045005000 从第一幅图即rectangular PSD中可以看出，单边带功率谱密度图的中心频率大致在1800Hz即载波频率处，第一零点大致在600Hz和3000Hz，即fc-Rb和fc+Rb处。

6 测试结果

6.1 调制解调正确性的测试

在前面我取用了一个特例sign=[1 0 0 1 0 0 1 1]来说明每一步的正确性。下面利用randint 产生定长度的二进制随机序列sign，多次运行结果都显示最终差分相干解调得出的信号Bn与基带数字信号bn一致，说明此次基于MATLAB设计的DPSK通信系统能够成功实现DPSK信号的调制与解调。

某几次运行结果的截图如下：

11

基带数字信号bn20-220-2差分信号dnDPSK调制信号mod20-2020040002004000200400过AWGN信道后的modz信号过BPF后的信号r220-210-10-2延时相乘后的信号demod10-1020040002004000200400过LPF后的信号mt经抽样判决后信号Bn20-202004000200400 基带数字信号bn20-220-2差分信号dnDPSK调制信号mod20-2020040002004000200400过AWGN信道后的modz信号过BPF后的信号r220-210-10-2延时相乘后的信号demod10-1020040002004000200400过LPF后的信号mt经抽样判决后信号Bn20-202004000200400 12

基带数字信号bn20-220-2差分信号dnDPSK调制信号mod20-2020040002004000200400过AWGN信道后的modz信号过BPF后的信号r220-210-10-2延时相乘后的信号demod10-1020040002004000200400过LPF后的信号mt经抽样判决后信号Bn20-2020040002004006.2 系统抗噪性能的测试 在前面调制解调正确性的测试中，设定信噪比SNR=20.下面通过检验在不同信噪比值下，系统调制解调的正确性，来简单评定系统的抗噪性能。 SNR=10时，正确解调 基带数字信号bn20-220-2差分信号dnDPSK调制信号mod20-2 020040002004000200400过AWGN信道后的modz信号过BPF后的信号r220-210-10-2延时相乘后的信号demod10-1020040002004000200400过LPF后的信号mt经抽样判决后信号Bn20-202004000200400

13

SNR=5时，正确解调 基带数字信号bn20-220-2差分信号dnDPSK调制信号mod20-2020040002004000200400过AWGN信道后的modz信号过BPF后的信号r220-210-10-2延时相乘后的信号demod10-1020040002004000200400过LPF后的信号mt经抽样判决后信号Bn20-202004000200400 SNR=2时，正确解调 基带数字信号bn20-220-2差分信号dnDPSK调制信号mod20-2 020040002004000200400过AWGN信道后的modz信号过BPF后的信号r220-210-10-2延时相乘后的信号demod10-1020040002004000200400过LPF后的信号mt经抽样判决后信号Bn20-202004000200400 SNR=1.5时，有误码出现 14

基带数字信号bn20-220-2差分信号dnDPSK调制信号mod20-2020040002004000200400过AWGN信道后的modz信号过BPF后的信号r220-210-10-2延时相乘后的信号demod10-1020040002004000200400过LPF后的信号mt经抽样判决后信号Bn20-202004000200400 所以，可以定性判断该系统的抗噪性能良好。

7 总结

经过两周的MATLAB课程设计，我收获多多。总结后主要有一下几点可谈

1) 良好编程习惯的培养

我们编程时要有良好的风格，源代码的逻辑简明清晰，易读易懂是好程序的重要标准。我的几点小收获是：标识符命名尽量含义鲜明，适当的注解，注意程序的视觉组织。

2) 信息筛查能力的训练

拿到一个课程设计题目后，应该从何下手，如何获取有效资源、筛查信息，是对个人能力的一种考验。经过这次试验后，我认为首先需要把理论知识弄透彻（尤其是课本上的），再阅读别人的设计报告，重复别人所做的工作是一种不错学习方法，但是需要懂得去粗取精，提取有效信息，否则会浪费大量时间。 3) 调试思想的建立

这次试验中我被很多问题难住，却常常不知如何下手找出错误，经老师指点后才懂得要按一定方法进行调试纠错，结果这些错误常常是一些拼写错误或符号不对称等小错。明确整体的架构，有针对地进行测试。调试方法多种多样，但调试思想的建立却是最为基础、根本。

8 心得体会

课程实验是培养学生综合运用所学知识,发现,提出,分析和解决实际问题,锻炼实践能力的重要环节,是对学生实际工作能力的具体训练和考察过程。

在这两周的时间里我切切实实的学到了许多知识，尤其是提高了在MATLAB中

15

的编程能力。回顾起此次通信原理课程实验，至今我仍感慨颇多，在整整两周的时间里，我的心情和课程设计的进展起起伏伏，一度怀疑自己的能力开始否定自己，直到最后把每处都弄通了、达到完美了，心情才真正舒畅开来。

但是，我确实学到很多很多的的东西，同时巩固了所学知识。整个设计的过程，可以说得是困难重重，这毕竟算是自己第一次独立完成，在设计的过程中发现了自己的不足之处，对软件不熟悉，动手能力比较欠缺。设计的过程很纠结，调试的过程也很痛苦，但是知识就是在这个过程中学习的，经验也是在这个过程中积累的。经历过困难才会有更多的收获。掌握正确的方法，保持认真的作风和平和的心态，这样才能真正的学到知识。

这次课程实验终于完成了，在设计中遇到了很多问题，多谢李晓峰老师的耐心讲解，从老师的身上我学到了很多实用的知识，也由衷地庆幸自己能遇到这样一位好老师，不仅专业知识扎实，而且品格高尚、谦和有礼，这就是我理想中的好老师。再次感谢李老师。同时，对给予我帮助的同学们也表示忠心的感谢。

附：完整代码

%%%%%%%%%%%% DPSK调制解调 %%%%%%%%%%%%%%%%%

%%%%%%%%%%%% 初始参量定义 %%%%%%%%%%%%% fc=1800; %载波频率

Rb=1200; %基带数据传输速率 Tb=1/Rb;

T=50; %每个Tb内进行的增样点数 M=8; %原始码字序列的长度 N=9; %差分码字序列的长度 sign=zeros(1,M);

sign=randint(1,N); %产生长度为M的随机二进制序列 SNR=1.5; %信噪比

%%%%%%%%%%%% 对sign进行T点增样后的基带数字信号 %%%%%%%%%% bn=[]; for n=1:M

if sign(n)==0; F=zeros(1,T); else

F=ones(1,T); end

16

bn=[bn,F]; end

%%%%%%%%%%%%% 差分编码得到差分码cfm %%%%%%%%%% cfm=[];

cfm(1)=1; %设cfm的初值为1 for i=1:M

cfm(i+1)=xor(sign(i),cfm(i)); end

%%%%%%%%%%%%% 对cfm进行T点增样 %%%%%%%%%%%%%% dn=[]; %差分信号 for n=1:N;

if cfm(n)==0;

A=zeros(1,T); %每个symbol进行T点增样 else

A=ones(1,T); end

dn=[dn A]; end

%%%%%%%%%%%%%%%% DPSK调制 %%%%%%%%%%%%%%%

t=Tb/T:Tb/T:Tb; %对Tb=1/Rb间隔进行T点划分 mod=cos(2*fc*pi*t); %调制后的信号 delta=0; %相位差 for n=1:M

if sign(n)==1;

delta=delta+pi; else

delta=delta; end

c2=cos(2*fc*pi*t+delta); mod=[mod c2]; end

%%%%%%%%%%% 调制信号过高斯白噪声的信道 %%%%%%% modz=awgn(mod,SNR); %过AWGN信道

%%%%%%%%%%% BPF %%%%%%%%%%%%%

modz_broad=[modz modz((N-1)*T+1:N*T)]; fs=Rb*T;

hBPF=fir1(100,2*pi*[fc-Rb fc+Rb]/fs); r0=filter(hBPF,1,modz_broad);

r=r0(T+1:(N+1)*T); %对r0进行截取

17

%%%%%%%%%%% 延时相乘 %%%%%%% r_shift=r(T+1:N*T); demod=[]; for i=1:M*T;

a=r_shift(i) .* r(i); demod=[demod a]; end

%%%%%%%%%%%%% LPF %%%%%%%%%%

demod_broad=[demod demod((M-1)*T+1:M*T)] hLPF=fir1(100,2*pi*[Rb]/fs); mt0=filter(hLPF,1,demod_broad);

mt=mt0(T+1:(M+1)*T); %对mt0进行截取

%%%%%%%%%%%%% 抽样判决 %%%%%%%%%% Bn=zeros(1,M*T); for i=1:M;

if sum( mt( (i-1)*T+1 :i*T ) ) > 0; Bn((i-1)*T+1:i*T)=0; else

Bn((i-1)*T+1:i*T)=1; end end

%%%%%%%%%%%% 画图 %%%%%%%%%%% %%%%%%%%%%%% figure(1) %%%%%%%%%% figure(1);

subplot(3,3,1); plot(bn);grid on;

axis([0 T*N -2 2]);title('基带数字信号bn');

figure(1);

subplot(3,3,2); plot(dn); grid on;

axis([0 T*N -2 2]);title('差分信号dn');

figure(1);

subplot(3,3,3); plot(mod); grid on;

axis([0 T*N -2 2]);title('DPSK调制信号mod');

18

figure(1);

subplot(3,3,4); plot(modz); grid on

axis([0 T*N -2 2]);

title('过AWGN信道后的modz信号');

figure(1);

subplot(3,3,5); plot(r);grid on

axis([0 T*N -2 2]); title('过BPF后的信号r')

figure(1);

subplot(3,3,6);

plot(demod);grid on axis([0 T*N -1 1]);

title('延时相乘后的信号demod')

figure(1);

subplot(3,3,7); plot(mt); grid on

axis([0 T*N -1 1]); title('过LPF后的信号mt');

figure(1);

subplot(3,3,8); plot(Bn); grid on

axis([0 T*N -2 2]);

title('经抽样判决后信号Bn')

%%%%%%%%%%%% figure(2) %%%%%%%%%% figure(2);

subplot(3,3,1);

plot(abs(fft(bn))); axis([0 T*N 0 50]); title('原始信号bn频谱');

figure(2);

subplot(3,3,2);

19

plot(abs(fft(mod))); axis([0 T*N 0 50]); title('DPSK信号mod频谱');

figure(2);

subplot(3,3,3);

plot(abs(fft(modz))); axis([0 T*N 0 50]);

title('加入白噪声的DPSK信号modz频谱');

figure(2);

subplot(3,3,4); plot(abs(fft(r))); axis([0 T*N 0 50]); title('过BPF后的信号r频谱');

figure(2);

subplot(3,3,5); plot(abs(fft(r))); axis([0 T*N 0 50]); title('信号r_shift频谱');

figure(2);

subplot(3,3,6);

plot(abs(fft(demod))); axis([0 T*N 0 50]); title('信号demod频谱');

figure(2);

subplot(3,3,7);

plot(abs(fft(mt))); axis([0 T*N 0 50]);

title('经低通滤波器后信号mt频谱');

figure(2);

subplot(3,3,8);

plot(abs(fft(Bn))); axis([0 T*N 0 50]);

title('经抽样判决后信号Bn频谱');

%%%%%%%%%%%%%%% PSD %%%%%%%%%%%%%%% window=boxcar(100); %矩形窗 window1=hamming(100); %海明窗

20

window2=blackman(100); %blackman窗 Nfft=N*T*Rb;

noverlap=0; %数据无重叠

range='onesided'; %频率间隔为[0 Fs/2]，只计算一半的频率 [Pxx,f]=pwelch(mod,window,noverlap,Nfft,fs,range); [Pxx1,f]=pwelch(mod,window1,noverlap,Nfft,fs,range); [Pxx2,f]=pwelch(mod,window2,noverlap,Nfft,fs,range);

figure(3) subplot(311)

plot(f,Pxx);axis([0 5000 0 0.001]); title('rectangular window PSD');

subplot(312)

plot(f,Pxx1);axis([0 5000 0 0.0007]); title('Hamming PSD');

subplot(313)

plot(f,Pxx2);axis([0 5000 0 0.0007]); title('Blackman PSD');

21

window2=blackman(100); %blackman窗 Nfft=N*T*Rb;

noverlap=0; %数据无重叠

range='onesided'; %频率间隔为[0 Fs/2]，只计算一半的频率 [Pxx,f]=pwelch(mod,window,noverlap,Nfft,fs,range); [Pxx1,f]=pwelch(mod,window1,noverlap,Nfft,fs,range); [Pxx2,f]=pwelch(mod,window2,noverlap,Nfft,fs,range);

figure(3) subplot(311)

plot(f,Pxx);axis([0 5000 0 0.001]); title('rectangular window PSD');

subplot(312)

plot(f,Pxx1);axis([0 5000 0 0.0007]); title('Hamming PSD');

subplot(313)

plot(f,Pxx2);axis([0 5000 0 0.0007]); title('Blackman PSD');

21

本文来源：https://www.bwwdw.com/article/7pxg.html