山东大学通信原理软件实验——MPSK 通信系统的设计与性能研究

本文为山东大学通信原理实验软件实验部分实验二的实验报告（2023年），因本人能力有限，仅供参考。
实验一、二完整版代码文件以及实验报告已上传至CSDN，感兴趣可以下载->山东大学通信原理实验软件部分-实验报告+代码-MATLAB-数字基带传输系统设计与性能研究-MPSK 通信系统的设计与性能研究

---

一、实验原理

1.1发送端理论分析

一组 M 载波相位调制信号波形的一般表达式为：
$${\mathrm{u}}_{\mathrm{m}}(\mathrm{t})=\mathrm{A}{\mathrm{g}}_{\mathrm{T}}(\mathrm{t})\cos (2\pi {\mathrm{f}}_{\mathrm{c}}\mathrm{t}+\frac{2\pi \mathrm{m}}{\mathrm{M}}),\mathrm{m}=0,1,\dots ,\mathrm{M}-1$$
其中，A 为信号幅度，M 为载波数目。$g_t(t)$为发送滤波器的脉冲波形，决定了传输信号的频谱特性。由于 PSK 信号对所有的 m 有相同的能量：
$${\epsilon }_m={\int }_{-\infty }^{\infty }{u}_m^2(t)dt={\epsilon }_s$$

其中 ${\epsilon }_s$ 为每个传输符号的能量。

所以在$g_T(t)=\sqrt{\frac{2}{T}},A=\sqrt{{\epsilon }_s}$ 时，信号波形可表示为：
$${\mathrm{u}}_{\mathrm{m}}(\mathrm{t})=\sqrt{\frac{2{\epsilon }_{\mathrm{s}}}{\mathrm{T}}}\cos (2\pi {\mathrm{f}}_{\mathrm{c}}\mathrm{t}+\frac{2\pi \mathrm{m}}{\mathrm{M}}),\mathrm{m}=0,1,\dots \dots \mathrm{M}-1$$

对式子中余弦函数进行展开，将其相角看作两个相角的和，可得：
$$\begin{array}{rlr}{\mathrm{u}}_{\mathrm{m}}(\mathrm{t})& =\sqrt{{\epsilon }_s}{\mathrm{g}}_{\mathrm{T}}(\mathrm{t})\cos \left(\frac{2\pi \mathrm{m}}{\mathrm{M}}\right)\cos (2\pi {\mathrm{f}}_{\mathrm{c}}\mathrm{t})-\sqrt{{\epsilon }_{\mathrm{s}}}{\mathrm{g}}_{\mathrm{T}}(\mathrm{t})\sin \left(\frac{2\pi \mathrm{m}}{\mathrm{M}}\right)\sin (2\pi {\mathrm{f}}_{\mathrm{c}}\mathrm{t})& \\ & =s_{\mathrm{mc}}{\phi }_1(\mathrm{t})+s_{\mathrm{ms}}{\phi }_2(\mathrm{t})& \end{array}$$

其中，$s_{\mathrm{mc}}=\sqrt{{\epsilon }_{\mathrm{s}}}\cos \left(\frac{2\pi \mathrm{m}}{\mathrm{M}}\right),s_{\mathrm{ms}}=\sqrt{{\epsilon }_{\mathrm{s}}}\sin \left(\frac{2\pi \mathrm{m}}{\mathrm{M}}\right),2{\epsilon }_{\mathrm{s}}={\epsilon }_g$ 。

${\phi }_1($t$)=\sqrt{\frac{2}{{\epsilon }_g}}{g}_T(t)cos(2\pi f_{c}t),{\phi }_2($t$)=-\sqrt{\frac{2}{{\epsilon }_g}}{g}_T(t)sin(2\pi f_{c}t)$为两正交基函数。
由于 PSK 信号 (除 2PSK 外)是二维信号，所以数字相位调制信号在几何上可以用二维向量表示：
$${\mathbf{s}}_{\mathbf{m}}=(\sqrt{{\epsilon }_s}\cos (\frac{2\pi \mathrm{m}}{\mathrm{M}}),\sqrt{{\epsilon }_s}\sin (\frac{2\pi \mathrm{m}}{\mathrm{M}}))$$

1.2接收机

接收机通常包含两部分：

信号解调器：将接收到的 r(t) 转换为 N 维矢量 r。实现方法有相关解调器和匹配滤波器型解调器。实验选择使用相关型解调器。

相关型解调器：
将接收信号 $r(t)=u_m(t)+n(t)$ 分别与 N 个 ${\phi }_i(t),i=1,2,\dots ,N$ 基函数作相关，用向量表示得：
$$\mathrm{r}={\mathrm{s}}_{\mathrm{m}}+\mathrm{n},{\mathrm{s}}_{\mathrm{m}}=({\mathrm{s}}_{\mathrm{m}}1,{\mathrm{s}}_{\mathrm{m}}2,\dots ,{\mathrm{s}}_{\mathrm{m}}\mathrm{N}),\mathrm{n}=({\mathrm{n}}_1,{\mathrm{n}}_2,\dots ,{\mathrm{n}}_{\mathrm{N}})$$

检测器：对矢量 r 进行分析并判断 M 信号中哪一个是发送的信号 $s_{\mathrm{m}}$ 。实现方法有 ML 检测器、MAP 检测器。实验中选择使用 ML 最佳检测器。
ML 最佳检测器：使错误概率最小的检测器。

判决依据：
1.最小欧式距离：求出接收到的信号向量分别与 M 个传输向量的欧式距离，选取对应的最小欧式距离的向量，该向量对应的符号为判决输出符号。
$$D(\mathbf{r},{\mathbf{s}}_{\mathbf{m}})=||\mathbf{r}-{\mathbf{s}}_{\mathbf{m}}|{|}^2$$

2.最大相关度量：将接收到的信号向量分别投影到 M 个传输信号向量上去，选取对应于最大投影的向量，该向量对应的符号即判决输出符号。
$$C(\mathbf{r},{\mathbf{s}}_m)=2\mathbf{r}\cdot {\mathbf{s}}_{\mathbf{m}}-||{\mathbf{s}}_m|{|}^2$$

将上述定义的距离度量展开：
$$D(\mathbf{r},{\mathbf{s}}_{\mathbf{m}})=|\mathbf{r}{|}^2-2\mathbf{r}\cdot {\mathbf{s}}_{\mathbf{m}}+|{\mathbf{s}}_{\mathbf{m}}{|}^2$$

其中，$|\mathbf{r}{|}^2$ 项对所有的判决度量是等价的，所以可以忽略。这样就得到相关度量。可以看出，距离度量越小，则相关度量越大。即相关度量与距离度量是完全等价的。

1.3接收端理论分析

信号在 AWGN 信道传输的时候, 会有噪声混杂进来, 此时在一个区间内收到的带通信号可以表示为：
$$r(t)=u_m(t)+n(t)=u_m(t)+n_c(t)cos(2\pi f_{c}t)-n_s(t)sin(2\pi f_{c}t)$$

其中 $n_c(t),n_s(t)$ 分别为加性噪声的同相分量和正交分量，是零均值互不相关的高斯随机过程。$E(n_c)=E(n_s)=0,E(n_c{n}_s)=0,{\sigma }^2(n_c)={\sigma }^2(n_s)=\frac{n_0}{2}$。
在相关型解调器中，r(t)分别与 ${\phi }_1$(t)、${\phi }_1$(t) 作相关，得到的输出为：
$$\mathbf{r}={\mathbf{s}}_{\mathbf{m}}+\mathbf{n}=\left(\sqrt{{\epsilon }_s}cos\left(\frac{2\pi m}{M}\right)+n_c,\sqrt{{\epsilon }_s}sin\left(\frac{2\pi m}{M}\right)+n_c\right)$$

选择最小欧式距离或者最大相关度量作为判决依据，得到判决输出 ${\mathbf{s}}_{\mathbf{m}}$'，可求差错概率。
下面讨论在采用最佳解调器和最佳检测器时，在 AWGN 信道上传输的 M 进制相位调制的错误概率。下面的推导得到的是理想情况下相位相干调制系统的性能。
假设发射信号 $u_0(t)$ 的相位 $\theta =0$，发射信号矢量为 $s_0=(\sqrt{{\epsilon }_s},0)$ 。
接收信号矢量的分量分别为：$r_1=\sqrt{{\epsilon }_s}+n_c,r_2=n_s$。因为$n_c$和$n_s$都是高斯随机变量，所以$r_1$和$r_2$也是高斯随机变量。且它们的均值分别为$\sqrt{{\epsilon }_{\mathrm{s}}}$ 和0，方差均为$\frac{N_0}{2}$ 。所以：
$$f_r(r_1,r_2)=\frac{1}{2\pi {\sigma }_r^2}{e}^{-[(r_1-\sqrt{{\epsilon }_s}{)}^2+r_2^2]/2{\sigma }_r^2}$$

进行变量代换：V$={\sqrt{\mathrm{r}}}_1^2+{\mathrm{r}}_2^2$ ，相干相位检测的度量为：${\Theta }_{\mathrm{r}}={\tan }^{-1}\frac{{\mathrm{r}}_2}{{\mathrm{r}}_1}$ ，带入到输出，得到关于随机变量 V 和 ${\Theta }_{\mathrm{r}}$ 的联合概率密度函数 $f_{V,{\Theta }_{\mathrm{r}}}(v,{\theta }_r)$ ，对 v 积分就可以得到 $f_{\mathrm{\Theta r}}({\theta }_r)$ ,即：
$$f_{{\Theta }_{\mathrm{r}}}({\theta }_r)=\frac{1}{2\pi }{e}^{-2{\rho }_{s}si{n}^2{\theta }_r}{\int }_0^{\infty }v{e}^{-(v-\sqrt{4{\rho }_s}cos{\theta }_r{)}^2/2}dv$$

当传输信号时，如果由于噪声干扰而使相位超出 $-\frac{\pi }{\mathrm{M}}\le {\Theta }_r\le \frac{\pi }{\mathrm{M}}$ 的范围，就会导致错误判决。由此可得码元错误概率为：
$$P_M=1-{\int }_{-\pi /M}^{\pi /M}{f}_{{\theta }_r}\text{(0)d}\Theta$$

一般情况下，$f_{{\theta }_r}$(θ) 的积分不能化为一个简单的函数，必须利用数值积分法来估计。但在 M=2 和 M=4 时，积分式可以写为简单函数的形式。对于二进制相位调制(M=2)，传输的信号为对极信号，因此其错误概率为：
$$P_{2PSK}=Q(\sqrt{\frac{2{\epsilon }_b}{N_0}})$$

M=4 的 QPSK 调制实际上可以看做是两个在正交载波上的 2PSK 调制系统，如果载波相位估计是准确的，则两对正交载波调制信号间就不会发生串扰现象，这时系统的比特错误概率等于 $P_M$。也即：QPSK 误比特率与 BPSK 的一样。
M=4 时，系统码元错误概率由下式决定：
$$P_4=1-P_c=1-(1-P_2{)}^2=2Q(\sqrt{\frac{2{\epsilon }_b}{N_0}})[1-\frac{1}{2}Q(\sqrt{\frac{2{\epsilon }_b}{N_0}}){]}^2$$
其中 $P_c$ 是二比特码元的正确判决概率。M>4 时，码元错误概率只能用数值积分法求出。此时可以求出大 M 值、大信噪比条件下 $f_{\mathrm{\Theta r}}({\theta }_r)$ 的近似值，进而可以得到错误概率的近似值。
$$\begin{array}{rl}{P}_{MPSK}& =2Q(\sqrt{\frac{2{\epsilon }_s}{N_0}}sin\frac{\pi }{M})\\ & =2Q(\sqrt{\frac{2(lo{g}_{2}M){\epsilon }_b}{N_0}}sin\frac{\pi }{M})\end{array}$$

1.4格雷码

格雷码（Gray code），又叫循环二进制码或反射二进制码。格雷码是一种 无权码，采用绝对编码方式，典型格雷码是一种具有反射特性和循环特性的单步自补码，它的循环、单步特性消除了随机取数时出现重大误差的可能，它的反射、自补特性使得求反非常方便。格雷码属于可靠性编码，是一种错误最小化的编码方式。它是一种数字排序系统，其中的所有相邻整数在它们的数字表示中只有一个数字不同。它在任意两个相邻的数之间转换时，只有一个数位发生变化。它大大地减少了由一个状态到下一个状态时逻辑的混淆。

---

二、系统框图

第一步：由均匀随机数发生器产生二进制信号序列。
第二步：通过 QPSK/8PSK 映射为 4/8 个信号点，采用格雷码映射。
第三步：加入高斯白噪声，产生一对相互独立的高斯噪声随机变量 $n_c(t)$ 和 $n_s(t)$ ，叠加在相互正交的两个信号支路上，产生相关接收机的两路输出 $r_c(t)$ 和 $r_s(t)$。
第四步：检测器判决。将接收机的两路输出 $r_c(t)$ 和 $r_{\mathrm{s}}(t)$ 按照最大投影准则或最小欧式距离准则判决成相应的信号点。
第五步：将检测器的输出与发送端的二进制信号序列进行比较，统计比特错误个数以及码元错误个数，并进行相应的 BER 和 SER 计算。
第六步：绘制星座图和误比特率曲线。

---

三、系统主函数

3.1系统代码流程图

3.2主程序代码

%  AWGN信道下QPSK和8PSK调制通信系统设计与性能研究
function main()
    Es=1;
    Eb=Es/2;
    L=input('请输入数据点数：');
    flag0=input('进行单指标噪声分析请输入0； 进行Monte Carlo仿真分析误比特率曲线请输入1\n');
    if flag0==0
        flag2=input('QPSK调制请输入0；    8PSK调制请输入1\n');
        if flag2==0
            M=4;
        else
            M=8;
        end
        [N,d]=GenerateSendSequence(L,M);
        disp('发送序列长度')
        disp(N)
        [gray_code,s]=Grayreflect(N,d,M);
        flag1=input('已知噪声方差请输入0；    已知信噪比请输入1\n');
        if flag1==0
            var=input('请输入噪声方差：');
        else
            SNR=input('请输入信噪比，单位dB：');
            var=SNR2var(SNR,Eb);
        end
        [rn]=ReciveSignal(s,var);
        flag3=input('选择最小距离法判决请输入0；     选择最大投影法判决请输入1\n');
        if flag3==0
            [Judge]=MinDistance(rn,M);
        else
            [Judge]=MaxProjection(rn,M);
        end
        [dnew]=Rebuild(Judge,M);
        [BER,SER]=ErrorRate(d,dnew,gray_code,Judge);
        fprintf('误比特率为%f\n',BER);
        fprintf('误符号率为%f',SER);

        x=0:1:size(d,2)-1;
        figure
        stem(x,d,MarkerSize=4);    axis([-1 size(d,2) -0.02 1.02])
        title('发送信号');
        xlabel('n')  
        ylabel('幅度')
        x=0:1:length(gray_code)-1;
        figure
        stem(x,gray_code,MarkerSize=4);    axis([-1 length(gray_code) -0.02 M-1+0.1])
        title('格雷码');
        xlabel('n')  
        ylabel('幅度')
        x=0:1:length(s(1,:))-1;
        figure
        stem(x,s(1,:),MarkerSize=4);    axis([-1 length(s(1,:)) -1.1 1.1])
        title('加噪声前发送信号的余弦分量');
        xlabel('n')  
        ylabel('幅度')
        x=0:1:length(rn(1,:))-1;
        figure
        stem(x,rn(1,:),MarkerSize=4);    axis([-1 length(rn(1,:)) -2.1 2.1])
        title('加噪声后发送信号的余弦分量');
        xlabel('n')  
        ylabel('幅度')
        x=0:1:size(dnew,2)-1;
        figure
        stem(x,dnew,MarkerSize=4);    axis([-1 size(dnew,2) -0.02 1.02])
        title('重建信号');
        xlabel('n')  
        ylabel('幅度')
        Planishphere(gray_code,rn,var,M);
    else
        flag2=input('QPSK调制请输入0； 8PSK调制请输入1;  QPSK和8PSK的对比请输入2\n');
        if flag2==0
            M=4;
            [N,d]=GenerateSendSequence(L,M);
            disp('发送序列长度')
            disp(N)
            flag=input('选择最小距离法判决请输入0,     选择最大投影法判决请输入1\n');
            Curve(N,d,M,flag);
        elseif flag2==1
            M=8;
            [N,d]=GenerateSendSequence(L,M);
            disp('发送序列长度')
            disp(N)
            flag=input('选择最小距离法判决请输入0,     选择最大投影法判决请输入1\n');
            Curve(N,d,M,flag);
        elseif flag2==2
            CurveAll(L);
        end
    end
end

本程序可以同时对 QPSK 和 8PSK 系统进行仿真，也即同时完成了实验内容（一）和实验内容（二）。本程序操作性好，交互性强，可以通过控制台输入对应的数字来选择进行不同类型的仿真，还可以自由选择调制类别、判决方法以及噪声量化方式。另外，程序的运行会输出 MPSK 通信系统仿真运行的各个阶段的图像，使系统运行过程更加清晰。

【思路】
首先从键盘输入要传输的数据点数 L，接着通过键盘输入 0 或 1 来分别选择进行单指标系统噪声分析还是进行 Monte Carlo 仿真分析。
若选择进行单指标系统噪声分析，则需继续选择系统的调制方式，输入 0 为 QPSK 调制，输入 1 为 8PSK 调制。选择好调制方式后，主函数会调用对应的子函数来生成随机的发送序列和补零操作，并在控制台输出发送序列的长度 N 。随即对发送序列进行格雷码映射。输入信道噪声的量化方式也是通过输入 0或 1 来分别选择输入噪声的方差或输入信噪比(dB)。为了方便计算，若输入的是信噪比，则会调用子函数来将信噪比转换为方差，随即在 AWGN 信道中传输调制信号，加入高斯白噪声。接下来则是选择判决方法，分别是最小距离判决发和最大投影判决，不同的判决方法会调用对应的子函数来执行判决。判决完成后会重构序列，对误比特数和误码数进行计数，并计算 BER 和 SER ，进行输出。最后会分别对发送序列、进行格雷码映射后的序列、加噪声前发送信号的余弦分量、加噪声后发送信号的余弦分量以及重建信号绘制图像，并画出眼图。
若进行 Monte Carlo 仿真分析，首先也是选择系统的调制方式，分为 QPSK调制、 8PSK 调制以及 QPSK 和 8PSK 的对比分析。选择好后会生成随机的发送序列并进行补零，在控制台输出发送序列的长度 N 。接着选择两种判决方式，随即调用 Monte Carlo 仿真子函数绘制理论误比特率/误符号率曲线以及仿真误比特率/误符号率曲线。

【功能验证】
点击运行后，可以在命令行窗口按提示输入实验所需的数值，如图 4 所示。此处以要传输的二进制比特数 L=1000，方差 ${\sigma }^2$ = 0.1 ，进行单指标噪声分析，调制方式为 QPSK ，使用最小距离判决法进行判决为例。

输入后经运行，得到的图像如表 2 所示。经验证，主函数完成了设计功能。

接下来验证 Monte Carlo 仿真分析的功能，这里以 L=1000， QPSK 调制 ，使用最小距离判决法进行判决为例。控制台运行结果如下图所示。

输入后经运行，得到的图像如下图所示。经验证，主函数完成了设计功能。

---

四、子函数实现

4.1发送序列产生子函数(GenerateSendSequence)

% 产生发送序列
function[N,d]=GenerateSendSequence(L,M)
    %L=1000;
    j=rand(1,L);%产生一行L列的随机序列，其值都在0-1内
    al=zeros(1,L);%产生一行L列的零序列
    for i=1:L  %将随机序列中小于0.5判为0；大于0.5判为1   
        if j(i)<0.5 
            al(i)=0;   
        else  
            al(i)=1;
        end
    end
    % 补零
    if M==4
        if (mod(L,2)~=0)% 当L为奇数时，在最后补一个0
            N=L+1;
            d=zeros(1,N);
            for j=1:L
                d(j)=al(j);
            end
        elseif(mod(L,2)==0)% 当L为偶数时，则不需要补0
            N=L;
            d=zeros(1,N);
            for j=1:L
                d(j)=al(j);
            end
        end
    elseif M==8
        if (mod(L,3)==0)% 当L为3的倍数时，不需要补0
            N=L;
            d=zeros(1,N);
            for j=1:L
                d(j)=al(j);
            end
        elseif(mod(L,3)==1)% 当L除3余1时，补两个0
            N=L+2;
            d=zeros(1,N);
            for j=1:L
                d(j)=al(j);
            end
        elseif(mod(L,3)==2)% 当L除3余2时，补一个0
            N=L+1;
            d=zeros(1,N);
            for j=1:L
                d(j)=al(j);
            end
        end
    elseif M==2
            N=L;
            d=zeros(1,N);
            for j=1:L
                d(j)=al(j);
            end
    end
end

1）编程思路：首先调用 rand 函数生成一行 L 列的随机序列，其值都在 0-1 内。对产生的随机序列进行判断，大于 0.5 判为 1，小于 0.5 判为 0，生成随机的长为 L 的二进制序列。接着根据系统的调制方式对序列进行补零。对于 QPSK 调制，由于两个比特表示一个符号，若 L 为偶数，则不需要补零；若 L 为奇数，则在最后补一个 0。对于 8PSK 调制，三个比特表示一个符号，则需要相应的将序列补为 3 的倍数。因为由于补充部分是确知的所以不会造成影响。
2）输入：要产生的序列的长度 L ，调制方式的载波数目 M 。
3）输出：发送序列 d，发送序列长度 N。
4）功能验证：在函数内赋值 L 为 49，M 为 4，画出 d 的 stem 图像，在函数内赋值 L 为 50，M 为 8，画出 d 的 stem 图像，如下图所示。经验证，函数完成了设计功能。

4.2格雷码映射子函数(Grayreflect)

%格雷码映射子函数
function [gray_code,s]=Grayreflect(N,d,M)
    %M=4;
    % 方法一：公式：f=|3x-y|可将相应的比特序列xy映射为格雷码
%     for i=1:2:N
%         gray_code((i+1)/2)=abs(d(i)*3-d(i+1));%将比特序列用格雷码映射转换成四进制数
%     end   
    
    % 方法二：一一映射 格雷码编码 00-00 01-01 10-11 11-10
    if M==4
        gray_code=zeros(1,N/2);
        for i=1:2:N-1% QPSK映射四进制格雷码
            if d(i)==0&&d(i+1)==0
                gray_code((i+1)/2)=0;
            elseif d(i)==0&&d(i+1)==1
                gray_code((i+1)/2)=1;
            elseif d(i)==1&&d(i+1)==1
                gray_code((i+1)/2)=2;
            elseif d(i)==1&&d(i+1)==0
                gray_code((i+1)/2)=3;
            end
        end
    %格雷码映射生成2维数组
    s=zeros(2,length(gray_code));
    s(1,:)=cos(2*pi*gray_code/M);
    s(2,:)=sin(2*pi*gray_code/M);
    elseif M==8
        gray_code=zeros(1,N/3);
        for i=1:3:N% 8PSK映射八进制格雷码
            gray_code((i+2)/3)=abs(d(i)*7-abs(d(i+1)*3-d(i+2)));
        end
    %格雷码映射生成2维数组
    s=zeros(2,length(gray_code));
    s(1,:)=cos(2*pi*gray_code/M);
    s(2,:)=sin(2*pi*gray_code/M);
    end 
end

1）编程思路：将发送序列映射为格雷码编码序列。QPSK 映射为 2 位格雷码，
8PSK 映射为 3 位格雷码。比特序列到格雷码的映射有两种方式，分别为公式法和一一映射的方法。公式法较为简单，但不够直观；一一映射的方法则较为复杂，但对应关系十分清晰。对于位数较高的格雷码映射来说，一一映射就显得繁琐了。在该函数中，2 位格雷码的映射使用一一映射的方法，3 位格雷码的映射则使用公式法。编码后，将序列映射为二维数组，包含同相分量和正交分量。
2）输入：发送序列长度 L、发送序列 d、调制方式的载波数目 M 。
3）输出：格雷码编码序列 gray_code、发送序列映射的二维数组 s。
4）功能验证：画出 gray_code 的 stem 图像，如下图所示，幅度为格雷码对应的数值。经验证，函数完成了设计功能。

4.3信噪比转方差子函数(SNR2var)

% 给定信噪比，噪声方差计算子函数
function [var]=SNR2var(SNR,Eb)
    N0=Eb/(10^(SNR/10));
    var=N0/2;
    disp("噪声方差为:")
    disp(var)
end

1）编程思路：根据公式 ${\sigma }^2=\frac{1}{2}\frac{{\mathrm{E}}_{\mathrm{b}}}{10^{SNR/10}}$ 将输入的信噪比(dB)转为方差。
2）输入：信噪比 SNR、每比特能量 $E_b$。
3）输出：噪声方差 var。
4）功能验证：输入信噪比，可在控制台观察到输出，如下图所示。经验证，函数完成了设计功能。

4.4接收信号生成子函数(ReciveSignal)

% 接收信号生成
function [rn]=ReciveSignal(s,var)%s是信道输入信号 %var是噪声的方差
    Len=length(s(1,:));% 根据格雷码映射数组s确定所需产生噪声序列长度Len
    % 方法一
    % 将噪声标准差与randn产生的随机数相乘得到高斯白噪声的同相和正交分量
    nc=randn(1,Len)*sqrt(var);
    ns=randn(1,Len)*sqrt(var);
    % 再分别与二维数组s的同相和正交分量相加，得到接收信号rn的同相和正交分量。
    rn=s;
    rn(1,:)=rn(1,:)+nc;
    rn(2,:)=rn(2,:)+ns;
    %方法二
%     noise=normrnd(0,var,[length(s(1,:)),length(s(1,:))]);
%     rn=s+noise;
end

1）编程思路：将噪声标准差与 randn 产生的随机数相乘得到高斯白噪声的同相和正交分量，再分别与二维数组 s 的同相和正交分量相加，得到接收信号 $r_n$ 的同相和正交分量。也可直接使用函数 normrnd 生成方差为 var 的噪声，如方法二。
2）输入：发送序列二维数组 s、噪声方差 var。
3）输出：接收信号 $r_n$。
4）功能验证：画出加噪声前发送信号的余弦分量和加噪声后信号（接收信号）的余弦分量如表 3。经验证，函数完成了设计功能。

4.5最小距离判决子函数(MinDistance)

% 最小距离判决
function [Judge]=MinDistance(rn,M)
    %生成理想的坐标
    ideal=zeros(2,M);% 2*M数组，存理想位置的横纵坐标
    m=0:M-1;
    ideal(1,:)=cos(2*pi*m/M);
    ideal(2,:)=sin(2*pi*m/M);% 理想的格雷码映射后的二维数组ideal
    Len=length(rn(1,:));% 得到需要判决的信号的长度
    Judge=zeros(1,Len);
    % 计算每一个接收符号的坐标分别到理想坐标的距离
    for a=1:Len
        for b=1:M
            x(b)=abs(rn(1,a)-ideal(1,b));%横坐标距离
            y(b)=abs(rn(2,a)-ideal(2,b));%纵坐标距离
            distance(b)=x(b)^2+y(b)^2;
        end
        mindistance=min(distance);
        Judge(a)=find(distance==mindistance)-1;%找到数组中距离最小的索引-1就是对应的格雷码
    end
end

1）编程思路：先生成理想的格雷码坐标，得到理想的格雷码映射后的二维数组 ideal。再根据 rn 得到需要判决的信号的长度。二维数组，第一行表示 x 轴，第二行表示 y 轴，通过 for 循环实现计算每一个接收符号的坐标分别到四个理想坐标的距离 distance （1 行 4 列的数组），再 min 得到最小距离在数组中对应的索引，因为 matlab 中数组索引从 1 开始，所以该索引-1 即为对应的格雷码标号。
2）输入：接收信号 rn、调制方式的载波数目 M。
3）输出：判决信号 Judge。

4.6最大投影点准则判决子函数(MaxProjection)

% 最大投影点准则判决
function [Judge]=MaxProjection(rn,M)
    % 生成理想的格雷码坐标
    ideal=zeros(2,M);% 2*M数组，存理想位置的横纵坐标
    m=0:M-1;
    ideal(1,:)=cos(2*pi*m/M);
    ideal(2,:)=sin(2*pi*m/M);% 理想的格雷码映射后的二维数组ideal
    Len=length(rn(1,:));% 得到需要判决的信号的长度
    Judge=zeros(1,Len);
    % 计算每一个接收符号的坐标分别到理想坐标的矢量积
    for i=1:Len
        jans=ideal'*rn(:,i); %转置再与rn的第i列相乘
        maxc=max(jans);
        Judge(i)=find(jans==maxc)-1;%找到最大向量积的位置再-1，得到对应的格雷码
    end
end

1）编程思路：同最小距离判决子函数。先生成理想的格雷码坐标，得到理想的格雷码映射后的二维数组 ideal。再根据 rn 得到需要判决的信号的长度。二维数组第一行表示 x 轴，第二行表示 y 轴，通过 for 循环实现计算每一个接收符号的坐标分别到四个理想坐标的矢量积，这里通过先将 ideal 转置得到 4*2 的数组再与 rn 的 2 行第 i 列相乘实现。再 max 得到最大矢量积在数组中对应的索引，因为 matlab 中数组索引从 1 开始，所以该索引-1 即为对应的格雷码标号。
2）输入：接收信号 rn、调制方式的载波数目 M。
3）输出：判决信号 Judge。

4.7 反格雷码子函数(Rebuild)

% 反格雷码子函数，重构二进制信号
function [dnew]=Rebuild(Judge,M)
    % 逐个对应恢复
    Len=length(Judge);
    if M==4
        for i=1:Len
            switch Judge(i)
                case 0
                  dnew(2*i-1)=0;
                  dnew(2*i)=0;
                case 1
                  dnew(2*i-1)=0;
                  dnew(2*i)=1;
                case 2
                  dnew(2*i-1)=1;
                  dnew(2*i)=1;
                case 3
                  dnew(2*i-1)=1;
                  dnew(2*i)=0;
            end
        end
    elseif M==8
        for i=1:Len
            switch Judge(i)
                case 0
                  dnew(3*i-2)=0;
                  dnew(3*i-1)=0;
                  dnew(3*i)=0;
                case 1
                  dnew(3*i-2)=0;
                  dnew(3*i-1)=0;
                  dnew(3*i)=1;
                case 2
                  dnew(3*i-2)=0;
                  dnew(3*i-1)=1;
                  dnew(3*i)=1;
                case 3
                  dnew(3*i-2)=0;
                  dnew(3*i-1)=1;
                  dnew(3*i)=0;
                case 4
                  dnew(3*i-2)=1;
                  dnew(3*i-1)=1;
                  dnew(3*i)=0;
                case 5
                  dnew(3*i-2)=1;
                  dnew(3*i-1)=1;
                  dnew(3*i)=1;
                case 6
                  dnew(3*i-2)=1;
                  dnew(3*i-1)=0;
                  dnew(3*i)=1;
                case 7
                  dnew(3*i-2)=1;
                  dnew(3*i-1)=0;
                  dnew(3*i)=0;
            end
        end
    end
end

1）编程思路：由 Judge 得到最小距离/最大相关度量所对应的格雷码的标号和恢复后的符号长度，使用格雷码一一对应的方法由 swith-case 语句实现反格雷码变换，重构得到发送序列。
2）输入：判决信号 Judge、调制方式的载波数目 M。
3）输出：重构的发送序列 dnew。
4）功能验证：以要传输的二进制比特数 L=50，方差 ${\sigma }^2$ = 0 ，进行单指标噪声分析，调制方式为 QPSK ，使用最小距离判决法进行判决为例，得到的发送序列和重构序列如表 4。经验证，函数完成了设计功能。

4.8误比特率和误符号率计算子函数(ErrorRate)

% 误比特率和误符号率计算子函数
function[BER,SER]=ErrorRate(d,dnew,gray_code,Judge)
    %用输入序列d与重建信号dnew得到BER
    %用gray_code与judge可以得到SER
    Len1=length(d); 
    Len2=length(gray_code); % 分别得到输入序列比特数Len1和输入序列的符号数Len2
    count1=0;
    count2=0; % 设置两个计数变量count1，count2
    for i=1:Len1 % for循环实现比较输入序列d和重构序列dnew
        if d(i)~=dnew(i)
            count1=count1+1; % 记录两者不相等的个数（误比特数）
        end
    end
    BER=count1/Len1; % 误比特数除以总比特个数Len1，得到误比特率BER
    % 同理对gray_code和Judge比较可以得到误符号率SER
    for i=1:Len2
        if gray_code(i)~=Judge(i)
            count2=count2+1;
        end
    end
    SER=count2/Len2;
end

1）编程思路：分别设置两个计数变量 count1，count2。通过 d 得到输入序列比特数 Len1，由 gray_code 得到输入序列的符号数 Len2 。for 循环实现比较输入序列 d 和重构序列 dnew，计算两者中不相等的个数除以总比特个数 Len1，得到误比特率 BER。同理 for 循环实现比较 gray_code 和 Judge，计算两者中不相等的个数除以总符号个数 Len2，得到误符号率 SER。
2）输入：发送序列 d、格雷码编码序列 gray_code、判决信号 Judge、重构的发送序列 dnew。
3）输出：误比特率 BER、误符号率 SER。
4）功能验证：输入噪声方差为 0.1，进行单指标噪声分析，调制方式为 QPSK ，使用最小距离判决法进行判决，得到的误比特率和误符号率如下图。经验证，函数完成了设计功能。

4.9Monte Carlo 仿真_求实际误比特率子函数(Curve_BER)

% Curve的子函数求误比特率
function [BER]=Curve_BER(d,s,var,flag,M)
    [rn]=ReciveSignal(s,var);
    if flag==0
        [Judge]=MinDistance(rn,M);
    else
        [Judge]=MaxProjection(rn,M);
    end
    [dnew]=Rebuild(Judge,M);
    Len1=length(d);
    count1=0;
    for i=1:Len1 % for循环实现比较输入序列d和重构序列dnew
        if d(i)~=dnew(i)
            count1=count1+1; % 记录两者不相等的个数（误比特数）
        end
    end
    BER=count1/Len1; % 误比特数除以总比特个数Len1，得到误比特率BER
end

1）编程思路：首先通过调用 ReceiveSignal 子函数获得接收信号 rn ，接着根据输入选用不同的判决方法，得到判决后的序列，接着调用 Rebuild 子函数得到重构序列 dnew。实际误比特率的计算同 ErrorRate 子函数，for 循环实现比较输入序列 d 和重构序列 dnew，计算两者中不相等的个数除以总比特个数 Len1，得到误比特率 BER。
2）输入：发送序列 d、噪声方差 var、判决方式选择标志变量 flag、发送序列二维数组 s、调制方式的载波数目 M。
3）输出：实际误比特率 BER。

4.10Monte Carlo 仿真_求实际误符号率子函数(Curve_SER)

% Curve的子函数求误码率
function [SER]=Curve_SER(s,gray_code,var,flag,M)
    [rn]=ReciveSignal(s,var);
    %用gray_code与judge可以得到SER
    if flag==0
        [Judge]=MinDistance(rn,M);
    else
        [Judge]=MaxProjection(rn,M);
    end
    Len2=length(gray_code); % 分别得到输入序列比特数Len1和输入序列的符号数Len2
    count2=0; % 设置两个计数变量count1，count2
     for i=1:Len2
        if gray_code(i)~=Judge(i)
            count2=count2+1;
        end
    end
    SER=count2/Len2;
end

1）编程思路：首先通过调用 ReceiveSignal 子函数获得接收信号 rn ，接着根据输入选用不同的判决方法，得到判决后的序列。实际误符号率的计算同 ErrorRate子函数，for 循环实现比较 gray_code 和 Judge，计算两者中不相等的个数除以总符号个数 Len2，得到误符号率 SER。
2）输入：格雷码编码序列 gray_code、噪声方差 var、判决方式选择标志变量 flag、发送序列二维数组 s、调制方式的载波数目 M。
3）输出：实际误符号率 SER。

4.11Monte Carlo 仿真子函数(Curve)

% 产生发送序列
% 仿真误比特率曲线和理论误比特率曲线绘制主函数
function []=Curve(N,d,M,flag)
    if M==4
        % QPSK仿真误比特率曲线和理论误比特率曲线
        [~,s]=Grayreflect(N,d,M);
        SNR=-5:1:11; % 设置SNR的范围为-5~11dB
        Eb=1/2; % 已知QPSK的比特能量为1/2
        Len=length(SNR);
        Theory2=zeros(1,Len);
        Real2=zeros(1,Len);
        for i=1:Len
	        N0=Eb/(10^(SNR(i)/10));
            var=N0/2; % 由SNR得var
            Theory2(i)=qfunc(sqrt(2*Eb/N0)); % 理论BER的计算
            Real2(i)=Curve_BER(d,s,var,flag,M);%调用得到实际BER 
        end
        figure;
        %用semilogy画出理论误比特率曲线 以10为基数的对数刻度与x轴的线性刻度
        Theorymap=semilogy(SNR,Theory2);
        %理论误比特率曲线用红色圆点表示
        Theorymap.Color='r';
        Theorymap.Marker='o';
        hold on;
        %实际误比特率曲线用蓝色星型表示
        Realmap=semilogy(SNR,Real2);
        Realmap.Color='b';
        Realmap.Marker='*';
        grid on;
        title('QPSK误比特率曲线');
        xlabel('信噪比SNR(Eb/N0)/dB');ylabel('误比特率');
        legend('QPSK误比特率理论值','QPSK误比特率实际值','Location','southwest');
    elseif M==8
        % 8PSK仿真误码率曲线和理论误码率曲线
        [gray_code,s]=Grayreflect(N,d,M);
        SNR=-5:1:11; % 设置SNR的范围为-5~11dB
        Eb=1/3;
        Len=length(SNR);
        Theory2=zeros(1,Len);
        Real2=zeros(1,Len);
        for i=1:Len
	        N0=Eb/(10^(SNR(i)/10));
            var=N0/2; % 由SNR得var
            Theory2(i)=2*qfunc(sqrt(6*Eb/N0)*sin(pi/8)); % 理论SER的计算
            Real2(i)=Curve_SER(s,gray_code,var,flag,M);%调用得到实际SER 
        end
        figure;
        %用semilogy画出理论误码率曲线 以10为基数的对数刻度与x轴的线性刻度
        Theorymap=semilogy(SNR,Theory2);
        Theorymap.LineStyle=':';
        Theorymap.Color='k';
        Theorymap.Marker='o';
        hold on;
        Realmap=semilogy(SNR,Real2);
        Realmap.LineStyle=':';
        Realmap.Color='m';
        Realmap.Marker='*';
        grid on;
        title('8PSK误码率曲线');
        xlabel('信噪比SNR(Eb/N0)/dB');ylabel('误码率');
        legend('8PSK误码率理论值','8PSK误码实际值','Location','southwest');
    end
end    

1）编程思路：分两种情况进行讨论。M=4 时绘制 QPSK 调制仿真误比特率曲线和理论误比特率曲线。首先设置 SNR 的范围为-5~11dB，接着调用信噪比转方差子函数将相应的信噪比转换为方差。接下来使用公式计算得到 QPSK 的理论误比特率；调用子函数计算得到仿真误比特率。 M=8 时绘制 8PSK 调制仿真误符号率曲线和理论误符号率曲线。同理，首先设置 SNR 的范围为-5~11dB，接着调用信噪比转方差子函数将相应的信噪比转换为方差。接下来使用公式计算得到 8PSK 的理论误符号率；调用子函数计算得到仿真符号特率。
2）输入：发送序列 d、调制方式的载波数目 M、判决方式选择标志变量 flag、发送序列长度 N。
3）输出：仿真误比特率曲线和理论误比特率曲线。
4）功能验证：设置要传输的二进制比特数 L=1000，进行 Monte Carlo 仿真，调制方式为 QPSK ，使用最小距离判决法进行判决，得到的仿真误比特率曲线和理论误比特率曲线如下图。经验证，函数完成了设计功能。

4.12Monte Carlo 仿真子函数_QPSK 和 8PSK 对比(CurveAll)

function CurveAll(L)
    M=4;
    flag=0;
    [N,d]=GenerateSendSequence(L,M);
    [gray_code,s]=Grayreflect(N,d,M);
    SNR=-5:1:11; % 设置SNR的范围为-5~11dB
    Eb=1/2; % 已知QPSK的比特能量为1/2
    Len=length(SNR);
    Theory=zeros(1,Len);
    Real=zeros(1,Len);
    for i=1:Len
	    N0=Eb/(10^(SNR(i)/10));
        var=N0/2; % 由SNR得var
        Theory(i)=2*qfunc(sqrt(2*Eb/N0))*(1-0.5*qfunc(sqrt(2*Eb/N0))); % 理论SER的计算
        Real(i)=Curve_SER(s,gray_code,var,flag,M);%调用得到实际SER 
    end
    
    M=8;
    [N,d]=GenerateSendSequence(L,M);
    [gray_code,s]=Grayreflect(N,d,M);
    SNR=-5:1:11; 
    Eb=1/3;
    Len=length(SNR);
    Theory2=zeros(1,Len);
    Real2=zeros(1,Len);
    for i=1:Len
        N0=Eb/(10^(SNR(i)/10));
        var=N0/2; % 由SNR得var
        Theory2(i)=2*qfunc(sqrt(6*Eb/N0)*sin(pi/8)); % 理论SER的计算
        Real2(i)=Curve_SER(s,gray_code,var,flag,M);%调用得到实际SER 
    end
    x=-5:1:11;
    figure;
    Theorymap=semilogy(x,Theory);
    Theorymap.Color='r';
    Theorymap.Marker='o';
    hold on;
    Realmap=semilogy(x,Real);
    Realmap.Color='b';
    Realmap.Marker='*';
    Theorymap=semilogy(x,Theory2);
    Theorymap.Color='k';
    Theorymap.Marker='o';
    Realmap=semilogy(x,Real2);
    Realmap.Color='m';
    Realmap.Marker='*';
    grid on;
    title('QPSK和8PSK误码率曲线');
    xlabel('信噪比SNR(Eb/N0)/dB');ylabel('误码率');
    legend('QPSK误码率理论值','QPSK误码率实际值','8PSK误码率理论值','8PSK误码实际值','Location','southwest');
    hold off
end

1）编程思路：同时绘制 QPSK 和 8PSK 调制仿真误符号率曲线和理论误符号率曲线，进行对比分析。方法与上面的子函数相同。
2）输入：要传输的数据点数 L。
3）输出：仿真误比特率曲线和理论误比特率曲线。
4）功能验证：设置要传输的二进制比特数 L=10000，进行 Monte Carlo 仿真， 进行两种调制方式的对比。得到的仿真误比特率曲线和理论误比特率曲线如下图所示。经验证，函数完成了设计功能。

---

五、实验内容

5.1AWGN 信道下QPSK 调制通信系统设计与性能研究

加粗内容为实验要求

1)最大投影点准则进行判决
a.画出噪声方差分别为 0、0.1、0.5、1.0 时，在检测器输入端 1000 个接收到的信号加噪声的样本（星座图）；

【方法】
从键盘输入数据点数 1000 点，选择进行单指标系统仿真以及 QPSK 调制，分别输入噪声方差 0、0.1、0.5、1.0，再选则最大投影点准则判决，即可得到仿真运行结果。程序给出了接收信号的样本（星座图），以及误比特率、误符号率仿真计算结果。

【结果与分析】
运行结果如表 5 所示。
结论：在信号点数为 1000 时，随着噪声方差从 0 到 1.0 不断增大，可以观
察到星座图上从噪声方差为 0 时分别集中在单位圆与坐标轴相交的四个点上（理想信号点），到随着噪声方差的增大，分布越来越离散。同时误比特率和误码率也随着噪声方差地增大不断增大。说明噪声的方差越大（功率越大），对接收信号的影响越大。随着噪声方差的增大，信号无法被正确判决，重构信号的错误概率越来越大。

b.在 AWGN 信道下，分别画出数据点（二进制比特数）为 1000、10000、 100000 时的 Monte Carlo 仿真误比特率曲线和理论误比特率曲线，比较差别，分析数据点的数量对仿真结果的影响（横坐标是 snr=Eb/No dB，格雷码映射）

【方法】
已知 QPSK 的比特能量为 1/2,设置 SNR 的范围为 -5~11dB,输入数据点数， 先构造误比特率理论和实际的空的长为 SNR 长度的序列，再根据公式 $N_0=$ $\frac{E_b}{10^{SNR/10}}$、${\sigma }^2=\frac{N_0}{2}$ 、$P_{QPSK}=Q(\sqrt{\frac{2{\epsilon }_b}{N_0}})$，由 for 循环实现理论误比特率曲线的 y 轴值，仿真误比特率曲线的 y 值通过调用误码率获取子函数得到，最后通过semilogy 语句和 hold on 实现仿真误比特率曲线和理论误比特率曲线在一张图上用不同颜色和不同点样式输出。

【结果与分析】
运行结果如表 6 所示。
结论：随着点数的增大，误比特率的实际值与理论值曲线越来越相近。对于信噪比与误比特率的关系，都是随着信噪比的增大，误比特率不断减小。

2）将检测器的判决准则改为最小距离法（星座图上符号间的距离）
A.画出数据点（二进制比特数）为 100000 时的 Monte Carlo 仿真误比特率曲线和理论误比特率曲线 （横坐标是 snr=Eb/No dB，格雷码映射，在 AWGN 信道下）；
B.比较与最大投影点准则进行判决结果的区别。
C.并理论分析证明两种判决方法是等价的。
D.理论证明 QPSK 误比特率与 BPSK 的一样。

【方法】
同最大投影判决法，区别是在选择判决方法时的输入。

【结果与分析】
运行结果如表 7 所示。
当数据点数为 100000 时，由最大投影准则和最小欧式距离准则时的 Monte Carlo 仿真误比特率曲线和理论误比特率曲线对比可得，在相同信噪比 SNR 范围下，二者的实际值和理论值曲线都几乎吻合。最大投影点准则与最小欧式距离准则等价。

理论分析：
根据 ML 最大似然准则，在 M 个信号上找寻最大的 $f(r|s_{\mathfrak{m}})$，通过求噪声的 PDF 后再求联合 PDF，可得：
$$\mathrm{f}(\mathbf{r}|s_m)=(\pi N_0{)}^{-N/2}{e}^{\frac{1}{N}\Vert \mathbf{r}-{\mathbf{s}}_{\mathbf{m}}{\Vert }^2}m=1,2,...,M$$

为了简化计算取 $f(r|s_{\mathfrak{m}})$ 的自然对数，可见它是一个单调函数，因此有：
$$\mathrm{lnf}(r|{\mathrm{s}}_{\mathrm{m}})=-\frac{\mathrm{N}}{2}\ln (\pi {\mathrm{N}}_0)-\frac{1}{{\mathrm{N}}_0}\sum _{\mathrm{k}=1}^{\mathrm{N}}({\mathrm{r}}_{\mathrm{k}}-{\mathrm{s}}_{\mathrm{mk}}{)}^2$$

寻找使 $\mathrm{lnf}(r|{\mathrm{s}}_{\mathrm{m}})$ 最大的信号 $s_{\mathfrak{m}}$，等价于寻找取得最小欧式距离的信号：
$$\mathbf{D}(\mathbf{r},{\mathbf{s}}_{\mathbf{m}})=\sum _{\mathbf{k}=1}^{\mathbf{N}}({\mathbf{r}}_{\mathbf{k}}-{\mathbf{s}}_{\mathbf{mk}}{)}^2$$

通过展开上式距离变量而得到：
$$\mathbf{D}(\mathbf{r},{\mathbf{s}}_{\mathbf{m}})=\Vert r{\Vert }^2-2r\cdot {\mathbf{s}}_{\mathbf{m}}+\Vert {\mathbf{s}}_{\mathbf{m}}{\Vert }^2$$

注意到选择使 $D^{\prime }(\mathbf{r},{\mathbf{s}}_{\mathbf{m}})$ 最小的信号 $s_{\mathrm{m}}$，等价于选择使变量 $C(\mathbf{r},{\mathbf{s}}_{\mathbf{m}})$最大的信号。

$$\mathbf{C}(\mathbf{r},{\mathbf{s}}_{\mathbf{m}})=2r\cdot {\mathbf{s}}_{\mathbf{m}}-\Vert {\mathbf{s}}_{\mathbf{m}}{\Vert }^2$$
其中 $\mathbf{r}\cdot {\mathbf{s}}_{\mathbf{m}}$ 项表示接收信号矢量 $\mathbf{r}$ 在 M 个可能发送信号矢量上的投影， 它是接收信号矢量与第 m 个发送信号之间相关性的度量。所以 $\mathbf{C}(\mathbf{r},{\mathbf{s}}_{\mathbf{m}})$ 叫相关度量。
综上可知最大投影准则与最小欧式距离准则是等价的。

理论证明 QPSK 误比特率与 BPSK 的一样：
已知 MPSK 的表达式为 $u_m$( t) = $\sqrt{\frac{2{\epsilon }_s}{T}}$cos( 2$\pi f_c$t+ $\frac{2\pi m}{M}$) , m= 0, 1, $\dots \dots$M- 1，因此 M=2 时，BPSK 表达式为：
$${\mathrm{u}}_{\mathrm{BPSK}}(\mathrm{t})=\{\begin{array}{l}\frac{\sqrt{2{\epsilon }_s}}{T}cos2\pi f_{c}t,m=0\\ -\frac{\sqrt{2{\epsilon }_s}}{T}cos2\pi f_{c}t,m=1\end{array}$$

M=4 时，QPSK 的表达式为：
$${\mathrm{u}}_{\mathrm{QPSK}}(\mathrm{t})=\pm \mathrm{A}{\mathrm{g}}_{\mathrm{T}}(\mathrm{t})\cos (2\pi {\mathrm{f}}_{\mathrm{c}}\mathrm{t})-{\mathrm{g}}_{\mathrm{T}}(\mathrm{t})\pm \mathrm{Asin}(2\pi {\mathrm{f}}_{\mathrm{c}}\mathrm{t})$$

可以看成两个正交的二相调制之差。当载波相位估计准确时，两对正交载波调制信号间不会发生串扰现象，所以 OPSK 系统的比特错误概率等于 BPSK 的错误概率为 $P_{bpsk}=Q\left[\sqrt{\frac{2{\epsilon }_b}{N_0}}\right]$。

5.2AWGN 信道下 8PSK 调制通信系统设计与性能研究

1）画出噪声方差分别为 0、0.1、0.5、1.0 时， 在检测器输入端 1000 个接收到的信号加噪声的样本（星座图）；

【方法】
从键盘输入数据点数 1000 点，选择进行单指标系统仿真以及 8PSK 调制，分别输入噪声方差 0、0.1、0.5、1.0，再选则最大投影点准则判决，即可得到仿真运行结果。程序给出了接收信号的样本（星座图），以及误比特率、误符号率仿真计算结果。

【结果与分析】
运行结果如表 8 所示。
结论：在信号点数为 1000 时，随着噪声方差从 0 到 1.0 不断增大，可以观察到星座图上信号点随着噪声方差的增大分布越来越离散。同时误比特率和误码率也随着噪声方差地增大不断增大。说明噪声的方差越大（功率越大），对接收信号的影响越大。随着噪声方差的增大，信号无法被正确判决，重构信号的错误概率越来越大。

2）检测器的判决准则选为最小距离法，格雷码映射，比较数据点为 100000 时8PSK 与 QPSK 的 Monte Carlo 仿真误符号率曲线、理论误符号率曲线，比较差别。 （横坐标是 snr=Eb/No）。（一张图上呈现 4 条曲线）

【方法】
输入数据点数 100000 点，选择进行 Monte Carlo 仿真中的 QPSK 和 8PSK的对比功能，分别将 QPSK 调制的仿真误符号滤和理论误符号率和 8PSK 的仿真误符号滤和理论误符号率绘制到同一张图上。

【结果与分析】
运行结果如下图所示。
可以看出，在相同的信噪比下，QPSK 调制的误符号率总是优于 8PSK 调制。也就是说，相对于 8PSK 系统，QPSK 系统的性能要更好一些。

3)理论分析 8PSK 性能为什么比 QPSK 差。
由于 8PSK 使用了更多的相位，它对噪声和失真更加敏感，从而导致性能相对较差。噪声会引起相位偏移和幅度衰减，而 8PSK 对相位的敏感性更高。当信道噪声较大时，8PSK 信号容易受到干扰和误差的影响，导致传输错误率增加。而且在 8PSK 中，不同相位之间的距离更近，使得相邻相位之间发生干扰的可能性增加。这种码间干扰会降低信号的可靠性，导致误比特率增加。
对于两种调制方式的误符号率，由实验原理部分的分析，M>4 时，错误概率的近似值如下式所示：$P_{MPSK}=2Q(\sqrt{\frac{2(lo{g}_{2}M){\epsilon }_b}{N_0}}sin\frac{\pi }{M})$。可得 8PSK 调制的误符号率为 $P_{8PSK}=2Q(\sqrt{\frac{6{\epsilon }_b}{N_0}}sin\frac{\pi }{8})$ ，而 QPSK 的误符号率为$P_{Qpsk}=2Q(\sqrt{\frac{2{\epsilon }_b}{N_0}})\left[1-\frac{1}{2}Q(\sqrt{\frac{2{\epsilon }_b}{N_0}})\right]$ 。由理论计算可以得到，对于相同的 $N_0$，QPSK 的误符号率更小一些。

---

正文完，实验总结部分在文件：山东大学通信原理实验软件部分-实验报告+代码-MATLAB-数字基带传输系统设计与性能研究-MPSK 通信系统的设计与性能研究`，包含问题与解决，实验心得等部分。

如有错误欢迎指正。
祝好！