一种航天测控链路复合调制信号的自动调制识别方法

本发明属于信号处理领域，特别是复合调制信号的信号调制制式识别设计。

背景技术：

1、复合调制(composite modulation,cm)指多用户基带数据采用不同的调制方式调制到统一载波上的多层调制方案，是提高航天测控通信系统传输效率、抗干扰能力和安全性的关键技术。复合调制识别(automatic composite modulation classification,acmc)在以认知无线电为基础的航天遥测、跟踪和控制(telemetry,tracking and command,tt&c)及通信系统中发挥着重要作用，如图1所示。目前航天通信对认知无线电的需求迅速增长，为实现下一代测控通信系统的频谱共享、自适应调制和灵活传输需求，必须对复合调制信号的调制类型进行准确识别。然而，复合调制方案引入了多层调制方法之间的非线性耦合，加上空间传输信道引入的信道噪声和大尺度频偏等损失因素，使得接收端在无先验知识的条件下对接收端复合调制制式的识别成为了巨大的挑战。因此，研究新的鲁棒高效的复合调制信号调制分类方法受到了广泛关注。

2、目前复合调制信号的调制分类方法研究有限，且因复合调制方案的复杂性和非线性，传统单调制信号的调制识别方法在引入到复合调制时面临困难。传统的单调制信号的调制识别方法大体上可以分为两类：基于似然函数的调制识别方法和基于特征提取的调制识别方法，因计算复杂度和特征稳定性限制此两种方法都不能取得令人满意的识别效果。

3、实际上，关于复合调制的识别问题罕有研究。传统的复合调制识别流程为“逐层识别并解调”，在此过程中，接收端需要首先使用现有方法识别复合调制信号的外层调制类型，进一步解调后识别内层。这种识别方法受诸多条件限制，如相位误差、定时误差、频率偏移等，外层的识别和解调在一定程度上会影响内层信号的识别精度。幅度谱和功率谱特征作为接收端信号变换域特征可被提取作为一种联合特征进行调制识别，但该种方法并非全无先验知识，而是在信号外层预解调的假设下完成的调制识别，使用场景受限。基于锁相环(phase lock loop,pll)及拟合优度的特征提取方法亦被提出应用于复合调制信号的识别，但当锁相环失锁时该方法失去有效性。高阶累积量(high order statistics,hos)作为信号统计特征可提取应用至信号盲识别，但其庞大的计算量对内存资源消耗很大，时效性不强。而且，毫米波雷达场景下的复合调制识别研究给了航天测控链路场景下的复合调制识别以启发，但因使用的复合调制体制有很大差异，特征改变，不可直接应用。因此探寻航天测控链路的复合调制识别方案的需求呼之欲出。

4、此外，随着深度学习方法及网络的蓬勃发展，通信信号的特征可以通过网络自动学习并输出分类结果。卷积神经网络(convolution nerual networks,cnn)、反馈神经网络(recurrent nerual networks,rnn)和对抗生成网络(generative adversarialnetworks,gan)及其变体，都被引入到调制识别当中来。然而，复合调制信号的特征难以提取，到目前不存在任何针对复合调制信号直接设计的深度学习网络。如何在不增加运算复杂度和特征提取难度的基础上开创一套崭新的高精度复合调制识别架构将是本发明关注并解决的一个重点问题。

5、本发明目前还未在国内外出版物上公开发表文章，未在国内外公开使用或者以其他方式为公众所知。

技术实现思路

1、本发明的目的在于在克服已有复合调制识别方法的不足，构造一套完整的针对航天测控统一载波体制复合调制信号的识别流程。该套估计流程的优势在于精度高，抗干扰性能好，可有效对抗航天链路引入的相位噪声、频率偏移和定时误差，且所需先验知识少，计算量和复杂度低。

2、本发明基于复合调制信号的二阶循环谱分析，提出了一种新的可提取特征量：循环爪印矩阵(cpp)，该矩阵可以唯一的表征不同种类的复合调制信号；考虑到该矩阵的维度较大，本发明引入了结合dct和dwt的矩阵降维方法将矩阵特征降维成一个1*18的特征向量，并引入序贯最小优化算法(sequential minimal optimization,smo)优化的支持向量机作为分类器输出最终的识别结果，整体识别方案如图2所示。

3、本发明技术方案为一种航天测控链路复合调制信号的自动调制识别方法，该方法包括：

4、步骤1：构造存储在接收端的不同信噪比下理想复合调制信号特征训练集

5、步骤1.1：生成理想复合调制信号集；

6、该信号集中包含十种统一载波体制下的复合调制信号，分别为pcm/bpsk/pm、pcm/qpsk/pm、pcm/bpsk1+bpsk2/pm、pcm/qpsk1+qpsk2/pm、pcm/bpsk+qpsk/pm、pcm/bpsk/fm、pcm/qpsk/fm、pcm/bpsk1+bpsk2/fm、pcm/qpsk1+qpsk2/fm、pcm/bpsk+qpsk/fm这些信号的调制指数kp＝1.2，理想复合调制信号无噪干扰、不考虑相位噪声和多普勒频移；这里给出理想复合调制信号的发射信号s(t)的模型：

7、

8、其中，是等效低通信号，fc是主载波频率，φ0是主载波的初始相位，代表取实部，

9、为：

10、

11、其中，a代表调制信号幅度；kpm,kfm分别是两种复合调制方式的调制指数，si(t)是内层调制信号，具体又分为bpsk、qpsk及多用户组合方式，其中表示内层调制信号的个数，即统一载波体制下的用户个数。

12、内层调制的信号模型如下：

13、

14、其中，

15、其中，是内层剩余载波调制的载波频率，φi,0是内层剩余载波调制的初始相位，i指的是内层调制信号的数目，i取到gi(t)代表矩形脉冲，ti代表第i路成型脉冲的持续时间，其中k代表的的是第k个用户信息码元；

16、步骤1.2：计算不同信噪比下理想复合调制信号的三维归一化循环谱

17、对不同信噪比下的理想离散复合调制信号采用二阶循环理论中的快速傅里叶变换累积方法求取二阶循环自相关函数；给定频谱频率f和循环频率ε，基于时域平滑的循环周期图表示为：

18、

19、其中，g′(t)是时域宽度为δt＝nts的归一化权重函数，n为每个傅里叶变换窗的采样点数，ts是时域采样周期；f1，f2定义了在fam方法中使用的滤波器中心频率，其中

20、f1＝f+ε/2，f2＝f+ε/2；rt(λ，f1),rt(λ,f2)是离散复合调制信号的复解调，这两个量由如下的公式得到：

21、

22、其中，ω(λ)表示一个时域采样长度为t＝n′ts的矩形窗函数，该函数的带宽和公式(5)二阶循环谱的频率分辨率δf保持一致；因此，循环自相关函数可由时域平滑循环周期量实现无偏估计，进一步表示为：

23、

24、由此求得的二维循环谱为具有非负幅值的三维频谱，该谱包含2n+1个循环频率，ε＝εp,p＝-n,-n+1,...,n和n+1个频谱频率f＝fq,q＝-n′/2,-n′/2+1,...,n′/2；更进一步，我们对得到的二维循环谱做归一化处理得到归一化循环谱表示为：

25、

26、步骤1.3：对得到的三维归一化循环谱做循环爪印矩阵映射；

27、由步骤1.2得到的三维归一化循环谱可转化成其俯视图；如图3所示。明显可察复合调制信号的循环谱俯视图集中在ε-f平面上对称分布的四块区域，其中每块区域的视图就像动物的爪印，故将一个由三维归一化循环谱得到的俯视图称为循环爪印图(cyclic paw print,cpp)。不同的复合调制信号得到的循环爪印图区别明显，因此该图可以作为唯一标识某种复合调制类型的特征图像。本发明将三维循环谱的俯视图用一个二维p-q矩阵表示，即有：

28、

29、中的元素是给定循环频率ε和频谱频率f时的归一化非负幅度值；进一步将中的元素采用16比特量化，由此将图片转化成16比特灰度矩阵，为：

30、

31、其中表示向下取整。因此根据以上公式，可以将复合调制信号的归一化三维循环谱转化成p-q的灰度矩阵，即cpp矩阵

32、步骤1.4：利用离散余弦变换(dct)和二维离散小波变换(dwt)对cpp特征矩阵降维，构建融合特征向量(hybrid feature vector,hfv)；

33、步骤2：对接收端实际接收到的某种或多种复合调制信号按照上述方法求取归一化循环谱(nscs)，得到cpp矩阵；

34、步骤3：对得到的cpp矩阵并行做dct和dwt+svd的处理，构建复合调制融合特征分量hfv，得到实际接收信号的测试集；

35、步骤4：采用序贯最小优化算法(sequential minimal opmization,smo)优化的支持向量机(support vector machine,svm)进行分类输出。

36、步骤4.1：接收端首先使用不同信噪比下理想复合调制信号的训练集训练支持向量机，引入smo算法后，主要针对的是支持向量机的训练过程进行优化；非线性svm分类器的决策函数为：

37、对应第k种复合调制

38、

39、不对应第k种复合调制

40、其中表示的最优非线性投影函数，为最优的偏差，二者均有待确定；

41、将k类复合调制识别的多分类问题简化为k个并行独立二分类的子问题，架构如图6所示；对于给定的融合特征向量有如下定义：

42、

43、传统的svm优化问题由公式(21)给出；

44、

45、

46、其中，表示第k个svm的k×1的拉格朗日乘子法向量，c说明了最小化错误分类样本数目的正惩罚项，k()表示一个非线性核函数；为了减少求解公式(21)的计算复杂度，引入smo启发式算法优化svm的目标函数，该算法将svm解决的qr问题分解成多个子问题，每个子问题包含两个独立的拉格朗日因子，其中至少有一个因子不满足kkt条件；设这两个拉格朗日因子为和其他为αk,k',k'＝ι3,ι4,...,ιk,并保持不变；由此，smo-svm分类方式中每个子问题涉及的目标函数表示为：

47、

48、其中：

49、

50、定义且为一个预设置的常数；在每次迭代中，和最初由最初的迭代得到的可行解分别由和表示，定义对应的优化解为分别为和公式(22)中给出了两个变量之间的限制关系，则可求得：

51、

52、公式(23)中的目标函数可以简化为单变量函数，进一步关于对其求微分令其等于0，可以得到：

53、

54、

55、

56、其中，定义为：

57、

58、

59、bk,k'表示对于第k个smo-svm决策面的第k'个拉格朗日乘子式的偏移量；

60、每次迭代优化解可以更新为：

61、

62、

63、其中，h和l分别表示公式(30)约束的上下界，对于k'＝1,2,...,k,有

64、

65、和

66、

67、两个拉格朗日因子确定好后，优化后的偏移量和可以通过计算得到：

68、

69、其中，表示初始迭代中得到的两个拉格朗日乘子式和的偏差；故当前更新后的优化决策面偏差为两偏移量均值：

70、

71、故经smo算法优化后的二分类决策面方程为：

72、

73、其中，表示对非线性决策面函数的估计值；在训练阶段，子问题中的拉格朗日因子至少一个需要违背公式(21)对应的kkt条件，提供优化空间；和的最优值可以由公式(30)和公式(31)确定；因此，通过反复迭代求解包含不同拉格朗因子对的子问题直至收敛，可以获得第k个smo-svm的最优拉格朗日因子向量和决策面最优偏移量在并行处理下，所有k个smo-svm的子分类器都能被训练，并且可以根据公式(35)和(36)获得最优非线性投影函数集和最优偏移量集

74、步骤4.2：接收端接收到航天测控链路中传输的复合调制信号，按照上述步骤1的方法构造融合特征向量hfv，进一步构造测试集将其作为接收端已训练好的多决策支持向量机的输入；对于第k个smo-svm，局部决策度量dk为：

75、

76、此时，实际接收信号的复合调制类型可以被识别为：

77、

78、进一步的，所述步骤1.4的具体方法为：

79、步骤1.4.1：提取二维dct变换后的系数矩阵低频部分作为hfv的第一部分；首先对cpp矩阵s(p,q)进行二维离散余弦变换，得到dct系数矩阵：

80、

81、其中，u＝0,1,...,p-1,v＝0,1,...,q-1，其中βu和βv表示两个相应的乘法因子，可以表示为：

82、

83、

84、经过二维dct变换得到的系数矩阵可以分成两个部分：低频部分和高频部分；其中低频系数主要集中在dct系数矩阵的左上角，幅值上明显大于其他分量，保留了更多的特征信息；因此选取系数矩阵的部分低频系数作为复合调制识别的特征；

85、步骤1.4.2：提取二维dwt变换后的两个维度上的ll子带，经过奇异值分解(singular value decomposition,svd)求取两个子带的奇异值作为hfv的第二部分；首先给出n×n的离散哈尔小波变换操作的公式：

86、

87、其中，为偶数，表示的单位矩阵；的子矩阵和分别表示离散小波变换的低通和高通算子；当对cpp矩阵做二维离散哈尔小波变换时，cpp矩阵的行和列将被分别处理，离散哈尔小波变换后得到的dwt矩阵b可以定义如下：

88、

89、cpp矩阵的dwt矩阵b能够被分解为四个子矩阵；位于b左上角的子矩阵c表示cpp矩阵的近似值，由先后对的列与行取均值得到；其他的三个子矩阵e、f和y是稀疏矩阵，分别阐明了的垂直、水平和对角差异。因为cpp矩阵的能量大部分聚集在c，故保留c作为的dwt简要表示；为提取不同维度上的dwt特征，采用递归的方法对矩阵c进行处理，表示为：

90、

91、其中，c1＝c；为“分解层级指数”，更进一步采用奇异值分解的方法对不同层级上的矩阵c进行矩阵降维，对应得到的奇异值将作为融合特征的第二部分；当分解层级为时，对矩阵cl奇异值分解的表达式为：

92、

93、其中，表示对角元素为的对角矩阵；vl＝min{p/2l,q/2l}；和是酉矩阵，矩阵的列分别为和的特征向量；σl表示包括了cl所有奇异值的对角矩阵，即σκ,κ＝1,2,...,vl。如果cl的秩为那么并且并且在本发明中，正奇异值被用于构建cm信号的融合特征向量；

94、步骤1.4.3：构建复合调制融合特征向量(hfv)训练集和测试集

95、本发明的融合特征向量由步骤1.4.1的dct系数矩阵d的低频部分和步骤1.4.2的dwt+svd方法得到的正奇异值组成；为了保证两部分分量对融合特征向量的影响相同，构建融合特征向量时须保持格式上的一致；经过dwt+svd处理得到的两组奇异值为1*6的向量，故考虑把dct低频系数部分也提取6个元素。以d中的行和列上的平方和作为系数幅值分布指标，其分布分别如图4和图5所示，可以看出位于d的第一行的前几个元素相较于其他的元素有明显更大的幅值；故提取其前六个元素作为特征组成部分构建复合调制信号的融合特征向量即有：

96、

97、其中，d(p,q)、σ1(p,q)和σ2(p,q)分别为d、σ1、σ2的第(p,q)个元素。

98、在有噪声的条件下，对于复合调制候选集中的第k种复合调制方式，即仍然可以采用上述的方法经过ρ次试验从相应的nscs提取cpp，构造器相应的融合特征向量集合，即有其中，表示第k个cm信号生成的混合特征向量k＝1,2,...,k；μ＝1,2,...,ρ；此外，复合调制测试集也能通过相同的方法构建，得到的训练集为

99、本发明相比于其他现有的航空测控链路中的复合调制识别方法，所需的先验知识少，精度高，抗干扰性能好，可有效对抗航天链路引入的相位噪声、频率偏移和定时误差，并有效降低了计算复杂度和内存量，展示出优越的性能。