CN107680120B - 基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法 - Google Patents

Abstract

本发明公开了一种基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，首先引入一个表征正常目标区域和表征噪声背景离群区域的遮挡二值指示向量，结合目标字典模板，建立红外小目标的稀疏表示模型；然后，通过显著性检测提取红外图像高频区域，作为粒子滤波模型状态转移过程中的先验信息，在此限制下进行目标状态预测和粒子采样；最后，建立了基于目标样本粒子稀疏重构误差差异性的粒子滤波观测模型，并结合在线模板更新策略实现对目标状态的估计与跟踪。本发明增强了随机粒子的状态估计能力，提升了粒子稀疏表征对小弱运动目标的适应能力和跟踪准确度。

Description

基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法

技术领域

本发明属于红外小目标跟踪技术，特别是一种基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法。

背景技术

红外小目标的捕获跟踪在红外精确制导、红外遥感与航天探测、导航等领域有着广泛应用。目标跟踪一直以来是计算机视觉和图像处理领域的研究热点之一，但在一些应用环境下，成像***需要在距离目标尽可能远的地方对目标进行摄像跟踪，这种情况下的目标在视场中常以小目标的形态出现。由于能量衰减和传感器噪声的影响，红外小目标呈现出亮度弱，面积小，像元数目少，缺乏纹理信息的特点。此外，当背景较为复杂时，红外小目标图像显示出较低的信噪比，目标很容易被噪声和背景杂波所淹没，因而在跟踪过程易发生飘移甚至跟踪失败。

对于红外小目标的跟踪，现有方法主要分为两类：即基于目标建模、定位理论的跟踪方法和基于滤波、数据关联理论的跟踪方法。基于目标建模、定位理论的跟踪方法最典型的是基于均值漂移理论(Meanshift)的跟踪方法，其利用各向同性的核函数对目标图像区域和候选目标图像区域的特征直方图分别进行加权处理,得到“核函数直方图”，然后对处理结果进行相似性度量,最后通过梯度下降法定位目标(1.杨一帆,田雁,杨帆,等.基于改进Mean-Shift方法的红外小目标跟踪,红外与激光工程.43(7)(2014)2164-2169.)。然而，红外小目标亮度弱，像元少，缺乏纹理信息，基于图像灰度空间利用目标直方图特征进行建模不一定能取得很好的效果。基于滤波、数据关联理论的跟踪方法最典型的是基于卡尔曼滤波或粒子滤波的跟踪方法，其基本思想是利用一组带有相关权值的随机样本以及基于这些样本的估算来表示动态样本的后验概率密度，将状态粒子决定的区域所对应的特征与参考模型的特征相比较，得出最佳的后验估计，运用最佳粒子方法确定目标的坐标，从而实现跟踪(2.胡洪涛,敬忠良,胡士强.基于辅助粒子滤波的红外小目标检测前跟踪算法,控制与决策.20(11)(2005)1208-1211.)。在以往的跟踪模型中，粒子滤波多用于大目标的跟踪，并且，对于粒子滤波，传统的粒子转移概率分布一般根据高斯分布建模，粒子数越多，跟踪的准确性越高。但是在实际应用中，红外小目标像素数一般只有几到几十量级，远少于一般目标，所以粒子经高斯分布转移后绝大多数落在非目标区域，获得极低的权值，导致很低的采样效率。因而在跟踪过程中，大量的计算时间被浪费在这些无效粒子上，目标很容易发生漂移并导致跟踪错误。

发明内容

本发明提供一种基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，能够对复杂背景下的红外小目标进行高效并且鲁棒的跟踪。

实现本发明目的的技术解决方案为：一种基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，步骤如下：

首先，输入待跟踪图像序列的首帧，手动框选首帧目标，根据选中的目标初始化稀疏表示模型中的目标字典，引入遮挡项，建立红外小目标的稀疏表示模型；

其次，根据图像首帧目标的选取框的位置、尺度和长宽比信息初始化目标状态参数，从第二帧图像开始，以前一帧目标粒子状态为基准在当前帧红外小目标图像的显著性图像中选取感兴趣的待搜索区域，计算其显著信息累计概率分布图，在该图的纵轴上采样一次，对应横轴位置即可得到一个状态粒子，且显著性越高的区域获得采样粒子的概率越高，则以纵轴上显著性最高处为起点向两侧均匀采样，从横轴上提取每个粒子的位置参数，并根据高斯分布更新其尺度和长宽比参数，即完成显著性先验信息限制下的粒子状态转移过程，采样N_s次得到N_s个状态粒子；

然后，先根据每一个采样粒子的状态参数从仿射空间中提取目标样本，将所有的目标样本代入稀疏表示模型计算稀疏重构系数和稀疏重构残差，并估计用于区分正常目标区域与目标模板离群区域的遮挡二值指示向量，记每个粒子的初始权重为1/N_s，基于稀疏重构残差对目标粒子进行观测并更新每个采样粒子的权重，再根据粒子权重进行重采样并更新目标状态，权重越大的粒子对应的目标状态在粒子滤波观测模型中的概率值也就越大，通过最大后验概率逼近来获得目标最优状态并输出当前时刻的跟踪结果；

接着，先根据稀疏表示模型中的目标稀疏重构系数对目标字典模板分配权重并衡量权重最大的项与当前帧目标的相似性，再根据稀疏表示模型中的遮挡二值指示向量计算当前帧目标中离群区域的占比，即遮挡率，当相似性大于阈值s_th并且离群区域占比小于阈值r_th时，用当前帧的目标模板替换目标字典中权重最小的项，并重新分配、归一化权重，从而对目标字典在线更新；

最后，输入待跟踪图像序列的下一帧，将重采样后的粒子代入状态转移模型以得到新的预测粒子，重复上述步骤完成对红外小目标的跟踪过程。

本发明与现有技术相比，其显著优点为：(1)结合稀疏表示模型对噪声和遮挡干扰不敏感的特点以及粒子滤波框架对非线性非高斯***良好的适应性，在贝叶斯推理框架下建立了基于目标粒子稀疏重构误差差异性的粒子滤波观测模型，能够对复杂背景下的红外小目标进行鲁棒的跟踪。(2)由于红外小目标的成像背景常常多云天空或者海面杂波的干扰，图像信噪比也较低，传统的稀疏表示模型对这类复杂背景下红外小目标的表征不够鲁棒，本发明在目标稀疏表示模型部分，在传统稀疏表示模型基础上引入区分正常目标区域与目标模板离群区域的遮挡二值指示向量，确保了在噪声和遮挡干扰下能对目标样本进行较好的描述，能克服复杂背景对红外小目标的局部遮挡干扰。(3)在粒子转移受限模型部分，根据经典蒙德卡洛法，将显著检测提取的红外图像高频信息转化为粒子滤波采样过程中状态转移的先验信息，从而提高了粒子采样过程中有效粒子的比率，减小跟踪方法的计算负荷，从总体上提升跟踪效率。(4)在目标字典模板更新机制部分，根据目标稀疏表示模型得到的重构系数和重构残差，采用在线方式对目标字典进行更新，能进一步适应背景杂波等因素引起的变化,从而提升跟踪方法的有效性和稳健性。

下面结合附图对本发明作进一步详细叙述。

附图说明

图1是本发明方法的流程图。

图2是红外小目标原图及其显著性图像：(a)红外小目标原图，(b)显著性图像ROI区域。

图3是采样过程示意图：(a)显著概率分布图，(b)显著概率累计分布图，(c)样本粒子示意图。

图4是基于稀疏表示的目标外观模型示意图。

图5是遮挡表征示意图。

图6是红外小目标跟踪结果示意图。

图7是本发明在不同信噪比下对多云小目标的跟踪结果图。

图8是本发明在不同信噪比下对多云小目标的跟踪结果的RMSE曲线图。

图9是本发明对多种复杂背景下红外小目标的跟踪结果图。

具体实施方式

结合图1，本发明基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，步骤如下：

首先，输入待跟踪图像序列的首帧，手动框选首帧目标，根据选中的目标初始化稀疏表示模型中的目标字典，引入遮挡项，建立红外小目标的稀疏表示模型；

其次，根据图像首帧目标的选取框的位置、尺度和长宽比信息初始化目标状态参数，从第二帧图像开始，以前一帧目标粒子状态为基准在当前帧红外小目标图像的显著性图像中选取感兴趣的待搜索区域，计算其显著信息累计概率分布图，在该图的纵轴上采样一次，对应横轴位置即可得到一个状态粒子，且显著性越高的区域获得采样粒子的概率越高，则以纵轴上显著性最高处为起点向两侧均匀采样，从横轴上提取每个粒子的位置参数，并根据高斯分布更新其尺度和长宽比参数，即完成显著性先验信息限制下的粒子状态转移过程，采样N_s次得到N_s个状态粒子；

然后，先根据每一个采样粒子的状态参数从仿射空间中提取目标样本，将所有的目标样本代入稀疏表示模型计算稀疏重构系数和稀疏重构残差，并估计用于区分正常目标区域与目标模板离群区域的遮挡二值指示向量，记每个粒子的初始权重为1/N_s，基于稀疏重构残差对目标粒子进行观测并更新每个采样粒子的权重，再根据粒子权重进行重采样并更新目标状态，权重越大的粒子对应的目标状态在粒子滤波观测模型中的概率值也就越大，通过最大后验概率逼近来获得目标最优状态并输出当前时刻的跟踪结果；

接着，先根据稀疏表示模型中的目标稀疏重构系数对目标字典模板分配权重并衡量权重最大的项与当前帧目标的相似性，再根据稀疏表示模型中的遮挡二值指示向量计算当前帧目标中离群区域的占比，即遮挡率，当相似性大于阈值s_th(根据经验取值，本发明取0.8)并且离群区域占比小于阈值r_th(根据经验取值，本发明取0.5)时，用当前帧的目标模板替换目标字典中权重最小的项，并重新分配、归一化权重，从而对目标字典在线更新；

最后，输入待跟踪图像序列的下一帧，将重采样后的粒子代入状态转移模型以得到新的预测粒子，重复上述步骤完成对红外小目标的跟踪过程。

结合图1，本发明实现上述基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法的具体步骤如下：

1.红外小目标稀疏表示模型的建立过程

在本发明方法中，基于目标模板向量与表征噪声、背景扰动的遮挡向量，利用稀疏表示建立红外小目标的外观模型，具体方法如下。

第一步，基于观测样本构造目标字典：根据图像稀疏编码理论，红外小目标图像子块可以在一组目标特征子空间中线性展开，这组目标特征子空间又可称为目标字典。通常情况下，红外小目标是以相对统一的、规律分布的幅度集中于一个较小区域的亮斑，并且在许多实际跟踪场景中，字典的条目常常是相似的。本发明方法中采用对原始图像采样的方式，直接利用观测样本来构造字典，即根据第一帧图像中目标的初始状态，通过添加一些微小的扰动，基于目标初始位置向上下左右方向移动一两个像素，截取图像中相应的目标子块来初始化字典中的目标原子，构成稀疏表示的目标字典T。这类方法便于结合模板更新策略，较好地抓住图像的结构特征。

第二步，引入遮挡项，建立目标稀疏表示模型：对于候选目标样本y，当不考虑噪声和背景扰动时，其基于目标字典T的稀疏表示方程为公式(1)。

y＝Tα+∈ (1)

其中，α＝(α₁，…，α_n)是重构系数向量，对应n个目标原子(n远小于字典原子的维数)，∈对应重构误差。将向量α看做一个待恢复的稀疏信号，对y的稀疏表示则可以通过一个L1正则化最小函数转化为一个最小误差重构模型为公式(2)，其中，

表示二范数，||·||₁表示一范数，参数λ用于控制重构系数向量α的稀疏度。

然而，在实际情况中红外小目标的运动过程常常面临噪声与背景扰动，比如，在对复杂云天背景下红外小目标进行跟踪研究时发现，当红外小目标在云区与非云区之间穿梭时，目标极易发生漂移甚至丢失。具体来说，在从原始红外图像中截取候选目标图像子块时，难免会引入一些图像噪声和背景，当红外小目标从云区运动到非云区时，背景区域的灰度值会相应发生变化，可以把这种因为图像噪声或背景产生的扰动作为遮挡来处理。那么，在原始稀疏表示模型的基础上，引入遮挡二值指示向量β对正常目标及遮挡区域进行显性区分，如图4所示，那么建立添加遮挡系数后的目标外观模型为公式(3)，即目标由一组目标模板、遮挡项以及残差项线性构成。

y＝Tα+Eβ+∈ (3)

其中，E是一个m阶单位矩阵(m＝w*h，w和h分别表示目标模板的宽高)，对应目标子块的m个像素的位置，β是一个二值列向量，其中非零项对应噪声和背景扰动存在的遮挡区域，并且根据E_i＝1可得Eβ＝β。

通常，将图像子块中的遮挡区域看做是无法纳入目标模板向量线性表示范围的离群区域，并且假设它们在一个相对较小的区域内连续分布，以此作为先验约束β系数。根据图割理论，如果将图像区域表示成一个图G＝(υ，ε)，υ是一组像素点构成的顶点集，ε为像素点的四邻域边界。目标的局部遮挡面积较小且具有空间连续性，这种分布信息可通过马尔可夫随机场MRF很好地表示，如图4所示，把Eβ项对应的目标遮挡区域的这种分布信息用伊辛模型的图像能量函数表示，则最终的目标稀疏表示模型建立为公式(4)：

其中，M_β(·)是一个函数，在计算稀疏重构残差时，排除候选目标子块中遮挡区域，即，对于ε＝y-Tα有

b_i(＞0)是β_i＝1的罚系数，对β进行约束，c_i，j(＞0)用于控制第i个与第j个像素的独立程度。

第三步，根据稀疏表示模型目标样本进行稀疏重构，计算稀疏重构系数、稀疏重构残差和用于区分正常目标区域与离群区域的遮挡二值指示向量。对本发明中目标稀疏表示模型总体求解的难点在于目标函数是非凸的，既包含离散部份又包含连续部分，因而采取一种交替求解的方式，通过对α和β交替估计的方式来估算模型的解。

首先，提取与系数α相关的项，作为一个非负约束下的非线性优化问题，利用加速近邻梯度方法对其进行求解，目标方程可以写成以下形式：

增加的项1_R+(α)定义为

其中，公式(6)式第一项是一个梯度Lipschitz连续的可微凸函数，第二项是一个非光滑的凸函数。只有β较稀疏时，即βE所表征噪声和背景扰动对应的离群区域面积较小时，才能保证该模型有效。因为当离群区域面积越来越大时，会导致

所表示的重建误差较小，使得求解偏离正确值，进而造成错误的跟踪结果。为了保证重构误差相对稳定，用a(β)对目标子块中正常目标区域和遮挡区域的占比进行约束，其中dim(β)表示β的维度，zeros(β)表示β中零项的个数，令a(β)＝min(a_th，dim(β)/zeros(β))。在本发明方法中利用加速梯度投影方法(3.P.Tseng,On accelerated proximal gradient methodsfor convex-concave optimization,Siam J.Optim.(2008).)估计稀疏表征系数α。对α的求解方法的具体迭代过程如方法A所示，根据加速梯度投影方法进行迭代，从而估计稀疏表征系数α，设定迭代终止条件为||α_K-α_*||＜∈，其中α_*是公式(6)的局部极小值，即至多经过

次迭代，输出稀疏系数α的当前最优解。

方法A

然后，提取与系数β相关的项得到式(8)，那么对稀疏表示模型中β的求解是一个组合优化问题，可以写成：

基于M_β(·)函数，对上式进行分解合并(过程见附录)，那么对β的求解可简化为：

argmin表示求令目标函数取最小值时的系数，，在每次迭代中令a(β)等于一个常数，上述能量函数简化为标准的一阶马尔可夫随机场，利用图割得到二值系数β。由于a与系数β是相关的，为了简化求解过程，在每次迭代中令a等于一个常数，那么上式可通过标准的一阶马尔可夫随机场建模，最终利用图论的方法解决。参数b用于控制系数β的稀疏程度，由上式易得，当

时，对应的β_i项更有可能为1，即目标项相应的线性空间剩余量较大时，认为对应区域可能存在噪声或背景扰动，如图5所示。因此参数b应当基于线性空间剩余来取值，本文中令b＝5σ²，σ根据

在线估计，初始化为一个较大的值，在本发明方法中取b＝0.5*[max_{i＝1，…，m}(y_i)]²。

2.显著先验信息限制下的粒子状态转移过程

本发明方法引入频域显著性检测目标状态预测和粒子采样过程进行限制，设计了一种针对红外小目标跟踪的粒子状态转移概率模型，在显著累计概率分布约束下完成粒子状态转移,获取新的目标样本。该模型显著提升了采样粒子中有效粒子的比率，提高了跟踪正确率，而且降低了稀疏表示阶段的方法计算负荷，节省了计算时间。第一步，获取红外小目标图像的显著性图像；第二步，以前一帧目标状态为基准在当前帧显著性图像中提取感兴趣区域；第三步，计算感兴趣区域内的显著概率分布及其累计概率分布；第四步，根据累计概率分布图进行粒子转移，获得新的采样粒子。

第一步，获取红外小目标图像的显著性图像。红外小目标的重要特征是其高频特性，所以常用高通滤波的方法进行小目标自动检测。相应的，对于红外小目标跟踪的粒子状态转移过程，红外图像高频值较大的区域出现目标的概率大于高频值较小的区域，以此为依据，引入傅里叶频谱变换，对t时刻的图像帧提取频域显著区域，其操作过程如下。

首先，利用傅里叶变换获取红外灰度分布图像的频域分布。对于t时刻的原始红外图像f_t(x，y)，其对应的傅里叶频谱为：

其中，W、H是原始红外图像的宽高，j是虚数单位。则对应的傅里叶变换幅度谱为：

A(f_t)(u，v)＝|F(f_t)| (14)

|·|表示求复数的模，在傅里叶频谱图像中，低频部分对应着空间图像中梯度较为平缓的背景区域，高频部分对应灰度变化更加剧烈的目标、噪点、陡峭边缘，由于图像的振幅谱随着频率增加很快衰减，因此可通过式(15)来获取修正后的傅里叶频谱F′(f_t)，从而达到抑制低频、增强高频成分的效果。

F′(f_t)(u，v)＝F(f_t)./A(f_t) (15)

那么利用以下基于傅里叶幅度谱的显著性检测方法，可得到原始红外图像f_t(x，y)的频域显著图像S(f_t)。

其中，

表示二维高斯滤波，

表示卷积操作。F^-1表示反傅里叶变换。

第二步，以前一帧目标状态为基准在当前帧显著性图像中提取感兴趣区域。在粒子滤波方法中，如果将t-1时刻转移粒子状态记为X(t|t-1)，对应的粒子位置状态为(x(t|t-1)，y(t|t-1))，简记为(x，y)。则将基于显著信息的粒子状态转移受限模型建立为：以(x，y)为窗口中心，从显著图像S(f_t)中取出一个(2B+1)×(2A+1)的感兴趣窗口S(y-B：y+B，x-A：x+A)，保证粒子中心不超出窗口范围，如图2所示。

第三步，计算感兴趣区域内的显著概率分布及其累计概率分布，设定阈值TH_s令：

其中，i＝x-A，…，x+A，j＝y-B，…，y+B，假使粒子的初步转移概率在窗口范围内服从均匀分布，根据当前帧的显著信息S′(i，j)，利用P_s′(i，j)＝S′(i，j)/sum(S′)来计算感兴趣区域内的显著概率分布，从而对粒子转移进行约束。然后，将概率分布(pdf)的形式P_s′写成累计概率分布(cdf)的形式P_cum，分别如图3(a)、(b)所示。

第四步，根据累计概率分布图进行粒子转移，获得新的采样粒子，即t时刻的样本粒子位置这样获得：根据经典蒙特卡洛理论，在显著信息累计概率分布图，即图3(b)的纵轴上以显著概率最大处为起点向两侧均匀取样，假如取到了0.36，那么那么纵轴取样处对应到图中横轴上的点就是依据上述概率分布取得的样本点，即在图3(b)的横轴上P_cum＝0.36的点就是依据上述概率分布取得的样本点。最终，t时刻的第k个(Pcum＝0.36)样本粒子的位置可以写为

其中i与j分别是其在图像中的横纵坐标。转移后t时刻该目标粒子的状态参数向量为X_t＝(i，j，l，τ)^T，若目标粒子数量为N_s(N_s基于目标尺寸大小和感兴趣区域大小取值，本发明中取50～200)，则经过N_s次采样得到N_s个目标状态粒子。此外，样本粒子的尺度参数和长宽比参数的状态转移分布依然根据高斯分布p((X_t)_s，r|(X_t-1)_s，r)＝N((X_t)_s，r；(X_t-1)_s，r，∑²)来建模，其中

表示各参数高斯方差，将转移后的尺度和长宽比记为(l，τ)，那么根据t时刻每一个目标粒子的状态参数向量X_t＝(i，j，l，τ)^T可获得对应目标模板的顶点坐标，然后根据仿射变换从t时刻的图像中取出相应的目标样本模板。

3.基于稀疏重构误差差异性的粒子观测和目标状态更新过程

粒子滤波提供了状态估计的总体搜索机制，其通过一组状态空间中随机粒子对概率密度函数进行近似,基于稀疏重构误差差异性对目标粒子进行观测并计算其权重，获得状态最小方差估计，进行目标状态更新，步骤如下：

第一步，根据t时刻每一个目标粒子的状态参数向量X_t获得对应目标模板的顶点坐标，然后根据仿射变换从t时刻的图像中取出相应的目标样本模板，将所有的目标样本代入稀疏表示模型计算稀疏重构系数和稀疏重构残差，并估计用于区分正常目标区域与目标模板离群区域的遮挡二值指示向量，记每个粒子的初始权重为1/N_s，基于稀疏重构残差对目标粒子进行观测并更新每个采样粒子的权重：

首先，定义目标在t时刻的状态变量和观测变量分别为X_t、Y_t。其中，X_t可通过起始点位置参量x、y和尺度以及长宽比参量s、r表示，即X_t＝(t_x，t_y，t_s，t_r)^T，第1帧到第t-1帧图像中目标的观测值表示为Y_1：t-1＝(Y₁，Y₂，…，Y_t-1}，在贝叶斯框架下，对目标信号的状态转移和观测模型建立问题转化对目标运动的状态转移概率和观测概率的计算，则目标状态的预测方程建立为：

p(X_t|Y_1：t-1)＝∫p(X_t|X_t-1)p(X_t-1|Y_1：t-1)dX_t-1 (10)

即可通过截止t-1时刻的概率密度p(X_t-1|Y_1：t-1)来估计t时刻目标状态变量X_t的先验概率密度p(X_t|Y_1：t-1)。由于红外小目标跟踪本质上是从不断涌现出的多个量测点中不断识别和确认目标运动状态。因此先根据(2)中显著先验信息限制下的粒子状态转移模型进行粒子状态预测和采样，计算每个采样粒子对应目标子块的顶点坐标，并根据仿射变换从t时刻的图像中取出相应的目标样本模板。然后根据(1)中的红外小目标稀疏表示模型计算稀疏重构系数和稀疏重构残差，并估计用于区分正常目标区域与目标模板离群区域的遮挡二值指示向量。

其次，根据贝叶斯推论，对目标t时刻的状态估计则转化计算后验概率，则目标状态更新方程建立为：

即可利用截止t时刻目标的观测值Y_1：t和截止t-1时刻的先验概率密度p(X_t|Y_1：t-1)估计其当前时刻的后验概率。p(Y_t|X_t)是一个似然概率，表明观察值为Y_t时状态变量X_t的可能取值，似然模型表明候选目标与目标模板之间的相似性，可直接利用稀疏重构误差建模，即利用目标字典中系数最大的字典对候选目标粒子重构，将重构后的能量残差

作为相似度衡量的标准。假设在第t帧时，所有重构残差、噪声服从独立分布，噪声为高斯白噪声，则在粒子滤波框架下的目标观测模型，某一候选目标粒子的似然函数可建立为

其中，t∈(i)为残差向量中的第i个分量，σ_t∈是高斯方差，取固定值。根据式(12)可知，稀疏重构的残差能量值越小，即候选目标粒子与目标字典中的原子越相似时，候选目标粒子所分配的权重越大；

第二步，根据粒子权重进行重采样并更新目标状态，权重越大的粒子对应的目标状态在粒子滤波观测模型中的概率值也就越大，最终通过最大后验概率逼近X_t＝argmax_xp(X|Y_1：t)来获得t时刻目标最优状态，输出当前帧的跟踪结果。

4.目标字典在线更新过程

在本发明方法中，根据稀疏表示模型的重构系数分配对应目标模板的权重，当遮挡率小于设定阈值且样本与目标字典的相似性大于设定阈值时，用当前帧跟踪目标替换字典中权重最小的项。

红外小目标跟踪的难点在于其成像背景一般为天空或海面，常常伴随云层和海面杂波等环境干扰，红外小目标的外观会随时间发生变化，目标可能被不同程度地遮挡(轻微遮挡或者全部遮挡)，这就导致目标字典中的一些模板不再具有对当前样本较好的表征能力。显然，当目标被局部遮挡且遮挡不严重时应当对目标字典进行更新，而当目标被全部遮挡时则不能更新目标字典。为了保证目标字典原子与实际目标图像子块的一致性，提出一种定量的评价方法，即计算遮挡率以及样本与目标字典的相似性，并建立总体的模板更新机制为:为字典中的目标模板分配权重，考察每个模板对于目标跟踪模型的可靠性，用当前帧的跟踪目标替换目标字典中可靠性较弱的模板，具体步骤如下：

(1)初始化当前字典中每个模板的权值ω_i＝1/n；

(2)y是最新选择的跟踪目标,α,β是式(8)的解，则更新权重为ω_i←ω_iexp(α_i)，i＝1，…，n；

(3)筛选目标字典中权重最大的模板T_m,其中，max＝max_1≤i≤nα_i；

(4)计算跟踪结果y与目标字典原子T_m的相似度s＝arccos<y，T_m>，其中arccos计算余弦相似度；

(5)根据遮挡二值指示向量β计算遮挡率r＝ones(β)/size(β)；

(6)当满足s＞s_th且r＜r_th时，利用当前的跟踪目标来更新字典中权重最小的目标原子，即令T_k＝y，其中k＝argmin_1≤i≤nω_i，s_th与r_th是预先设置的阈值；

(7)对目标字典第k项重新分配权重令ω_k＝median(ω)，标准化权重向量ω，令sum(ω)＝1；对ω设置阈值令max(ω)＜ω_th防止斜交；

(8)标准化更新的目标字典T。

5.输入待跟踪图像序列的下一帧，将重采样后的粒子代入状态转移模型以得到新的预测粒子，重复权利要求1至权利要求5的步骤完成对红外小目标的跟踪过程。

本发明可通过以下实施例进一步说明。

为说明本发明方法在红外小目标跟踪中的有效性和鲁棒性，对六组不同背景下的红外小目标图像序列进行了测试，并对基于加速梯度逼近的L1跟踪方法(L1APG)，基于连续异常概率模型的视觉跟踪方法(PCOM)，基于增量学习的鲁棒性视觉跟踪方法(IVT)，基于低秩稀疏表示的红外小目标跟踪方法(RLRST)以及本发明方法(Ours)的跟踪结果进行比较，分别对跟踪准确率(通过均方根误差衡量)和跟踪成功率两项指标进行统计计算，分析本发明方法与其他方法相比存在的优势。

测试数据集1中的红外小目标以抛物线轨迹运动，拍摄于噪声干扰的多云天空背景下，其详细信息见表1。这组图像可以等分为三组子图像序列，对应三种不同的信噪比，用于测试本发明方法在不同信噪比下的鲁棒性。在这组试验中，粒子数设置为50，目标字典中模板数设置为10，采样过程中感兴趣区域大小设置为9×9，根据无优化的MATLAB代码仿真结果，本发明方法的平均计算消耗时间为0.0584秒每帧。

表1测试数据集详细信息

对测试数据集1中小目标的跟踪结果如图7所示。对于每一帧图像，将帧号打印在图像的左上角，将目标字典子空间和遮挡二值指示向量对应的遮挡表征模板显示在图片下方，跟踪目标用白色框标出。实验处理结果显示，第176帧中的云区遮挡和噪点对应位置的值都是1，这说明当存在遮挡干扰时，本发明方法时能够有效地将目标模板中的正常目标区域提取出来用于稀疏表示。此外，根据显示在每幅图像下方的目标子空间，可以看到信噪比SNR＝3时目标字典中的模板更新频率明显小于信噪比SNR＝1时目标字典中的模板更新频率。实际上，对于SNR＝6的第一组子图像序列，尽管其信噪比很高，在基于连续异常概率模型的视觉跟踪方法和基于增量学习的鲁棒性视觉跟踪方法的试验中，当目标从非云区运动到云区时，由于图像背景的灰度特征发生突变，目标也迅速产生了漂移。与此形成对比的是，本发明方法在三组信噪比下都保持了比较稳定的跟踪。也就是说，本发明方法在噪声较强或者背景发生突变的情况下依然具有很强的鲁棒性，这在一定程度上印证了方法的抗噪声和抗遮挡干扰能力。

为衡量本发明方法的跟踪准确率，引入均方根误差(RMSE)作为定量标准，RMSE的值越小，跟踪准确性越高，其定义为

其中，

和

分别表示t时刻目标位置的真值和本发明方法的估计值，M表示蒙德卡洛实验的次数。在本次试验中，M设为100。图8是测试数据集1的RMSE曲线图。其中，x轴表示时间，y轴对应着某时刻下不同跟踪方法的均方根误差值，x轴的第1-80、81-160、161-240分别对应着测试数据集1的三组信噪比分别为SNR＝6、SNR＝3、SNR＝1的图像子序列。在图8中，基于加速梯度逼近的L1跟踪方法，基于连续异常概率模型的视觉跟踪方法，基于增量学习的鲁棒性视觉跟踪方法，基于低秩稀疏表示的红外小目标跟踪方法以及本发明方法的RMSE曲线分别用空心上三角、空心五边形、实心下三角，实心圆，空心圆标记。需要说明的是，在这些远离x轴的标记点所对应时刻的图像帧中，目标均发生不同程度的漂移，并且，偏离x轴越远，跟踪误差越大。

所以，根据图8可以得出以下结论，当信噪比较高时(SNR＝6或SNR＝3)，基于加速梯度逼近的L1跟踪方法、基于低秩稀疏表示的红外小目标跟踪方法以及本发明方法的RMSE曲线均比较贴近x轴，即都能够保持比较稳定的跟踪。然而，当信噪比降低到SNR＝1时，基于低秩稀疏表示的红外小目标跟踪方法和基于加速梯度逼近的L1跟踪方法的RMSE曲线上远离x轴的标记点迅速增多，也就是说，跟踪结果中偏离真值的帧数迅速增加，这两种方法的跟踪性能迅速恶化，相反的，本发明方法的RMSE曲线依然贴近x轴，即依然保持了鲁棒的跟踪。此外，基于连续异常概率模型的视觉跟踪方法和基于增量学习的鲁棒性视觉跟踪方法在三组信噪比不同的子图像序列中，RMSE曲线均迅速偏离x轴，这两种方法对于多云背景干扰下的红外小目标跟踪完全失败。

测试数据集1-6用于测试本发明方法在不同成像背景下对红外小目标的跟踪结果，图7给出了测试数据集1的三种信噪比下的跟踪结果图，图9给出了测试数据集2-6的五种背景下的跟踪结果图，表2展示了本发明方法以及几种广泛应用的跟踪方法对这六组测试数据集的跟踪成功率。其中，测试数据集2和序列3拍摄于晴朗天空背景下，测试数据集4和测试数据集6分别拍摄于噪声和遮挡干扰的多云天空和海天交界背景下，测试数据集5中的目标由极远处向成像***靠近，因此红外小目标尺寸缓慢增大。这五组测试数据集的详细信息均在表1中给出。为了定量地衡量本发明方法及几种比较方法在多种成像背景下对红外小目标的跟踪能力，利用跟踪正确率(TSR)作为标准对各方法进行了评估。具体来说，当测试图像帧中目标跟踪结果的均方根误差小于某一阈值时，认为该帧的跟踪结果是可接受的，那么跟踪正确率定义为：

其中，NT表示每一组测试数据集中存在目标的帧数，

表示每一组测试数据集的目标跟踪结果中可接受的帧数，本实验中D₀取3。那么，对于表1中的6组测试数据集，上述几种方法跟踪正确率如表2所示。

表2几种方法的跟踪正确率比较

根据表2，对于所有测试数据集，本发明方法保持了95％以上的跟踪正确率，而基于连续异常概率模型的视觉跟踪方法和基于增量学习的鲁棒性视觉跟踪方法在测试数据集4-6中的跟踪正确率均低于30％。此外，对于测试数据集1的第三组图像子序列，本发明方法的跟踪正确率是基于低秩稀疏表示的红外小目标跟踪方法的6倍，是基于加速梯度逼近的L1跟踪方法的3.41倍，也就是说本发明方法在低信噪比条件下对于红外小目标的跟踪鲁棒性相对其他方法具有明显优势。并且，根据上述数据，当存在遮挡和背景杂波干扰时，本发明方法在多种成像背景下总体具有最优的表现。

另外，根据本发明所建立的目标外观模型，结合M_β(·)函数，在计算候选目标粒子y在目标字典T上的稀疏重构残差时，需排除候选目标子块中无法纳入目标模板线性表示范围离群区域，即对于目标子块中的某项，若判定其被遮挡，则相应有

等于零，该过程通过将目标模板的每一项乘以一个系数(1-β_i)来实现，多项式分解合并过程如下：

Claims (5)

1.一种基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，其特征在于步骤如下：

首先，输入待跟踪图像序列的首帧，手动框选首帧目标，根据选中的目标初始化稀疏表示模型中的目标字典，引入遮挡项，建立红外小目标的稀疏表示模型；

其次，根据图像首帧目标的选取框的位置、尺度和长宽比信息初始化目标状态参数，从第二帧图像开始，以前一帧目标粒子状态为基准在当前帧红外小目标图像的显著性图像中选取感兴趣的待搜索区域，计算其显著信息累计概率分布图，在该图的纵轴上采样一次，对应横轴位置即可得到一个状态粒子，且显著性越高的区域获得采样粒子的概率越高，则以纵轴上显著性最高处为起点向两侧均匀采样，从横轴上提取每个粒子的位置参数，并根据高斯分布更新其尺度和长宽比参数，即完成显著性先验信息限制下的粒子状态转移过程，采样N_s次得到N_s个状态粒子；

然后，先根据每一个采样粒子的状态参数从仿射空间中提取目标样本，将所有的目标样本代入稀疏表示模型计算稀疏重构系数和稀疏重构残差，并估计用于区分正常目标区域与目标模板离群区域的遮挡二值指示向量，记每个粒子的初始权重为1/N_s，基于稀疏重构残差对目标粒子进行观测并更新每个采样粒子的权重，再根据粒子权重进行重采样并更新目标状态，权重越大的粒子对应的目标状态在粒子滤波观测模型中的概率值也就越大，通过最大后验概率逼近来获得目标最优状态并输出当前时刻的跟踪结果；

所述遮挡二值指示向量是指当红外小目标从云区运动到非云区时，背景区域的灰度值会相应发生变化，把这种因为图像噪声或背景产生的扰动作为遮挡来处理，在原始稀疏表示模型的基础上，引入遮挡二值指示向量β对正常目标及遮挡区域进行显性区分；

接着，先根据稀疏表示模型中的目标稀疏重构系数对目标字典模板分配权重并衡量权重最大的项与当前帧目标的相似性，再根据稀疏表示模型中的遮挡二值指示向量计算当前帧目标中离群区域的占比，即遮挡率，当相似性大于阈值s_th并且离群区域占比小于阈值r_th时，用当前帧的目标模板替换目标字典中权重最小的项，并重新分配、归一化权重，从而对目标字典在线更新；

最后，输入待跟踪图像序列的下一帧，将重采样后的粒子代入状态转移模型以得到新的预测粒子，重复上述步骤完成对红外小目标的跟踪过程。

2.根据权利要求1所述的基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，其特征在于红外小目标稀疏表示模型的建立过程如下：

第一步，基于观测样本构造目标字典：采用对原始图像采样的方式，直接利用观测样本来构造字典，即根据第一帧图像中目标的初始状态，通过添加一些微小的扰动，基于目标初始位置向上下左右方向移动一两个像素，截取图像中相应的目标子块来初始化字典中的目标原子，构成稀疏表示的目标字典T；

第二步，引入遮挡项，建立目标稀疏表示模型：对于候选目标样本y，当不考虑噪声和背景扰动时，其基于目标字典T的稀疏表示方程如下：

y＝Tα+∈ (1)

其中，α＝(α₁，…，α_n)是重构系数向量，对应n个目标原子，n远小于字典原子的维数，∈对应重构误差，将向量α看做一个待恢复的稀疏信号，对y的稀疏表示则通过一个L1正则化最小函数转化为一个最小误差重构模型为公式(2)，其中，

表示二范数，||·||₁表示一范数，参数λ用于控制重构系数向量α的稀疏度；

建立添加遮挡系数后的目标外观模型为公式(3)，即目标由一组目标模板、遮挡项以及残差项线性构成：

y＝Tα+Eβ+∈ (3)

其中，E是一个m阶单位矩阵，m＝w*h，w和h分别表示目标模板的宽高，对应目标子块的m个像素的位置，β是一个二值列向量，其中非零项对应噪声和背景扰动存在的遮挡区域，并且根据E_i＝1可得Eβ＝β；

根据图割理论，图像区域表示成一个图G＝(υ，ε)，υ是一组像素点构成的顶点集，ε为像素点的四邻域边界，目标的局部遮挡面积较小且具有空间连续性，这种分布信息通过马尔可夫随机场MRF很好地表示，把Eβ项对应的目标遮挡区域的这种分布信息用伊辛模型的图像能量函数表示，则最终的目标稀疏表示模型建立为公式(4)：

其中，M_β(·)是一个函数，在计算稀疏重构残差时，排除候选目标子块中遮挡区域，即，对于ε＝y-Tα有

b_i(＞0)是β_i＝1的罚系数，对β进行约束，c_i，j(＞0)用于控制第i个与第j个像素的独立程度；

第三步，根据稀疏表示模型目标样本进行稀疏重构，计算稀疏重构系数、稀疏重构残差和用于区分正常目标区域与离群区域的遮挡二值指示向量，即采取一种交替求解的方式，通过对α和β交替估计的方式来估算模型的解：

首先，提取与系数α相关的项，作为一个非负约束下的非线性优化问题，利用加速近邻梯度方法对其进行求解，目标方程写成以下形式：

增加的项1_R+(α)定义为

其中，公式(6)第一项是一个梯度Lipschitz连续的可微凸函数，第二项是一个非光滑的凸函数；

为了保证重构误差相对稳定，用a(β)对目标子块中正常目标区域和遮挡区域的占比进行约束，其中dim(β)表示β的维度，zeros(β)表示β中零项的个数，令a(β)＝min(a_th，dim(β)/zeros(β))；根据加速梯度投影方法进行迭代，从而估计稀疏表征系数α，设定迭代终止条件为|α_K-α_*||＜∈，其中α_*是公式(6)的局部极小值，即至多经过

次迭代，输出稀疏系数α的当前最优解；

然后，提取与系数β相关的项得到公式(8)，那么对稀疏表示模型中β的求解如下：

基于M_β(·)函数，对上式进行分解合并，那么对β的求解可简化为：

argmin表示求令目标函数取最小值时的系数，在每次迭代中令a(β)等于一个常数，上述能量函数简化为标准的一阶马尔可夫随机场，利用图割得到二值系数β。

3.根据权利要求1所述的基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，其特征在于显著先验信息限制下的粒子状态转移过程如下：

第一步，获取红外小目标图像的显著性图像，即对t时刻的图像帧提取频域显著区域，其操作过程如下：

首先，利用傅里叶变换获取红外灰度分布图像的频域分布，对于t时刻的原始红外图像f_L(x，y)，其对应的傅里叶频谱为：

其中，W、H是原始红外图像的宽高，j是虚数单位，则对应的傅里叶变换幅度谱为：

Λ(f_t)(u，v)＝|F(f_t)| (14)

|·|表示求复数的模，在傅里叶频谱图像中，低频部分对应着空间图像中梯度较为平缓的背景区域，高频部分对应灰度变化更加剧烈的目标、噪点、陡峭边缘，由于图像的振幅谱随着频率增加很快衰减，因此通过公式(15)来获取修正后的傅里叶频谱F′(f_t)，从而达到抑制低频、增强高频成分的效果：

F′(f_t)(u，v)＝F(f_t)./Λ(f_t) (15)

那么利用以下基于傅里叶幅度谱的显著性检测得到原始红外图像f_L(x，y)的频域显著图像S(f_L)：

其中，

表示二维高斯滤波，

表示卷积操作，F^-1表示反傅里叶变换；

第二步，以前一帧目标状态为基准在当前帧显著性图像中提取感兴趣区域，即在粒子滤波方法中，如果将t-1时刻转移粒子状态记为X(t|t-1)，对应的粒子位置状态为(x(ι|ι-1)，y(ι|ι-1))，简记为(x，y)，则将基于显著信息的粒子状态转移受限模型建立为：以(x，y)为窗口中心，从显著图像S(f_L)中取出一个(2B+1)×(2A+1)的感兴趣窗口S(y-B：y+B，x-A：x+A)，保证粒子中心不超出窗口范围；

第三步，计算感兴趣区域内的显著概率分布及其累计概率分布，设定阈值TH_s令：

其中，i＝x-A，…，x+A，j＝y-B，…，y+B，假使粒子的初步转移概率在窗口范围内服从均匀分布，根据当前帧的显著信息S′(i，j)，利用P_s′(i，j)＝S′(i，j)/sum(S′)来计算感兴趣区域内的显著概率分布，从而对粒子转移进行约束，然后，将概率分布(pdf)的形式P_s′写成累计概率分布(cdf)的形式P_cum；

第四步，根据累计概率分布图进行粒子转移，获得新的采样粒子，即t时刻的样本粒子位置这样获得：根据经典蒙特卡洛理论，在显著信息累计概率分布图的纵轴上以显著概率最大处为起点向两侧均匀取样，那么纵轴取样处对应到图中横轴上的点就是依据上述概率分布取得的样本点；最终，t时刻的第k个样本粒子的位置写为

其中i与j分别是其在图像中的横纵坐标；同时，样本粒子的尺度参数和长宽比参数的状态转移分布依然根据高斯分布p((X_t)_s，r|(X_t-1)_s，r)＝N((X_t)_s，r；(X_t-1)_s，r，∑²)来建模，其中

表示各参数高斯方差，将转移后的尺度和长宽比记为(l，τ)，结合之前的位置参数(i，j)，转移后t时刻该目标粒子的状态参数向量为X_t＝(i，j，l，τ)^T，若目标粒子数量为N_s，则经过N_s次采样得到N_s个目标状态粒子。

4.根据权利要求1所述的基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，其特征在于基于稀疏重构误差差异性的粒子观测和目标状态更新过程如下：

第一步，根据t时刻每一个目标粒子的状态参数向量X_t可获得对应目标模板的顶点坐标，然后根据仿射变换从t时刻的图像中取出相应的目标样本模板，将所有的目标样本代入稀疏表示模型计算稀疏重构系数和稀疏重构残差，并估计用于区分正常目标区域与目标模板离群区域的遮挡二值指示向量，记每个粒子的初始权重为1/N_s，基于稀疏重构残差对目标粒子进行观测并更新每个采样粒子的权重：

(1)定义目标在t时刻的状态变量和观测变量分别为X_t、Y_t，其中，X_t可通过起始点位置参量x、y和尺度以及长宽比参量s、r表示，即X_t＝(t_x，t_y，t_s，t_r)^T，第1帧到第t-1帧图像中目标的观测值表示为Y_1：t-1＝{Y₁，Y₂，…，Y_t-1}，在贝叶斯框架下，对目标信号的状态转移和观测模型建立问题转化对目标运动的状态转移概率和观测概率的计算，则目标状态的预测方程建立为：

p(X_t|Y_1：t-1)＝∫p(X_t|X_t-1)p(X_t-1|Y_1：t-1)dX_t-1 (10)

即可通过截止t-1时刻的概率密度p(X_t-1|Y_1：t-1)来估计t时刻目标状态变量X_t的先验概率密度p(X_t|Y_1：t-1)；

(2)根据贝叶斯推论，对目标t时刻的状态估计则转化计算后验概率，则目标状态更新方程建立为：

即可利用截止t时刻目标的观测值Y_1：t和截止t-1时刻的先验概率密度p(X_t|Y_1：t-1)估计其当前时刻的后验概率；p(Y_t|X_t)是一个似然概率，表明观察值为Y_t时状态变量X_t的可能取值，似然模型表明候选目标与目标模板之间的相似性，可直接利用稀疏重构误差建模，即利用目标字典中系数最大的字典对候选目标粒子重构，将重构后的能量残差

作为相似度衡量的标准；假设在第t帧时，所有重构残差、噪声服从独立分布，噪声为高斯白噪声，则在粒子滤波框架下的目标观测模型，某一候选目标粒子的似然函数可建立为

其中，tc(i)为残差向量中的第i个分量，σ_t∈是高斯方差，取固定值；根据公式(12)可知，稀疏重构的残差能量值越小，即候选目标粒子与目标字典中的原子越相似时，候选目标粒子所分配的权重越大；

第二步，根据粒子权重进行重采样并更新目标状态，权重越大的粒子对应的目标状态在粒子滤波观测模型中的概率值也就越大；最终，通过最大后验概率逼近X_t＝argmax_xp(X|Y_1：t)来获得t时刻目标最优状态，输出当前帧的跟踪结果。

5.根据权利要求1所述的基于稀疏表示和转移受限粒子滤波的红外小目标跟踪方法，其特征在于目标字典在线更新过程如下：

(1)初始化当前字典中每个模板的权值ω_i＝1/n；

(2)y是最新选择的跟踪目标，α，β是公式(4)的解，则更新权重为

ω_i←ω_iexp(α_i)，i＝1，…，n；

(3)筛选目标字典中权重最大的模板T_m，其中，max＝max_1≤i≤nα_i；

(4)计算跟踪结果y与目标字典原子T_m的相似度s＝arccos<y，T_m>，其中

arccos计算余弦相似度；

(5)根据遮挡二值指示向量β计算遮挡率r＝ones(β)/size(β)；

(6)当满足s＞s_th且r＜r_th时，利用当前的跟踪目标来更新字典中权重最小的目标原子，即令T_k＝y，其中k＝argmin_1≤i≤nω_i，s_th与r_th是预先设置的阈值；

(7)对目标字典第k项重新分配权重令ω_k＝median(ω)，标准化权重向量ω，令sum(ω)＝1；对ω设置阈值令max(ω)＜ω_th防止斜交；

(8)标准化更新的目标字典T。

Images