多跳频信号频率跟踪与二维波达方向实时估计算法

多跳频信号频率跟踪与二维波达方向实时估计算法
张东伟;郭英;张坤峰;齐子森;韩立峰;尚耀波
【摘要】为了实时提取跳频(FH)通信参数以及为通信对抗提供所需信息，该文提出一种多跳频信号频率跟踪和2维波达方向实时估计算法。

首先建立跳频信号的L 型阵列接收数据模型，并推导证明了自回归滑动平均(ARMA)模型对L型阵列数据的适用性，然后采用粒子滤波思想对阵列流型矩阵和频率进行实时估计。

再基于频率估计值建立ARMA模型实时检测跳时刻，并结合流型矩阵估计值实现无需参数配对的2维波达方向(2D-DOA)准确估计。

新方法通过设计合理的粒子生成以及权值更新方式，使流型矩阵与频率估计值能够迅速收敛至稳定状态。

最后蒙特卡罗仿真结果验证了该算法的有效性。

%In order to extract Frequency-Hopping (FH) communication parameters and provide the necessary information for the communication countermeasure, an online estimation algorithm of 2D-DOA and frequency tracking for multiple FH signals is proposed in this paper. Firstly, the data model of the L-array for FH signals is built and the applicability of Auto Regresive Moving Average (ARMA) model to L-array data is proved. Then, the particle filtering is introduced to conduct the online estimation of manifold matrix and the frequency, and the ARMA model is built based on the frequency estimates, depending on which, the online detection of hop timing is obtained. After that , the precise estimation of 2D-DOA can be gained via manifold matrix estimates and without parameter matching. With the rational method of particle generation and the weight updating, the new method makes the estimates of manifold matrix and the frequency reach to the stable value promptly.
Finally the the Monte-Carlo simulation results show the effectiveness of
the proposed algorithm.
【期刊名称】《电子与信息学报》
【年(卷),期】2016(038)009
【总页数】8页(P2377-2384)
【关键词】信号处理;跳频;频率估计;波达方向;粒子滤波;自回归滑动平均模型
【作者】张东伟;郭英;张坤峰;齐子森;韩立峰;尚耀波
【作者单位】空军工程大学装备发展与运用研究中心西安 710051; 空军工程大学信息与导航学院西安 710077;空军工程大学信息与导航学院西安 710077;空军工程大学信息与导航学院西安 710077;空军工程大学信息与导航学院西安 710077;空军工程大学装备发展与运用研究中心西安 710051;空军工程大学装备发展与运
用研究中心西安 710051
【正文语种】中文
【中图分类】TN911.7
跳频(Frequency Hopping, FH)通信具有优越的抗干扰、抗截获和多址组网等能力，已逐渐成为军事通信的主导技术，并在越来越多的信息系统中得到广泛应用[1]。

FH信号跟踪和参数估计一直是通信侦察和通信对抗等领域的重要研究内容。

针对FH信号参数盲估计，已有研究大多基于批处理方式，即需要获得充分的观测数据之后才能从中提取频率和方向等所需信息，因此具有一定的滞后性，难以应用于通信对抗等对侦察实时性要求较高的场合。

如在对FH信号的跟踪式干扰中，就要求侦察系统能够实时跟踪FH信号频率，并迅速估计信号波达方向(Direction Of
Arrival, DOA)以引导干扰分系统实施干扰。

针对FH信号的实时估计问题，文献[8]最早提出了一种基于粒子滤波的单/多通道FH信号频率跟踪方法。

文献[9]在文献[8]基础上建立了FH信号的相位统计模型，在单通道时取得了更好的跟踪性能。

但文献[8,9]所提方法仅能处理一个信号，且
在多通道时未能给出DOA信息，限制了其应用范围。

文献[10]提出了一种基于时域自回归滑动平均(Auto Regressive Moving Average, ARMA)模型的跳时刻实
时检测方法，首次实现了对多个FH信号的频率跟踪，但该方法的信噪比适应能力较弱。

文献[11]将粒子滤波思想与ARMA模型法相结合，提高了文献[10]的信噪
比适应能力，同时可获取DOA信息。

文献[12]采用贝叶斯稀疏学习实现跳频跟踪，其信噪比适应能力略优于文献[10]，但进行跟踪时所需积累的样本数多于文献[10]和文献[11]，且计算复杂度较高。

综上可知，目前可同时跟踪多个FH信号并能获取DOA信息的方法较少，且尚未见到能实时估计FH信号2维波达方向(2D-DOA)的成果报道。

而2D-DOA对FH 网台分选、信号识别和干扰等任务均具重要作用。

针对该问题，本文在文献[11]基础上，利用两个正交布置的均匀线阵(ULA)构造L型阵列，提出一种多FH信号频
率跟踪和2D-DOA实时估计算法，并通过设计更为合理的粒子采样和权值更新方式，提高了文献[11]的频率跟踪性能(跳时刻检测概率和频率估计精度)。

所提算法
对同步、异步网台都能适应，且在估计2D-DOA时无需参数配对。

阵列结构如图1，假定L阵交点为坐标原点，每边是阵元间距分别为, 的元ULA，两ULA分别与+X轴和+Y轴重合(记为ULA1, ULA2)。

若空间有()个FH信号同时入射到该阵列上，其2D-DOA分别为和，并假设在某
个频率驻留时间(Frequency Resident Time, FRT)内各信号频率分别为，定义，则ULA1和ULA2快拍数据可表示为
其中，, 包含了第个信号的幅度和初相信息；, 分别为ULA1和ULA2的流型矩
阵；, 分别为ULA1, ULA2对第个信号的导向矢量，表达式为
观察式(3)，式(4)可见，, 由X轴和Y轴方向的空间相移因子(Spatial Phase Shift Factor, SPSF) , 唯一确定。

只要能够估计各信号的SPSF，就可恢复流型矩阵，进
而估计2D-DOA和频率。

将两子阵数据进行拼接，写成矩阵形式，有
式中为整个L型阵列的流型矩阵。

若某信号在某时刻，其频率从跳变至，则在跳频发生后阵列流型矩阵会改变，快拍数据也不再满足跳变前的形式。

多FH信号频率跟踪与2D-DOA实时估计就是基于，实现对这些跳频发生时刻的实时检测，并在跳频发生后迅速估计各信号的频率和2D-DOA。

文献[10]证明了当个FH信号同时入射到ULA时，各FRT内连续+1个时刻的快拍数据之间符合阶ARMA模型，即
其中，为噪声矢量，, 为[11]
并且，是式(8)所示一元方程的系数。

令式(8)中，，并在两边同时乘以，再将这个等式累加，整理后满足
同理，在式(8)两边乘以并累加，有
可见，对于ULA1, ULA2，其快拍数据均满足阶ARMA模型。

由于ARMA模型系数由各信号频率唯一确定，而信号频率具有唯一性，因此。

根据式(5)，有式(11)
关系成立。

即L型阵列仍满足ARMA模型。

假设各通道及各时刻的噪声均相互独立，服从均值为0，方差为的高斯分布，即，则式(11)的分布函数为
依据前面连续个时刻的快拍数据结合式(12)，可得当前时刻的快拍数据预测值为
时刻实际快拍与之间的偏差也满足高斯分布，即
若在该时刻有一个或多个FH信号发生跳频，将显著偏离，因此可利用实际快拍数
据与其预测值之间的偏差来判断是否发生了跳频，即
其中，表示假设“发生跳频”，表示假设“未发生跳频”；为判决门限，计算方式可参考文献[10]。

对于侦察方，利用快拍数据仅有可能得到阵列对各信号的SPSF。

文献[11]引入粒子滤波[13]估计SPSF，但设计的粒子采样方法未充分利用阵列信息，导致需经多次粒子更新才能收敛至稳定估计值，故在信噪比较低或跳周期较短时，其性能迅速下降。

本节提出一种新的采用粒子滤波的SPSF估计方法。

(1)空间相移因子的后验概率密度形式：以ULA1为例，未发生跳频时，假设已采集到组快拍，则这组样本服从独立高斯分布，即
假设先验分布服从均值为0，协方差为的高斯分布[14]，即
其中，表示观测信噪比(实际应用中可取无穷大)。

若SPSF和噪声方差分别服从内的均匀分布和参数为的逆伽马分布，将流型矩阵简记为，根据贝叶斯定理，未知参数, , 关于的联合后验概率为
式中，, 。

对式(18)关于和积分，并省略常数项，化简后得到关于的后验概率密度函数为[11]
式中，表示矩阵求迹； ; 为协方差矩阵估计值，即
同理可得的后验概率，其与式(19)的形式相同。

(2)空间相移因子的序贯估计
步骤1 粒子初始化。

设定粒子维数、数目，产生初始粒子群并将所有粒子权值置为；
步骤2 粒子及权值更新。

在时刻更新粒子并根据后验概率密度式(19)计算粒子权值；
步骤3 SPSF估计。

对粒子加权平均求得SPSF；
步骤4 重采样。

根据需要，对粒子进行重采样；
步骤5 时刻，转到步骤2。

(a)粒子初始化：利用ULA1的快拍构造如式(21)Toeplitz矩阵
式中，, 。

由于，则式(21)可写为
式中，
为噪声矩阵(形式见文献[15])。

显然的秩等于[15]。

由于为Vandermonde阵，因此可保证的秩也为。

此时对进行特征值分解得到噪声子空间，并构造关于空间相移因子的MUSIC谱，有
同理，构造关于空间相移因子的MUSIC谱：
其中，为ULA2对应Toeplitz矩阵的噪声子空间。

对式(24)和式(25)进行搜索，提取前个峰值，作为和的估计值, ；再在以, 为均值，为方差的正态分布函数上分别采样个样本得到个维粒子矢量，记为
经过以上分析可见，相比于文献[11]，本文采取的粒子生成方法充分利用了阵列信息，使得初始粒子更加接近SPSF真实值，利于后续粒子更新环节SPSF估计值迅速收敛至稳定状态。

(b)粒子更新与空间相移因子估计：假设时刻粒子状态为, ，当接收到时刻快拍后，由于其与该FRT内所积累的历史样本具有相同参数，故只需在这组粒子上附加随
机高斯扰动(为扰动方差)以补偿参数估计值可能存在的误差，即
新候选粒子的接受概率为，其中
式中为新候选粒子生成的流型矩阵；为的广义正交投影矩阵。

用更新的粒子生成流型矩阵并依次代入式(19)，选择使式(19)取最大值的作为最优粒子，归一化粒子权值计算方法为
则时刻ULA1对各信号的SPSF估计值为
设置有效粒子数阈值，若有效抽样尺度，则进行重采样，这样能够自适应地根据粒
子情况决定是否需要重采样，避免了文献[11]在每个时刻都进行重采样所带来的粒子枯竭和计算负担问题。

ULA2粒子更新、权值计算和相移因子估计方式与式(28)~式(32)同理。

假设SPSF估计值为, ，并定义其对应的流型矩阵分别为, ，由式(16)可知当前FRT 内信号波形估计值为
同时，相邻时刻波形之间满足
因此各信号频率估计值为
其中，表示时刻信号波形估计值的第个元素，表示取相位运算。

由于在任一时刻，估计得到的和，其个元素是随机排序的，而式(33)~式(35)估计
结果与SPSF在排序上具有一一对应关系，因此需首先对和进行配对，使两者排序一致，以保证信源和频率估计值的唯一性。

假设ULA1, ULA2的SPSF估计值分别为, ，频率估计值分别为, 。

由于频率具有
唯一性，且在同等条件下ULA1, ULA2的频率估计结果相对稳定，因此可通过调
整中各频率排序使其与方差最小，实现频率配对，即
式中表示与配对后的结果。

再基于频率与SPSF之间的对应关系调整各元素排序，最终实现SPSF以及频率的正确配对。

经过上述配对过程，将与进行平均得到频率估计值，则俯仰角和方位角估计值为
由式(37)、式(38)可见，俯仰角和方位角自动配对。

将代入式(7)计算系数，建立ARMA模型用于下一时刻的跳频检测。

在每个FRT 内，随着快拍数的积累，SPSF估计精度会有所提高，频率和2D-DOA估计精度
也会逐步改善。

(1)算法流程总结：
步骤1 分别利用ULA1, ULA2在初始时刻的快拍数据构造各自的Toeplitz矩阵和，根据式(24)~式(27)生成两组(每组个)维粒子；
步骤2 每接收一组新的快拍，分别(并行)对两组粒子及其权值进行更新，并计算SPSF估计值。

(a)依据式(28)生成新粒子；
(b)由式(29)计算新粒子接收概率；
(c)由式(19)结合式(30)~式(32)计算更新后各粒子的归一化权值并估计各信号的SPSF；
(d)计算，根据需要，进行粒子重采样；
步骤3 利用SPSF估计值恢复流型矩阵，根据式(33)~式(35)估计频率，并与SPSF估计值进行配对；
步骤4 利用式(7)计算系数建立ARMA模型，通过式(37)，式(38)估计各信号的
2D- DOA；
步骤5 通过式(13)和式(15)进行跳频检测。

(a)若未检测到频率跳变，回到步骤2；
(b)若检测到频率跳变，先回到步骤1再转入步骤2，即开始新SPSF内的参数估计。

(2)运算复杂度分析：以复数乘次数来考察运算复杂度，忽略峰值提取(可通过数值比较实现)、反正弦、反正切函数以及标量幂运算计算量，矩阵求逆采用原始定义。

粒子初始化的复杂度为，其中表示, 的搜索范围；粒子及权值更新复杂度为 (表示
此FRT内已积累的快拍数，表示由某粒子生成广义正交投影矩阵的计算量)；SPSF估计的运算量为; 和粒子重采样(多项式重采样)运算量基本可忽略；频率和
2D-DOA估计运算量为; ARMA模型建立与跳频检测运算量为；故算法总运算量
约为，若本时刻发生跳频，则时刻运算量还应加上。

由总的运算复杂度表达式可知，减小阵元数和粒子数均可降低复杂度，然而，和的选取都对参数的实时准确估计有一定影响。

当信噪比不是很低(如0 dB)时，在本
文仿真条件下，综合考虑复杂度以及参数估计性能，可取, 。

实验1 对异步网台信号实时估计的正确性
空间存在两个异步FH信号，其2D-DOA分别为和, , 均为频率400 MHz信号波
长的一半；信噪比12 dB, SPSF抖动参数；观测时间内信号1的跳频次序为，其FRT分别为18, 30, 30个采样点；信号2的跳频次序为，其FRT分别为30, 30, 18个采样点。

图2(a)是频率的实时跟踪结果，可见：算法分别在第19, 31, 49, 61个样点处检测到跳频发生，与真实跳频情况一致；算法能够实时准确估计多个异步FH信号频率，且频率估计值能迅速收敛至真实值，验证了新粒子生成方法的有效性。

图2(b)和
图2(c)是2D-DOA的实时估计偏差，可以看出，算法能够实现对2D-DOA的实
时精确估计，在各FRT稳定期，估计误差小于。

相比于频率偏差，2D-DOA偏差的震荡幅度略大，这主要是由误差累积导致的，但估计误差已能满足绝大多数场合应用需要。

实验2 对同步网台信号实时估计的正确性
观测窗包含完整的3跳信号，两信号同时在第31, 61采样点处发生跳频，其余条
件与实验1相同。

频率跟踪结果如图3(a)所示。

可见，算法分别于第31和61样点处检测到跳频发生，与真实跳频情况一致；算法能够实时准确估计多个同步FH信号频率，且频率估计值能迅速收敛至真实值。

对比图3(a)和图2(a)发现，总体上看，同步组网时，频率估计精度更高，这主要由于在同等条件下(采样率、跳周期均相同)，同步组网时各FRT相对更长(异步组网时，两信号到达时刻随机，驻留时间长度随机)，更利于算法收敛。

图3(b)和图3(c)是2D-DOA的实时估计偏差，可见，本文方法能够对同步FH信号的2D-DOA进行实时精确估计；在各FRT的稳定期，估计误差小于
实验3 验证算法对跳时刻的检测性能
采用正确检测概率评价跳时刻检测性能。

两信号跳频集与实验1、实验2相同，观测样点数设为32，信号1在第17个采样点发生跳频，信号2在观测时间内频率
保持不变。

每次试验，两信号的跳频频率均在跳频集中随机选择。

正确检测的定义为：检测到跳频发生次数为1，且位置为第17或第18采样点。

不同信噪比下进
行100次独立实验，本文方法与文献[11]的正确检测概率如图4所示。

由图可见，本文方法的跳时刻检测性能优于文献[11]，主要原因在于：一是采取新的粒子生成方法使初始粒子分布于真实空间相移因子周围的概率显著增大，易于快速收敛，导致频率估计值和求解的ARMA模型系数更准确；二是L型阵列相比于ULA可获得信源更多的有用信息。

实验4 频率和2D-DOA的估计性能
仿真条件同实验3，发生跳频(第1，第2，第17，第18)的采样点不参与统计，
各参数偏差的定义如式(39)~式(41)：
其中，表示第k次试验中第n个信号在第l个采样点的频率真实值，为对应估计值；K为试验次数，N为信号数，L为样点总数，为跳频带宽；, 分别为第n个信号俯
仰角和方位角真实值；和分别表示第k次试验中第n个信号在第l个采样点的俯仰角和方位角估计值。

文献[11]与本文算法的参数估计性能曲线如图5，图6所示。

图5表明，本文方法对频率的估计性能优良(信噪比10 dB以上时，估计偏差小于跳频带宽的1%)，且优于文献[11]，原因在于L型阵列可供利用的有效信息更多，且对两ULA的频率估计值进行了平均。

图6表明，算法对2D-DOA的估计精度
随着信噪比的提高逐渐变好：信噪比10 dB以上时，方位角和俯仰角的估计偏差
均小于。

7 结束语
在电子对抗实时情报支援等对侦察实时性要求很高的场合，需快速估计空间多个跳
频信号的技术参数，本文基于L型阵列提出一种对同步和异步网台均适用的多跳频信号频率跟踪和2维波达方向实时估计算法。

该方法频率跟踪速度快，2D-DOA估计准确，且无需参数配对。

理论分析和仿真实验表明，所提方法可以很好地解决FH信号多参数的实时估计问题。

如何使所提方法适用于信号源数目随机变化以及更低信噪比环境有待进一步研究。

[1] ZHAO Lifan, WANG Lu, BI Guoan, et al. Robust frequency- hopping spectrum estimation based on sparse Bayesian method[J]. IEEE Transactions on Wireless Communications, 2015, 14(2): 781-793.
[2] 张东伟, 郭英, 齐子森, 等. 多跳频信号波达方向与极化状态联合估计算法[J]. 电子与信息学报, 2015, 37(7): 1695-1701. doi: 10.11999/JEIT141315. ZHANG Dongwei, GUO Ying, QI Zisen,et al. Joint estimation algorithm of direction of arrival and polarization for multiple frequency-hopping signals[J]. Journal of Electronics & Information Technology, 2015, 37(7): 1695-1701. doi: 10.11999/JEIT141315.
[3] FU Weihong, HEI Yongqiang, and LI Xiaohui. UBSS and blind parameters estimation algorithms for synchronous orthogonal FH signals [J]. Journal of Systems Engineering and Electronics, 2014, 25(6): 911-920. doi: 10.1109/JSEE.2014. 00105.
[4] 赵新明, 金艳, 姬红兵. 稳定分布噪声下基于Merid滤波的跳频信号参数估计[J]. 电子与信息学报, 2014, 36(8): 1878-1883. doi:
10.3724/SP.J.1146.2013.01436.
ZHAO Xinming, JIN Yan, and JI Hongbing. Parameter estimation of frequency-hopping signals based on Merid filter in a stable noise environment[J]. Journal of Electronics & Information Technology, 2014,
36(8): 1878-1883. doi: 10.3724/ SP.J.1146.2013.01436.
[5] ANGELOSANTE D, GIANNAKIS G B, and SIDIROPOULOS N D. Sparse parametric models for robust nonstationary signal analysis: Leveraging the power of sparse regression[J].IEEE Signal Processing Magazine, 2013, 30(6): 64-73.
[6] SHA Zhichao, HUANG Zhitao, ZHOUYiyu, et al. Frequency-hopping signals sorting based on underdetermined blind source separation[J]. IET Communications, 2013, 14(7): 1456-1464.
[7] FU Kuoching and CHEN Yungfang. Blind iterative maximum likelihood-based frequency and transition time estimation for frequency hopping systems[J]. IET Communications, 2013, 7(9): 883-892.
[8] VALYAKIS A, TSAKONAS E E, SIDIROPOULOSN D, et al. Stochastic modelling and particle filtering algorithms for tracking a frequency-hopped signal[J]. IEEE Transactions on Signal Processing, 2009, 57(8): 3108-3118.
[9] 李世鹏, 赵正予, 黄硕, 等. 一种基于粒子滤波的跳频信号频率跟踪技术[J]. 华中科技大学学报(自然科学版), 2011, 39(9): 33-37. doi:
10.13245/j.hust.2011.09.009.
LI Shipeng, ZHAO Zhengyu, HUANG Shuo, et al. Frequency tracking technology for frequency-hopping signals using particle filters[J]. Journal of Huazhong University of Science and Technology (Natural Science Edition), 2011, 39(9): 33-37. doi: 10.13245/j.hust.2011.09.009.
[10] LIU Zhangmeng, HUANG Zhitao, and ZHOU Yiyu. Hopping instants detection and frequency tracking of FH signals with single or multiple
channels[J]. IET Communications, 2012, 6(1): 84-89.
[11] SHA Zhichao, LIU Zhangmeng, HUANG Zhitao, et al. Online hop timing detection and frequency estimation of multiple FH signals[J]. ETRI Journal, 2013, 35(5): 748-756.
[12] 王丰华, 沙志超, 刘章孟, 等. 基于贝叶斯稀疏学习的多跳频信号频率跟踪方法[J]. 电子与信息学报, 2013, 35(6): 1395-1399. doi:
10.3724/SP.J.1146.2012.01493.
WANG Fenghua, SHA Zhichao, LIU Zhangmeng,et al.A frequency tracking method for multiple frequency-hopping signals based on sparse bayesian learning[J]. Journal of Electronics & Information Technology, 2013, 35(6): 1395-1399. doi: 10.3724 /SP.J.1146.2012.01493.
[13] LAROCQUE J R, REILLY J P, and NG W. Particle filters for tracking an unknown number of sources[J]. IEEE Transactions on Signal Processing, 2002, 50(12): 2926-2937.
[14] ANDRIEU C and DOUCET A. Joint Bayesian model selection and estimation of noisy sinusoids via reversible jump MCMC[J]. IEEE Transactions on Signal Processing, 1999, 47(10): 2667-2676.
[15] 毛维平, 李国林, 谢鑫. 单快拍数据预处理的相干信号波达方向估计[J]. 上海交通大学学报, 2014, 48(10): 1362-1367.
MAO Weiping, LI Guolin, and XIE Xin.DOA estimation of coherent signals with preprocessed single snapshot[J]. Journal of Shanghai Jiaotong University, 2014, 48(10): 1362-1367.
张东伟：男，1987年生，博士，讲师，研究方向为通信信号处理及阵列信号处理.
郭英：女，1961年生，博士，教授，博士生导师，研究方向为通信信号处理、自适应信号处理等.
张坤峰：男，1989年生，博士生，研究方向为通信信号处理及阵列信号处理.。