中国指挥与控制学会会刊 
在电磁环境日益复杂的今天,低截获概率(Low Probability of Intercept, LPI)雷达得益于自身功率低,带宽大等属性[1],使得侦察接收机难以对其进行有效截获和分析,因此,电子侦察仍面临巨大的挑战且有强烈的现实需求。雷达信号的脉内有意调制类型主要分为频率调制、相位调制,对于相位调制信号,尤其是多相编码信号,具有低旁瓣电平和高距离分辨率特性,且兼具类线性调频信号的多普勒容限性质[2],因此,被广泛应用。相位编码信号通过改变码元的相位达到调制的目的,相比二相编码的相位只有0和π,多相码的相位更加复杂,这使得在电子侦察中对截获的信号进行正确的调制类型识别和调制参数估计十分困难。
早期的雷达信号脉内调制类型识别主要是基于对信号的变换域提取人工特征,文献[3]对维格纳-威尔分布(Wigner-Ville Distribution, WVD)和CWD(Choi-Williams Distribution)时频图像提取瞬时频率和伪Zernike矩等特征,对5种多相码及其LFM信号、Costas码和BPSK信号进行分类识别;文献[4]对模糊函数(Ambiguity Function, AF)、WVD时频图像做Radon变换、分数阶自相关(Fractional Autocorrelation, FA)等,提取5个特征,并分别设置阈值进行调制类型识别,但其多相码调制参数设置单一;文献[5]将LPI雷达信号基于可视图提取的特征输入支持向量机(Support Vector Machine, SVM)进行分类,虽然在低信噪比(Signal-to-Noise Ratio, SNR)下效果斐然,但该特征无法区分相同调制类型的不同调制参数。近几年,深度学习发展十分迅速,因其在图像识别方面效果显著而被引入该领域。文献[6⇓-8]将雷达信号预处理后的CWD图像和模糊函数矩阵送入卷积神经网络(Convolutional Neural Network, CNN)进行特征提取和调制类型识别;文献[9]首先利用SSD(Single-shot Multi-box Detector)网络识别调制类型,再通过一个额外的分类器将易混淆的多相码信号做互相关达到分类的目的,该文献将深度学习和传统方法相结合,是本文主要的对比文章;文献[10]首先通过同步压缩傅里叶变换(Fourier-based Synchrosqueezing Transform, FSST)进行时频转换,分别将三种像素大小不同的时频图像通过CNN进行特征提取,最后,通过SVM分类器对融合后的特征进行调制类型识别,相比其他方法,该方法网络效果更好。
由于多相码信号具备相位跳变的特点,同时兼具类线性调频的特征,是一种非平稳信号。时频分析方法是处理非平稳信号的有效手段,因此,多相码的参数估计方法主要是基于其时频特性,通过提取时频分布的不同特征进行参数估计。文献[11]对拉东魏格纳变换(Radon-Wigner Transform,RWT)和三角函数特征计算调制参数;文献[12-13]分别通过拉东模糊变换(RAT)和分数阶自相关提取特征进行调制参数估计;文献[14]通过改进的维格纳霍夫变换 (Wigner-Hough Transform,WHT)对多相码信号进行检测和参数估计;文献[15-16]分别通过循环自相关函数和循环相关熵谱对多相码信号的载频进行估计。
1 多相码时频脊线特征
1.1 多相码信号及其定义
雷达相位调制信号脉内部分的数学模型为复信号形式,如下所示:
s(t)=A*
其中,A为脉冲幅度,fc为载频,ϕk为第k个码元的相位。
多相码的编码长度也称为码元个数或压缩比,用Nc表示,每个码元内的载频周期数为cpp[1]。常用的5种常规多相码信号包括Frank码、P1码、P2码、P3码、P4码,每个码元的时宽相同,相位调制规律如表1 所示。其中,Frank码、P1码、P2码将脉宽T分为M段,每段分为M个码元,编码长度Nc与步进频率数M存在如式(2)所示的关系。由于P2码为了满足自相关低旁瓣的需求,其M一般取正整数中的偶数,Frank码和P1码的M取正整数。P3码和P4码的编码长度Nc则没有该限制,取值为正整数即可。5种多相码信号的具体调制模型如表1 所示。
Nc=M2
表1 多相码信号相位调制方式 |
| 多相码 | 相位变化规律 |
|---|---|
| Frank码 | ϕi,j= (i-1)(j-1) |
| P1码 | ϕi,j= [M-(2j-1)][(j-1)M+(i-1)] |
| P2码 | ϕi,j= [2i-1-M][2j-1-M] |
| P3码 | ϕi= (i-1)2 |
| P4码 | ϕi= -π(i-1) |
1.2 时频脊线特征
时频分析是常用的雷达信号分析方法。相比时域和频域分析,时频分析可以反映信号频率成分随时间变化的规律,适合于非平稳信号分析。短时傅里叶变换(Short-time Fourier Transform,STFT)中,因为窗长度的选择无法同时保持较高的频率分辨率和时间分辨率,WVD会产生交叉项的干扰,所以,本文提取特征时采用了CWD时频变换的时频分析方法。表1 中的5种多相码信号的CWD等高线如图1 所示。
通过多相码信号的CWD图像,发现Frank码和P3码的时频脊线个数为2,P1码、P2码、P4码的时频脊线个数为1;P2码时频脊线斜率为负数,其他四类的斜率为正数。此外,Frank码、P1码、P2码的时频脊线在时频图像中可分为M段,每段内部斜率小,段与段之间斜率大,这在等高线图体现为时频脊线的轮廓呈“阶梯”状。而P3码和P4码的时频图像没有该现象,十分平滑。这是因为Frank码、P1码、P2码的调制方式与P3码、P4码不同,其相位分为M段,段内相位按照特定规律小幅有序变化,而段间(指上一段最后一个码元与下一段第一个码元)的相位差相差过大,这种相位的大幅无规律变化引起时频脊线段间斜率增加;而P3码、P4码的相位没有分段,相邻码元的相位变化值小,因此,其时频脊线无“阶梯”特征。这5种多相码的离散相位值如图2 所示。
根据以上分析可知,不同多相码具有不同的特征,具体见表2 。
表2 常规多相码特征汇总 |
| 调制类型 | 主时频脊 线个数 | 时频脊 线斜率 | 是否存在 “阶梯特征” |
|---|---|---|---|
| Frank码 | 2 | 正 | 是 |
| P1码 | 1 | 正 | 是 |
| P2码 | 1 | 负 | 是 |
| P3码 | 2 | 正 | 否 |
| P4码 | 1 | 正 | 否 |
1.3 “阶梯”特征提取算法
本文提取的多相码特征包括时频脊线个数、时频脊线斜率以及与“阶梯”特征有关的参数。提取步骤具体如下:
Step1:判断主时频脊线个数A,若A为1,直接对脊线按照下述方法提取特征;若A为2,则对每条主脊线依次提取特征;
Step2:对主时频脊线进行直线拟合,将直线斜率α作为时频脊线斜率的估计值,并保存其拟合的均方根误差E;
Step3:提取主时频脊线的轮廓;
Step4:计算轮廓的上沿部分与下沿部分间距的方差V,并视其为特征之一;
Step5:对轮廓进行平滑拟合,取上下边沿之一,对平滑拟合的曲线进行求导,若时频图像存在“阶梯”特征,则理论上该导数有M个零值段,且每个零值段长度相等,若时频图像不存在“阶梯”特征,则实际中求出的零值段数量少且长度短。因误差导致实际每个零值段长度不一致,分别计算每个零值段的长度,根据长度对其进行K-means聚类,聚类中心设为3,删除零值段数量最少的一类。将剩余零值段的个数N1作为“阶梯”个数的估计值,零值段的均值L1作为“阶梯”平均长度的估计值;
Step6:平滑拟合后,若存在“阶梯”特征,则视每段内斜率为0,段间斜率大于0,提取斜率不为0的段落,对其个数求和得到N2,对其长度求均值得到L2;
Step7:挪动Step2中拟合得到的直线,平移该直线使其分别与时频脊线轮廓相外切,计算上下切线的距离作为一个特征DO; 平移该直线使其分别与时频脊线相内切,计算上下切线的距离DI,将DO与DI的比值R1视为一个特征;计算两条外切线内部的时频脊线面积占比R2;
Step8:若A为1,提取的特征参数可直接应用;若A为2,则所提取特征除A外各有2个值,对其取平均值的有α、L1、L2、DO、R1、R2,取二者之和的有E、V、N1、N2。
本文共提取特征11个,采用的分类器为SVM,SVM是一种常用的机器学习分类器,因其对样本数量需求低且运算量小而得到广泛应用。将训练集提取的特征输入SVM进行模型训练,再将测试集特征输入学习到的模型进行分类,达到调制类型识别的目的。
Step7中的部分时频脊线内外切线如图3 所示。
2 特征提取及参数估计流程
本文提出的多相码信号时频脊线特征用于调制类型识别,其中,直线拟合斜率α能直接用于多相码的脉内调制参数估计。
2.1 带宽及编码长度估计
文献[9]对连续波(Continuous Wave, CW) LPI 雷达信号进行调制类型识别,通过RWT对脉宽等参数进行估计;本文针对的是脉冲信号,脉冲宽度(Pulse Width, PW)在脉间信号分选时已计算得出。文献[9]通过文献[11]提出的方法进行调制参数估计。如图4 所示,文献[11]将连续波雷达信号做RWT后求出时频谱线的夹角α及相隔1个时频谱线的距离d,再通过三角函数关系计算出带宽B和单个多相码调制周期T,最终根据公式(3)求出编码长度。因为WVD会产生交叉项干扰,导致在其时频图中相邻时频脊线中间会产生实际不存在的虚假脊线,所以,可通过计算相隔1个时频谱线的距离来得到单个多相码调制周期T。
Nc=B*T
式中,Nc为编码长度,B为频带宽度,T为脉冲宽度。
2.2 载频和码元内载频周期数估计
本文对多相码时域信号做时频变换,时频分布的峰值对应的纵坐标所代表的频率视为多相码载频的估计值,根据估计的带宽B和载频fc,通过公式(4)对码元内载频周期数cpp进行估计。
cpp=fc/B
2.3 调制类型识别及参数估计流程
3 仿真实验分析
3.1 仿真参数设置
本文通过生成仿真脉内信号进行实验,采样率fs设为1.2 GHz,每个样本的载频取值范围为[fs/5,fs/4],脉内调制参数如表3 所示,每个样本的载频、码元内载频周期数和编码长度均在指定范围内随机产生,则所有样本的长度范围均在[288, 1 280] 区间内。训练集每类信号共产生1 000个样本,每个样本的信噪比均在[20 dB] 范围内根据生成的随机数不定向取值,训练集共包括5 000个样本;测试集信噪比从-2 dB以2 dB的步进值依次增加至20 dB,每类信号在每个信噪比下生成200个测试样本,即测试集共包括12 000个样本。通过调用时频分析工具箱中的相关函数对脉内信号做CWD,采用默认参数,即频率窗长度为每个样本信号长度的1/4,时间窗长度为信号长度的1/10。
表3 多相码信号仿真参数 |
| 调制类型 | 参数值 |
|---|---|
| Frank码、P1码 | 码元内载频周期数cpp:[2, 3, 4] |
| 步进频率数M:[6, 7, 8] | |
| P2码 | 码元内载频周期数cpp:[2, 3, 4] |
| 步进频率数M:[6, 8] | |
| P3码、P4码 | 码元内载频周期数cpp:[2, 3, 4] |
| 编码长度Nc:[36, 49, 64] |
3.2 实验结果及分析
本文对5种多相码信号调制类型进行识别,并对每个样本的调制参数进行估计。1)参数估计部分,对带宽B、编码长度Nc和载频fc和cpp进行了估计;2)调制类型识别部分,将本文提出的方法与文献[9]采用的互相关法进行对比,两种算法的识别结果都是基于本文根据上节所定参数仿真的脉内信号。此外,还将本文方法与深度学习方法进行了对比。
本文对误差的定义包括相对误差(Relative Error, RE)和均方根误差(Root Mean Square Error, RMSE),如式(5)和式(6)所示:
RE= *100%
RMSE=
式中,RE为单个样本参数估计值的相对误差,RMSE为一类样本参数估计值的均方根误差,N为样本数量,δ为待估计参数的真实值, 为估计值。
3.2.1 参数估计
首先,按照提取的参数α对带宽B和编码长度Nc进行估计,然后,分别通过CWD和WVD时频图像对载频fc进行估计,再通过公式(4)对cpp进行估计。数据采用的是3.1节提到的全部测试集样本,共12 000个。B和fc计算了估计值与真实值的均方根误差和相对误差平均值。
参数估计及对比结果如图7 所示。
随着信噪比的提升,B和fc的估计误差逐渐减小, Nc和cpp的正确估计率随之上升。当信噪比为6 dB时,B的估计误差、Nc和cpp的正确估计率基本收敛。在不同的信噪比下,通过WVD对fc进行估计的误差都明显小于通过CWD进行估计的误差,通过WVD图像对fc进行估计的方法更耐噪声干扰,其误差始终保持在较低的程度,因此,本文采用WVD对fc进行估计。
3.2.2 调制类型识别
实验中,将本文方法与文献[9]采用的互相关方法进行了对比,其中,做互相关之前,仿真参考信号时所用的调制参数为上一节的估计结果,识别效果如图8 所示。
随着信噪比的提升,本文方法的平均识别率逐步提高,最终接近100%,在信噪比为2 dB时平均识别率达到90.7%,信噪比为6 dB时平均识别率达到95%;互相关方法的平均识别率始终在80%~85%区间内,没有随信噪比的提升有所改善。随着信噪比的提升,各类调制参数的估计误差逐渐减小,但是,与其相对应的互相关法的平均识别率却没有明显提升,说明互相关法的识别效果不随着调制参数的估计精度提高而线性增长。两种方法在信噪比为2 dB的混淆矩阵如图9 所示。
在信噪比为2 dB时,互相关法的平均识别率为82%,P4码的识别率结果最好,P1码易被错分为P4码,Frank码和P3码易混淆;本文方法的平均识别率为90.7%,P3和P4码的成功识别率相对较低,该信噪比下,本文方法对P2码的成功识别率为100%,且其查全率和查准率均为100%。
3.3 对比实验
由于深度学习处于快速发展阶段,深度学习的方法被广泛应用在各种实际问题中,LPI雷达信号识别问题也不例外。在电子战中,敌我双方的接收机与发射机为非合作方式,因此,在模型训练中,想要获得大量的真实样本十分困难,而深度学习方法却对训练样本数量有要求,因此,本次实验对比了训练样本数量较少时本文方法与深度学习方法的平均识别率。此外,如何能在更短时间内对截获的信号进行识别分析尤为重要,因此,将本文方法与深度学习方法的计算量进行对比,即对比总消耗时长。
本次对比试验的训练集每类样本由10个依次递增到50个,步进值为10,信噪比范围与3.1节中训练集保持一致;测试集样本的信噪比为10 dB,每类多相码有200个样本,共有1 000个测试样本,所有信号调制参数如表3 所示。深度学习网络采用Resnet101网络,Batch size设为32,最大Epoch为10,初始学习率为0.000 1,本文方法的总耗时包括特征提取、训练集测试集标签生成以及SVM模型训练和预测的耗时,深度学习的总耗时包括网络训练及预测耗时,二者计时起点均为预处理后的时频图像,计时终点为测试集预测结束。仿真结果如图10 所示。
当训练集每类信号样本数只有20个时,本文方法的平均识别率为91.6%,明显高于深度学习方法的56.2%,且耗时更少;当训练集每类信号样本数为50个时,两种方法的平均识别率基本相同,但深度学习的耗时约为本文方法耗时的4倍多,因此,在小样本和实时性处理的情况下,本文方法优于深度学习方法。
4 结束语
针对LPI雷达多相码信号易混淆这一问题,本文提出了基于CWD时频脊线的特征提取方法,用于调制类型识别及参数估计。针对时频脊线的“阶梯”现象,本文通过对脊线进行拟合提取特征,通过SVM对特征进行分类,达到调制类型识别的目的,同时,根据特征对调制参数进行估计。相比已有文献,本文选取的多相码参数有随机性,覆盖范围广,更具代表性。参数估计方面,本文分别对带宽、编码长度、载频和码元内载频周期数进行有效估计。调制类型识别方面,在不同信噪比下,将本文方法与经典的互相关法进行了对比,本文方法优于互相关法,并得出互相关法识别率不随着估计精度增加而增长的结论。此外,针对电子战中对实时性和小样本的现实需求,对比了本文方法和深度学习方法在训练样本数不同时的调制类型结果和计算量,在这两个方面本文方法的结果均优于深度学习方法。但本文提取的时频脊线“阶梯”特征在低信噪比下易被噪声干扰,如何在保持现有优势的基础上对序列或图像进行有效降噪是笔者下一步关注的方向。