基于时频能量谱纹理特征的跳频调制方式识别

李红光，郭英，眭萍，齐子森

（空军工程大学信息与导航学院，陕西 西安 710077）

1 引言

跳频（FH,frequency hopping）信号调制方式作为FH 通信的一维重要参数，能够为FH 组网通信的网台分选、情报分析和干扰引导等任务提供重要依据。现有传统特征提取识别方法通过提取时频能量、瞬时波形、星座图和高阶统计量等特征实现调制方式的分类识别。文献[1]通过提取信号8 阶以上矩的甚高阶统计量特征，实现了幅相调制信号的调制方式识别，该算法识别信号种类多，计算特征维数较高，运算量大。文献[2]通过计算分数低阶循环谱和频率谱截面的5 个相关系数组成特征集，利用判决树分类器实现Alpha 稳定分布噪声下的数字调制方式识别，该算法能有效抑制非高斯噪声影响，计算复杂度较低。文献[3]利用线性调制信号的星座图对称性，通过改进KS（Kolmogorov-Smirnov）分类器，实现调制方式的快速识别，提高了对噪声不确定性的稳健性。文献[4]通过构造4 阶和6 阶高阶累积量和循环谱特征，利用神经网络分类决策算法，实现5 类信号的调制方式识别，该算法在低信噪比条件下，识别成功率较高，但可适用的信号调制方式种类较少。文献[5]提出了基于顺序统计量和降阶统计量特征的自动调制分类器，利用反向传播神经网络分类算法实现在已知或未知信道场景中的调制方式分类。文献[6]提取不同阶数和共轭项数的多个循环累积量平稳特征组成特征向量，实现了部分数字调制信号的盲识别，但所需观测样本较多，计算量较大。文献[7]利用信号高阶循环累积量序列在傅里叶域的稀疏性，提取循环累积量稀疏特征，在较少观测样本条件下，实现信号的调制方式识别。该算法计算量小，但当信噪比较低时，识别率下降明显。文献[8]提出了一种基于长短期记忆的调制方式自动分类模型，该模型从训练数据的时域幅度和相位信息中学习提取特征，对不同调制符号速率的信号具有较高的特征表示能力。文献[9]提出了一种基于抗噪声处理和深度稀疏滤波卷积神经网络的高频信号调制识别算法，该算法首先计算信号的循环频谱作为训练样本，利用稀疏滤波准则，以非监督的方式，逐层预先训练网络，提高了算法的泛化能力。但该算法借鉴深度学习[10-12]的方法，所需训练样本数大，训练时间长。上述方法均是针对定频信号提出的，并不适用于载频伪随机变化的FH 信号，由于FH 信号时频能量谱能够反映出不同调制方式的FH 信号在时域和频域的能量分布情况，因此本文提出了一种基于时频能量谱纹理特征的FH 信号调制方式识别算法。该算法主要由时频能量谱预处理、纹理特征提取和支持向量机（SVM,support vector machine）分类器3 个部分组成。时频能量谱预处理部分通过二维维纳滤波去除背景噪声，再提取每跳FH 信号的时频能量谱并转化为灰度图；纹理特征提取部分通过计算时频灰度图的直方图统计特征和灰度共生矩阵特征组成多维特征向量，作为FH 各调制方式的特征集；SVM 分类器先优化得出最优参数组合，再对特征向量进行训练分类。

2 算法原理

假设接收机接收到M跳FH 信号，可表示为

其中，w(t)是FH 信号与噪声的叠加，sm(t)是第m跳信号，m∈[1,2,…,M]，v(t)是均值为0、方差为σ2的高斯白噪声。假设FH 信号的跳周期为Ts，第m跳对应载频为fm，φTs为起始跳的持续时间，不失一般性，观测到的第一跳信号应该是不完整周期的，因此0＜φ＜ 1。则s m(t)可表示为其中，t′=t-φTs，a是s m(t)的幅度，φ0是s m(t)的初相，rect(t)是单位矩形脉冲，φ(t′)是s m(t)的相位函数，同时决定了FH 信号的调制方式。本文利用不同调制方式的时频能量谱分布差异特点，通过提取时频灰度图的纹理特征，实现FH 信号调制方式分类识别，其计算流程如图1 所示。

图1 计算流程

本文算法的核心思想是将FH 信号调制识别问题转换为时频图识别问题。由图1 可知，数据预处理部分完成对FH 信号时频图的提取去噪和转换分割，特征提取部分计算时频灰度图的直方图统计特征和灰度共生矩阵特征并组成多维特征向量，SVM 训练识别部分完成对纹理特征向量的训练和分类。

2.1 FH 信号时频变换

目前，提取FH 信号时频图的方法主要有短时傅里叶变换（STFT,short-time Fourier transform）、伪Wigner-Ville 分布（PWVD,pseudo Wigner-Ville distribution）、平滑伪Wigner-Ville 分布（SPWVD,smoothed pseudo Wigner-Ville distribution）和Cohen分布（CWD,Choi-Williams distribution）。STFT 的时间和频率分辨率互相制约，时频分辨率低。PWVD是在Wigner-Ville 分布上加时间窗，其模糊域核函数φ(τ,υ)=h(τ)，只能消除模糊平面τ方向上的交叉项影响，υ方向的交叉项影响严重。SPWVD 同时具有时间窗和频率窗，在模糊域的核函数φ(τ,υ)=h(τ)g(υ)，可以同时消除模糊平面τ和υ方向的交叉项，而且h(τ)和g(υ)相互独立，可根据实际FH 信号单独设计时间和频率窗长，从而可获得较高的时频分辨率。CWD 在模糊域的核函数其时间聚集性和频率聚集性不能单独控制，而是由参数σ统一控制。当σ较大时，时频聚集性好，交叉项抑制能力差；当σ较小时，时频聚集性差，交叉项抑制能力好。这导致时频聚集性和交叉项抑制相互制约，而且当τ=0或υ=0时，φ(τ,υ)=0，表明CWD 核函数对模糊平面的τ和υ轴均没有滤波效应。为了对比4 种方法提取FH 信号时频图效果，取一段含有白噪声的FH信号，信噪比为0 dB，采样数据为1 500 个采样点，跳周期为200 个采样点，调制方式为QPSK，采样频率为200 kHz，频率集为[10,45,20,30,5,15,25,30,50]kHz，PWVD 的h(τ)为127 长矩形窗，SPWVD 的h(τ)和g(υ)分别为127 和63 长的Gauss窗，CWD 参数σ=10，各方法时频图如图2 所示。

由图2 可知，对于FH 信号，在相同信噪比条件下，SPWVD 方法的交叉项抑制和时频分辨率更好，因此本文采用SPWVD 方法提取FH 信号时频图。

2.2 时频图去噪与转换

图2 FH 信号时频图

虽然SPWVD 具有一定抗噪性能，但在低信噪比条件下，时频图中仍然存在很多背景白噪声[13]，严重影响纹理特征提取。二维维纳滤波[14]能够根据图像的局部方差自适应调整滤波器，有效降低白噪声影响。因此本文采用二维维纳滤波方法去除时频图噪声。

假设矩阵Hα×β表示α×β的二维图像，H(i,j)表示第i行第j列的像素，i=1,…α,，j=1,…β,。在Hα×β中选取a×b的图像域ηa×b，1＜a＜α，1＜b＜β，则图像域ηa×b均值μ和方差σ2由式(3)和式(4)计算得出，最后滤波恢复图像Pα×β可由式(5)计算得到。

为了提取每跳信号的时频灰度图，本文采用连通域检测算法提取FH 信号每一跳的时频能量谱，再通过式(6)将时频能量谱转化为时频灰度图。

其中，P′(γ,λ)为每跳时频能量谱的第γ行第λ列的值，Q(γ,λ)为转换后时频灰度图的第γ行第λ列的值，round为就近取整运算，L为灰度级，本文中L=256，即Q(γ,λ)=[0,1,…,255]。图3 是4 种不同调制方式的FH 信号单跳时频灰度图。

从图3 可知，不同调制方式的FH 信号时频灰度图具有明显的纹理差异，本文通过提取纹理特征实现FH 调制方式的分类识别。

2.3 纹理特征提取

直方图统计特征用来表示灰度图中每个灰度级出现的频率，而本文时频灰度图的直方图统计特征即可表征每种调制FH 信号的能量谱值大小和概率。假设单跳时频灰度图尺寸为μ ν×，像素点总数n=μ×ν，灰度为l的像素点个数为nl，l=[0,1,…,255]，则直方图统计，其中可提取的6 维直方图特征分别介绍如下。

图3 FH 信号单跳时频灰度图

灰度共生矩阵能够反映图像在方向、幅值和局部邻域的纹理信息分布情况，本文时频灰度共生矩阵即可表征每种调制FH 信号的能量谱在特定方向上的空间分布和能量大小变化，以及局部能量谱的差异对比。假设单跳时频灰度图中位置(x1,y1)的像素为l1，位置(x2,y2)的像素为l2，像素l1和l2的像素对位置关系为(Δx,Δy)，则单跳时频灰度共生矩阵为

其中，l1,l2=[0,1,…,255]，num (l1,l2,Δx,Δy)为位置关系为(Δx,Δy)的像素对(l1,l2)出现的总次数。

同理，可计算出(1,-1)、(1,0)和(-1,-1)这3种位置关系的共生矩阵特征参数[V11,V12,…,V22]，则最终组成本文用来表征不同调制FH 信号时频能量谱纹理特征的22 维特征向量[V1,V2,…,V22]。

2.4 SVM 参数优化

SVM 分类器参数优化[16]是提高SVM 分类模型分类性能和运算效率的重要途径，其核心问题是如何确定合适的初始惩罚参数和核函数参数。其中，惩罚参数决定了SVM 最优分类平面[17]对错误样本的容忍程度。核函数参数决定了样本间隔的尺度范围。通过选择合适的惩罚参数和核函数参数控制SVM 分类器的经验风险和计算复杂度[18]。

在分类过程中，当2 类样本不能被最优平面分类时，SVM 通过引入松弛因子ξ允许错误样本存在，并加入惩罚项C，调整最优平面对错误样本的容忍度，平衡了SVM 的经验风险和复杂度，此时非线性分类问题[19]计算式为

其中，k为样本类别数，ui为第i类样本，vi为分类结果判决输出，αi为拉格朗日乘子。

SVM 为了解决非线性分类问题，通过核函数将原始的k维非线性空间转化到高维线性空间，而且只需要在输入时做点积运算，降低了运算复杂度。

根据文献[19-20]，SVM 在高维空间的分类平面满足其中，φ(ui)为核函数，w为分类平面系数，为分类间隔。在线性分类问题中，寻找最优分类面就是最大化分类间隔，即最小化函数。在非线性空间中引入松弛因子后，由式(9)可知，当0＜ξi＜ 1时，样本ui正确分类；当ξi＞ 1时，样本ui被错误分类。因此，加入惩罚项，可得到高维空间的最小化目标函数为

由于径向基（RBF,radial basis function）核函数不受样本数量及特征维数的限制，而且需要确定的参数更少，复杂度更小，因此本文的SVM 分类器采用RBF 核函数，其计算式为

其中，γ为核函数参数[21]，决定着核函数作用样本间隔的尺度范围，当γ较大时，φ(u i,uj)值较小，导致核函数只能影响样本间隔为的小范围内的样本数据。所以确定合适的γ能有效控制SVM 的经验风险及模型复杂度，提高分类性能。

传统的网格搜索算法（GS,grid search）[22]理论上只要不断扩大搜索空间范围，不断缩短搜索步长，即可找出全局最优解。但在实际应用过程中，优化参数数量的增加、搜索空间范围的扩大、搜索步长的不断减小，会导致搜索算法计算量增加，时效性降低。为了能够得到合适的SVM 初始化参数，降低SVM 的经验风险，控制模型复杂度，本文改进了网格搜索参数优化算法，其具体算法流程如下。

输入α，β，lCstep，lγstep，δ，K，其中，α和β控制C和γ的取值范围，lCstep和lγstep分别为C和γ的初始网格搜索间隔，δ为过拟合临界值，K为交叉验证折数

步骤1初始化参数空间。

步骤2设置搜索网格lCstep,lγstep←λl，其中λ≥ 1，进行网格粗搜索，得到参数优化结果为Ci1，γi1，CVi1。

步骤3增加网格生成密度，

步骤4根据C1i与δ关系，进行局部细搜索。当Ci1＜δ时，在最优等高线区域进行细搜索计算，得出最优参数组合Ci1，i2γ，CVi2。当C1i≥δ时，沿最优等高线边界反向细搜索，具体如下。

步骤5输出最优参数组合：Cbest←Ci2，γbest←γi2，CVbest←CVi2。

本文参数优化算法首先在初始参数空间(C i,γi)上，增加初始步长为λl，减少网格密度，在较短时间内完成首次遍历搜索，得到最初SVM 的参数组合(C1i,γ1i)和评估结果的等高线图。然后在对初次惩罚参数C进行判断时，减少步长为，增加网格密度，进行精细搜索。

3 实验仿真与分析

3.1 信噪比对本文算法分类性能的影响

将本文算法与文献[2]、文献[4]和文献[7]的算法进行对比，分析4 种算法的分类识别性能。文献[2]计算每跳FH 信号的5 个低阶循环谱相关系数ρ1、ρ2、ρ3、ρ4和ρ5组成特征集，利用决策树分类器进行分类。文献[4]先计算信号的4 阶累积量C42、6阶累积量C61和谱相关系数，并将3 种特征组成特征向量，通过BP 神经网络分类器进行训练分类。文献[7]为了减少运算量，利用高阶累积量的稀疏性，通过正交匹配追踪重构算法，计算四阶累积量和6阶累积量Cs63，将和组成特征向量，利用SVM 分类器进行训练分类。

实验中FH 信号跳周期为1 500 个采样点，采样率为1.5 MHz，频率集为[460,100,150,100,225,75,300]kHz，调制方式有BPSK、QPSK、SDPSK、QASK、16QAM 和GMSK 共6 种，初始相位设为0，巴克码长为随机[150,200,250,300]个采样点，加性噪声为高斯白噪声，均值为0，样本长度为9 000 个采样点，每个样本提取5 个完整跳周期的时频灰度图，每个时频灰度图大小为1 500×64，每种调制信号的信噪比为-10～10 dB，每隔2 dB 采集400 个样本作为训练集，共4 400 个训练样本。每种调制信号的信噪比为-10～10 dB，每隔2 dB 采集200 个样本作为测试集，共2 200 个测试样本。将6 种调制信号的测试样本随机排列，信噪比以2 dB 间隔从-10 dB递增至6 dB，每个信噪比下均做20 次分类实验，结果取平均值。本文算法与文献[2]、文献[4]和文献[7]算法的仿真对比结果如图4 所示。

图4 不同信噪比下FH 调制识别正确率

由图4 可知，4 种算法的识别正确率均与信噪比有关，当信噪比大于0 dB 时，本文算法的识别正确率达到94%左右；当信噪比大于-4 dB 时，本文算法的识别正确率大于87%，随着信噪比继续降低，4 种算法的识别正确率均不同程度下降。但是对于QASK、QPSK 和GMSK 这3 种调制方式，当信噪比为-6 dB 时，本文算法的识别正确率仍然大于82%。文献[2]算法由于提取特征集的单一性，当信噪比较低时，相关系数特征波动较大，导致分类器判决阈值设置不合理，识别率下降明显。文献[4]和文献[7]算法均是提取信号的高阶累积量特征，由于FH 信号载频的伪随机变化对特征的稳定性和表征能力有较大影响，导致算法在低信噪比条件下，分类性能下降明显。本文算法在信号预处理阶段增加了二维维纳滤波，降低了时频图背景噪声，提高了算法抗噪性。

当信噪比为-4 dB 时，6 种调制方式的FH 信号混合识别的结果如表1 所示，其整体平均识别正确率为91.4%。从表1 中可知，本文算法对于6 种调制方式信号的识别正确率均较高。识别正确率最低的是BPSK 信号，算法识别误差主要发生在纹理特征相似的调制方式之间，例如SDPSK 和BPSK 调制信号，由于2 种调制信号的时频图在背景噪声的影响下，一些细微纹理特征变得模糊甚至丢失，导致提取的特征表征能力不足，分类正确率下降。

表1 信噪比为-4 dB 时6 种调制FH 信号识别结果

3.2 训练样本数对本文算法分类性能的影响

实验中FH 信号的参数设置与3.1 节实验相同，当信噪比为-2 dB 时，每种调制信号采集800 个训练样本，100 个测试样本。将6 种调制信号的测试样本随机排列，训练样本数以100 步长从100 递增至800，每个训练样本数条件下均做20次分类实验，结果取平均值作为识别正确率。本文算法与文献[4]和文献[7]算法的识别结果如图5 所示。

图5 不同训练样本数下FH 调制识别正确率

由图5 可知，在100 个训练样本条件下，本文算法的识别正确率约为88%，明显高于文献[4]和文献[7]算法。随着训练样本数增加，本文算法的识别正确率增加较少，而文献[4]和文献[7]的识别正确率分别增加了8%和13%，说明本文算法提取的特征集具有更强的稳定性和表征能力，能够在少量训练样本条件下，得到较好的分类结果，有效降低了算法的计算复杂度。

3.3 二维维纳滤波对本文算法分类性能的影响

清晰的时频图是纹理特征提取的关键，本文采用二维维纳滤波方法对SPWVD 时频变换后的时频图进行去噪处理。图6 是当信噪比为-8 dB 时，滤波前和滤波后的GMSK 调制FH 信号时频图。从图6中可知，在经过二维维纳滤波之后，时频图背景噪声被明显消除，保证了FH 信号时频图在低信噪比条件下的清晰度。

图7 是二维维纳滤波对本文算法分类性能的影响，实验中6 种调制信号的测试样本随机排列，信噪比以2 dB 步长从-10 dB 递增至10 dB，每个信噪比下均做20 次分类实验，结果取平均值。从图7中可知，当信噪比小于2 dB 时，经过二维维纳滤波后的SVM 分类正确率明显高于未经过二维维纳滤波的SVM 分类正确率。当信噪比为-4 dB 时，经过二维维纳滤波的SVM 分类正确率约为90%，而未经过二维维纳滤波的SVM 分类正确率约为70%，说明二维维纳滤波有效提高了低信噪比条件下本文算法的分类性能。

图6 二维维纳滤波前和滤波后的时频图

图7 二维维纳滤波对本文算法分类性能的影响

3.4 SVM 参数优化实验

本实验采用GS 算法和本文参数优化算法分别对SVM 的参数进行优化，惩罚参数C和RBF 核函数参数γ的寻优范围均设为，增加搜索步长λ为1.5，利用3.1 节中0 dB 信噪比条件下采集的训练和测试样本进行训练分类实验。经过参数优化后(C,γ)对应的训练集准确率等高线图如图8 所示。

图8 训练集等高线图

从图8(a)中可知，经过GS 参数优化后的SVM惩罚参数C=1、RBF 核函数参数γ=0.015 625，平均分类正确率为73.1%。从图8(b)中可知，经过本文算法参数优化后的SVM 惩罚参数C=1.515 3、核函数参数γ=0.020 382，平均分类正确率达到了92.4%，相比传统的SVM 分类器提高了19%。

4 结束语

本文提出了一种基于时频能量谱纹理特征的跳频调制方式识别算法。该算法提取每跳时频能量灰度图的直方图统计特征和灰度共生矩阵特征组成22 维特征向量，通过参数优化后SVM 分类器进行了训练、分类和识别。实验结果表明，本文算法的识别正确率高于其他2 种算法，尤其在低信噪比条件下识别性能较好，说明本文算法所提取的特征集具有较高的稳定性、充分性和可分性。如何降低多维特征向量的冗余度，提高算法的时效性是下一步研究重点。