在军事领域，面对复杂人工红外干扰的现代战场环境，红外成像导引头的空中目标跟踪技术是影响红外空空导弹抗干扰与精确制导打击能力的关键技术之一。随着红外对抗技术的不断发展，红外空空导弹面临着复杂新型人工红外干扰、机载红外干扰装备等战场干扰环境，尤其在由远及近的过程中，目标发生转弯、快速形变等剧烈运动，人工干扰的红外辐射特性导致导引头视场中目标严重遮挡，红外空空导弹的作战能力面临严峻挑战。

近些年，随着Henriques等^[1-2]提出的相关滤波跟踪算法在可见光图像目标跟踪领域的广泛应用，其高效、精确的跟踪性能也吸引了众多学者进行改进并应用于红外目标跟踪。KCF^[1]引入了多通道的方向梯度直方图(HOG)特征，体现了图像的边缘梯度信息。文献[3]提出一种融合深度特征的空间正则化判别相关滤波(DeepSRDCF)跟踪算法，将卷积特征考虑进来，在判别相关滤波跟踪框架下使用卷积神经网络(CNN)卷积层中的激励。RPT算法^[4]提出一种基于分块模型的方法，核心思想是通过判断各个分块间的相对位置是否达到设定阈值，进而估计目标尺度。Wang等^[5]将相关滤波算法和支持向量机(SVM)相结合，并提出利用平均峰值相关能量(APCE)指标为依据来自适应地更新滤波模型，能够在一定程度上有效缓解模型的漂移。Ma等^[6]提出了一种长时相关跟踪(LCT)算法，在跟踪失败时，通过检测机制对目标进行重定位，使跟踪器恢复到正确的跟踪状态，从而实现长时跟踪。目前的深度学习技术不仅局限于特征提取部分，而且一些算法^[7-10]还结合特征提取与分类器，使用CNN替代相关滤波的整个过程，从而建立端到端输出的跟踪框架，在公开数据集的评测中效果显著。但会极大地影响算法的执行速度，难以满足在线跟踪场景下的实时性需求。

目前主流的核相关滤波跟踪算法及其改进算法，均在空间域提取目标图像特征，没有挖掘目标在频率域中的深层信息。本文提出一种基于二维频域Gabor特征的相关滤波空中红外目标跟踪方法(简称GF-KCF)，该算法与传统相关滤波算法不同之处在于：①通过引入二维Gabor形状纹理特征克服了红外图像的缺陷；②由于二维频域Gabor特征的尺度变化特性，引入了一种利用频域特征进行目标尺度信息估计的方法，克服了传统相关滤波器在空域提取特征的局限性。③引入了高置信分块的跟踪模型解决目标部分遮挡的问题，同时提出了遮挡区域的预测机制，提高跟踪鲁棒性。

1 基于频域GF特征的跟踪模型

已有的基于KCF跟踪算法框架下的目标尺度估计方法，在本质上均是利用空域特征信息^[10]分析和提取目标尺度变化信息，存在实时性较差、对目标变化适应性较弱等缺点。频域GF特征既具有目标边缘纹理特征信息，同时在频域下具有尺度变换方面的特性，因此首先分析和推导了频域GF特征的尺度缩放特性、平移尺度不变性和旋转尺度特性，进而提出一种基于二维频域GF特征的目标尺度信息估计方法。

1.1 频域GF理论

空域中的标准二维Gabor函数一般表达式为：

(1)

式中: g(x, y, σ, β)为二维Gaussian函数; η(x, y, ξ₀, ν₀)为一个二维调制函数; ξ₀, ν₀为函数在频域的坐标。频率为(ξ₀, θ₀)的复平面波的传播方向为椭圆Gaussian函数的短轴方向。将频域下的直角坐标(ξ₀, ν₀)转换为极坐标(ω₀, θ), 并细化二维Gabor滤波器的模型参数, Gabor滤波器的函数模型可进一步得到

(2)

式中：ω₀为复平面波的尺度；θ为复平面波的传播方向；Ω为阶跃带宽, x′，y′分别为x，y旋转角度θ后的坐标。对于特定的Ω, k也是固定常数。(2)式即为二维Gabor滤波器函数模型, 选择特定的尺度、方向的Gabor滤波器函数分别提取相应的目标特征。应用Gabor滤波器组对图像进行滤波处理, 其空域响应值q_s为图像函数f(x, y)和Gabor滤波器函数ϕ(x, y, ω₀, θ)的卷积, 即

(3)

根据Fourier变换的卷积定理:2个函数在空域的卷积等于函数Fourier变换的乘积。因此, 可以将二维函数和二维Gabor滤波器函数转换到频域, 然后将2个函数直接相乘。等价于在频域下对二维函数进行Gabor滤波, 得到频域GF特征Q_s(ξ, ν, ω₀, θ), 即

(4)

式中：函数Q_s(ξ, ν, ω₀, θ)，F(ξ, ν)，Ψ(ξ, ν, ω₀, θ)分别为q_s(x, y, ω₀, θ)，f(x, y)，ϕ(x, y, ω₀, θ)的二维Fourier变换表达式。在使用不同尺度Gabor滤波器提取纹理特征时, 较大尺度的Gabor滤波器提取出的纹理较为细密, 目标边缘纹理特征较为准确, 但是易受到原图中噪声的干扰; 相反, 较小尺度的Gabor滤波器提取的纹理较为粗糙, 只有目标的大致边缘轮廓, 对目标精确定位较为困难, 但是可以抑制图中噪声。

由于Gabor特征本身具有非正交性, 这些频域响应值包含了大量的冗余信息, 需要将响应值中最能体现目标的特征提取出来。根据二维Gabor滤波器参数选择方案, 对于频域Gabor滤波响应矩阵位置

(5)

式中：J_m(ξ, ν)即为所有频域Gabor滤波响应矩阵位置处的响应最大值；i(i=0, 1, …, K-1)，j(j=0, 1, …, S-1)代表Gabor滤波器的不同尺度和方向, 每一个目标图像块经(5)式滤波后得到一组频域GF特征Q_s(ξ, ν, ω_i, θ_j)。对频域响应矩阵的所有位置点进行筛选、提取、融合成频域响应矩阵, 即

(6)

式中：m, n分别为响应融合矩阵的行数和列数。选择出对应频率下的最大响应值即该最大值对应的尺度和方向滤波响应最大, 即该尺度和方向下的纹理较为丰富, 更易于提取目标的边缘纹理信息。

1.2 频域GF特征的尺度特性分析 1.2.1 尺度缩放特性

由1.1节可知, Q_s(ξ, ν, ω₀, θ)和q_s(x, y, ω₀, θ)分别为频域和空域下的Gabor滤波响应矩阵, 即

(7)

假设m、n分别为x方向和y方向的尺度因子。为了推导出目标在空域发生尺度变化对于频域的影响, 令mx=τ₁，ny=τ₂, 同时考虑到符号m，n的组合情况, 即mn>0和mn < 0, 最终整理可得到

(8)

由(8)式可以看出, 在相位谱中, 尺度缩放因子m、n在ξ、ν方向上相互独立且变化方向相反。在幅度谱中, 尺度因子m、n相互耦合, 即2个方向的尺度缩放共同影响GF特征在频域幅值。同时, 由能量守恒原理, 二维尺度缩放后的幅值变为原来的1/mn倍。

1.2.2 尺度平移不变性

空域GF特征q_s(x, y, ω₀, θ)发生相应的平移变化, 其频域表达式推导过程如下

(9)

令x-x₀=τ₁、y-y₀=τ₂, 带入(10)式可得

(10)

式中，x₀、y₀的任意符号情况均满足等式(10)。分析可知, 二维函数在空域发生平移, 在频率域中只发生相位变化, 而其频谱的幅值响应不变。

1.2.3 尺度旋转特性

根据二维Fourier变换, 令 arctan(ξ/ν), 可以得到ξ=ωcosφ, ν=ωsinφ。则空域的GF特征q_s(x, y, ω₀, θ)的Fourier变换极坐标表达式为

(11)

式中：μ为复平面的频率；φ为复平面波的传播方向。以GF特征中任意一点作为旋转坐标系中心, 在空域旋转角度δ后的GF特征为

(12)

式中：x，y经过坐标旋转变换分别得到x′，y′。对等式(12)两边进行二维Fourier变换到频域, 同时推导和建立x与x′，y与y′的等式关系, 可得

(13)

因此, 在空域发生旋转变化的GF特征, 在频域中的的幅值响应不变, 其相位发生对应的旋转变化。

1.3 频域GF特征尺度信息估计方法

结合1.2节分析, 本文提出一种频域尺度信息估计方法, 分别提取径向、角向特征向量, 并结合帧间特征变化情况进行目标尺度估计, 具体包括旋转角度估计、尺度信息估计、跟踪框比例更新3个部分。

1) 旋转角度估计。由二维幅度谱关于频域原点中心对称, 提取当前帧跟踪框幅度谱0°~180°的角向特征值, 组成角度估计特征向量, 如图 1所示。

将得到的特征向量进行归一化处理, 进行循环移位r, 并与上一帧特征向量ϕ^t-1={ϕ₁^t-1, …, ϕ_i^t-1, …, ϕ_N^t-1}进行SAD模板匹配

(14)

响应值最小的位置r即为目标旋转角度估计值θ_r。

2) 尺度估计。当目标发生尺度变化时, 由(8)式可得, 振幅谱原点只发生振幅值的相应尺度变化, 其频谱相位不变。由于目标是刚体, 假设相邻2帧之间水平和垂直方向尺度变化率相等。因此, 首先对尺度变化率m、n进行粗估计: 。

为了获得更精确的目标尺度m，n参数估值, 提取上一帧频谱图 0°的频谱离散曲线作为径向特征模板T^t-1, 根据旋转角度估计值θ_r提取当前帧对应θ_r的频谱离散曲线作为径向特征向量T^t。根据尺度粗估计得到目标参数的初值, 使用非线性最小二乘拟合的算法通过多次迭代来求解。采用ASDF估计进行参数估计的细化过程

(15)

高斯牛顿迭代法求解为

(16)

式中：F=(F₁, …, F_N)^T；J_r为雅可比矩阵。根据(16)式的迭代表达式对参数向量θ进行多次迭代运算, 当‖θ^(s+1)-θ^(s)‖低于预设的门限值, 迭代终止。

3) 跟踪框更新。计算出的m，n即为x，y方向的2个尺度比例因子。根据上一帧跟踪框尺度w，h估计出当前帧目标跟踪框尺度w′，h′, 定义对角线长度为, 对角线与横轴夹角为θ_l, 因此

(17)

由于目标旋转会改变跟踪框的长宽比, 仍需要考虑旋转对跟踪框尺度的影响。根据目标旋转角度更新跟踪框尺度

(18)

2 高置信度分块的抗干扰跟踪策略

通常, 高置信度分块具有2个属性:高置信跟踪、与目标物体黏连。从边界框x=[x, y, w, h]∈R⁴采样图像分块, 观察前一帧中的z_1:t-1={z₁, z₂, …, z_t-1}, 推导分块置信度概率密度函数, 分块的跟踪可靠性p(z_t|x_t)可以表示为

(19)

式中：p_t(z_t|x_t)表示有效跟踪的分块置信度；p_o(z_t|x_t)表示分块在被跟踪对象上的可靠性。本文使用粒子滤波器来推断出高置信度的分块是否位于被跟踪的物体上, 通过相关滤波算法跟踪高置信分块。

为了估计有效分块的跟踪可靠性, 采用峰值旁瓣比(PSR)作为置信度。基于相关滤波跟踪算法估计目标位置

(20)

(21)

由于本文算法选择相关滤波作为基本跟踪算法框架, R(X)表示相关滤波器响应图, Φ表示峰值周围的旁瓣面积。μ_Φ和σ_Φ分别表示不包括旁瓣区域Φ的响应图R的均值和标准差。将s(X)视为分块算法测量其是否正确跟踪的置信度, 跟踪可靠性如(21)式所示, 其中是衡量似然函数贡献的系数, 本文中取λ=2。

同时, 利用运动信息计算分块位于被跟踪对象上的概率。跟踪前景和背景分块, 并记录每个分块的相对轨迹:V_t=[v_t-k+1^T, …, v_t^T]^T∈R^2k, 其中v_t=Ψ₂(x_t-x_t-1)表示相对移动矢量, Ψ₂=[E_2×2, 0]∈R^2×4表示用于选择初始状态下的位置向量的选择性矩阵。由于分块的位移可能对应于不同的对象, 因此记录了连续k帧的相对运动矢量信息, 以使其跟踪更具有鲁棒性, 同时采用2-范数来测量轨迹之间的距离。通过将相似性得分函数表述为以下方式来衡量分块与其标记组之间的相似性

(22)

式中：y_i∈{+1, -1}是指示x_i是否位于目标上的标签；Ω_t⁺是包含正样本分块索引的集合；Ω_t^-是包含负样本分块索引的集合；N⁺和N^-分别是响应集合的数量。计算分块在运动目标上的概率如(23)式所示, 其中μ是平衡在目标上的概率所贡献的系数, 本文中取μ=1。

(23)

3 GF-KCF空中红外目标跟踪算法

在相关滤波跟踪算法的框架下, 本文提出一种基于频域GF特征的KCF空中红外目标跟踪算法(简称GF-KCF)。同时结合标准导引头全过程抗干扰跟踪架构, 根据算法功能和目标当前状态选择相应的跟踪策略, 实现全过程抗干扰跟踪。算法整体架构如图 2所示。

首先判断目标部分区域是否被遮挡, 即是否为干扰态。若目标无遮挡, 则进入全局跟踪模式, 并对目标进行频域尺度估计; 若目标发生遮挡, 则切换至局部跟踪模式, 采用高置信分块模型进行抗部分遮挡跟踪。GF-KCF算法的实现可以分为构建频域Gabor滤波器组、全局/局部跟踪、跟踪框和模板更新3个部分。下面对每个部分的原理进行详细描述。

1) 构建频域Gabor滤波器组。初始帧由高斯窗函数和调制函数得到频域Gabor滤波器, 根据1.1小节Gabor滤波器参数选择方案, 变换滤波器尺度和方向得到频域滤波器组, 用于提取频域GF特征。

2) 全局/局部跟踪。通过上一帧目标响应最大值判断目标是否发生遮挡, 将GF-KCF算法分为全局跟踪和局部跟踪。

(1) 全局跟踪:如图 3所示, 在目标无干扰遮挡情况下, 算法采用全局跟踪器。训练阶段完成以后, 对于第t帧, 以上一帧目标位置为中心、大小为M×N的图像块z上进行频域GF特征提取、融合, 新的目标位置即为响应图f(z)最大响应位置。根据第(t-1)帧跟踪框尺度大小, 提取第t帧目标频域GF特征的尺度特征向量。基于相邻2帧的频谱振幅特征向量估计目标当前帧旋转角度θ_r。设原频域GF响应矩阵为Q(x, y), 旋转后的响应矩阵为Q_r(x, y)=Q(x′, y′), 并由此重新计算当前帧频域特征向量。基于相邻2帧的特征向量, 计算目标尺度比例因子m，n, 完成一次目标尺度估计。

(2) 局部跟踪:如图 4所示, 在目标发生部分遮挡情况下, 算法采用高置信分块模型进行局部抗干扰跟踪。检测各分块目标位置, 更新高置信分块模板并重采样低置信分块, 根据所有高置信分块估计目标位置和尺度。将粒子滤波器的归一化粒子权重w_t⁽ⁱ⁾视为跟踪结果的置信度, 采用所有高置信的分块来对目标的中心进行投票, 如(24)式所示

(24)

最终的跟踪目标状态可以估计为

(25)

3) 跟踪框和模板更新。为了更好适应目标外观和旋转变化, 根据(18)式重构跟踪框长宽比和尺度w″，h″, 采用线性插值形式训练并更新第t帧滤波器系数 和目标模板x, γ表示学习率, t表示帧数。

(26)

4 实验结果与分析

为了验证GF-KCF目标跟踪算法的有效性和准确性, 选择了7种主流的跟踪算法进行对比:KCF、DSST、SAMF、LCT、DPCF、ECO、CCOT。实验基于实验室的仿真图像序列数据集进行测试, 算法运行PC配置为Inter Xeon E5-1620 CPU@3.60 GHz, RAM为32 GB。实验平台为MATLAB R2018A。

4.1 实验测试条件

本节实验基于战场态势对抗仿真实验平台, 生成红外空战仿真图像数据集进行算法测试, 输出的仿真图像为256×256的16位、单通道灰度图像。该仿真测试数据集的生成条件可以设置包括初始发射条件、目标机动方式、干扰投放策略3个维度的对抗条件参数, 并在限定范围内涵盖了所有对抗条件参数并进行量化, 主要量化参数共包含:1)初始对抗态势共有7个参数:目标高度、载机高度、载机速度、水平进入角、发射距离、综合离轴角; 2)目标机动包含4种类型:左机动、右机动、跃升、俯冲; 3)红外干扰投放策略:总弹数、组数、弹间隔、组间隔。

因此, 本文设置仿真数字平台的仿真图像序列生成条件如下:1)导弹发射距离为7 000 m；2)目标高度及载机高度均为6 000 m；3)目标机动类型为左机动、右机动、跃升、俯冲；4)红外点源干扰弹总数为24枚；5)投弹组数分别为24，12，6；6)组间隔分别为0.8 s，1.5 s，弹间隔为0.1 s；7)水平进入角10°, 45°, 120°, 150°。

4.2 定性分析

本节采用GF-KCF目标跟踪算法进行飞机目标跟踪, 分别选取目标飞机俯冲、跃升、左机动、右机动4个机动类型的红外仿真数据集, 以及红外实测数据集进行测试。部分测试结果如图 6所示。

对于发射距离7 000 m的红外仿真图像序列, 如图 6所示, 目标在序列全过程投放红外点源干扰, 本文算法相比于KCF算法体现出明显的抗部分遮挡效果, 并且对于目标机动有较好的适应性, 更加精确地跟踪目标; 同时, 对于向跃升、俯冲态势的图像序列, 目标同时会发生横向或纵向的旋转变化, 这符合飞机目标的红外成像规律。

4.3 定量分析

跟踪算法的性能好坏通常使用中心位置误差精度(CLE)和重叠精度(OS)来评判。以跟踪算法得到的目标中心位置(x_T, y_T)与人工标注的目标中心位置(x_S, y_S)的欧氏距离来评估跟踪算法精度, 即

(27)

精确度图显示欧氏距离在阈值距离之内的帧数占总帧数的百分比。成功率图显示跟踪目标区域和真实目标区域重合程度, 重叠率公式如(28)式所示, 其中r_t表示跟踪目标区域范围, r_a表示真实目标区域范围。

(28)

本文采用一次性评估方法(one-pass evaluation, OPE)测试本文算法的性能优劣, 即在红外图像数据集上运行跟踪器, 并计算得到精确度和成功率。本文实验测试和统计了107条图像序列的精确度一次性评价和成功率一次性评价, 其中包括96条仿真图像数据和11条外场实测数据。同时, 分别与KCF算法、DSST算法、SAMF算法、LCT算法、DPCF算法、ECO算法、CCOT算法进行对比。测试算法的平均精确度和平均成功率比较分别如图 7、图 8所示。

表 2为每种算法在投放红外点源干扰的仿真数据下的误差阈值为20pixel时的精确度以及跟踪速度结果对比。

由表 2可以看出，本节算法的平均精确度较ECO算法高3.3%，帧率高57.76；与其他算法的精确度相比较，本节算法比CCOT算法高2.2%，比LCT算法高8.9%；与其他算法的帧率比较，本节算法低于KCF算法3 212.26 frame/s，低于DSST算法252.89 frame/s，远高于CCOT算法。

从图 7和图 8中可以看出，本文算法精确度优于其他跟踪算法，成功率仅次于ECO算法，略胜于CCOT算法。其主要原因是本文算法可以稳定跟踪目标仍在视场中的部位区域，这部分区域的分块模板没有被干扰污染，所以平均精确度较高。而由于在图像序列后半部分，即进入末端攻击距离时，目标尺度变化较为剧烈，不能准确的估计干扰投放时目标的尺度变化，因此跟踪框仅集中在目标部分可信区域。

5 结论

本文提出了一种基于GF-KCF的空中红外目标跟踪方法。对于红外运动目标图像序列本算法主要基于二维频域Gabor特征融合进行核相关滤波跟踪。构造二维频域Gabor滤波器组，用于提取、融合目标图像块的频域边缘纹理特征，得到频域响应融合矩阵，有效地抑制背景噪声和突出目标纹理信息。核相关滤波算法对频域响应融合矩阵进行相关滤波和快速检测，得到目标响应矩阵的最大位置即为目标所在位置坐标，将新的目标区域用于训练、更新目标模板，完成运动目标图像序列连续跟踪。