本发明涉及一种视频多目标跟踪方法，具体是一种基于高斯混合概率假设密度(GM-PHD)滤波的动态视频多目标跟踪方法。
背景技术：
多目标跟踪(MTT)技术已在智能监控，机器人自动导航，生物学等不同领域得到广泛应用。MTT的目的是从包含有杂波的量测集合中估计出多目标的数目和状态。由于量测中存在杂波以及检测的不确定性等原因，精确跟踪数目变化的多目标仍然是一个难题。传统的多目标跟踪应用“量测-航迹”数据关联技术，旨在把多目标跟踪问题解耦为多个单目标的跟踪问题。Reid等提出的多假设方法(MHT)和Bar-Shalom等提出的联合概率数据关联方法(JPDA)是最具代表性和有效性的两种关联方法。由于方法的组合特性，数据关联方法的共同缺点是计算量庞大，因而限制了该类算法的实际应用。
近年来，基于随机有限集(RFS)的多目标跟踪方法取得了巨大突破，并得到了广泛的关注。基于RFS理论框架方法将不确定的量测与多目标状态表示成随机集合，进而将多目标跟踪问题表述为多维贝叶斯滤波器，从而避免了“量测-航迹”的数据关联。其中，Mahler提出的概率假设密度(PHD)滤波及后来学者提出的PHD实现方法最具代表性。PHD使用多目标随机集的后验概率密度的一阶统计量(概率假设密度，即PHD)近似代替多目标的后验概率密度，简化多维贝叶斯滤波递推公式中的积分运算，使得多维贝叶斯滤波的实现成为可能。但是PHD滤波器存在目标数估计不准确的问题。为了解决这个问题，Mahler又提出了基数概率假设密度(CardinalizedProbabilityHypothesisDensity，CPHD)滤波器。CPHD通过同时传播状态随机集的强度函数和基数分布得到准确的目标数估计。PHD和CPHD滤波器虽然能够处理多目标的新生、孵化和死亡问题，但仍存在一个缺点，即标准PHD和CPHD滤波器假定新生目标的强度函数为已知。然而，一般情况下新生目标的强度函数并不容易获得。因此，该假设制约了PHD滤波器的广泛应用。
经对现有技术的文献检索发现，B.RISTIC等人在2012年的IEEETransactionsonAerospaceandElectronicSystems(国际电气和电子工程师协会航空电子系统学报)第48卷第2期上发表的“AdaptivetargetbirthintensityforPHDandCPHDfilters(PHD和CPHD滤波器中的新生目标自适应强度函数)”中提出一种基于粒子PHD滤波器的自适应新生目标强度方法。该方法由各帧量测设计新生目标强度函数，把量测空间分为量测可测子空间p与量测不可测子空间v，对于p由似然函数构造新生目标强度，但对于v仍然需要已知先验信息，没有根本解决新生目标强度函数的构造问题。MichaelBeard等人提出新生目标强度分布为均匀分布时的高斯混合PHD滤波器，但仍然假设新生目标强度分布为已知的。欧阳诚等人针对文献，提出一种归一化因子修正方法，用以改进原算法的航迹归一化失衡问题，但该方法仍需已知部分先验信息，而且该方法在滤波后仍然包含大量杂波，弱化了PHD滤波器去除杂波的主要优点。Wang等人提出使用全部量测的概率似然比测试(SPRT)方法进行新生目标的检测，进而构造出生目标强度函数。但是由于其使用全部量测的组合优化，不可避免的带来巨大的计算负担。
技术实现要素：
本发明针对上述现有技术中存在的不足，提供一种基于迭代随机采样一致性(RANSAC)算法的自适应新生目标强度估计的GM-PHD多目标跟踪方法。本发明提出一种可以从量测集中估计出新生目标位置的迭代RANSAC(I-RANSAC)算法，基于估计出的新生目标位置信息，构建出新生目标的强度函数，从而解决了原始PHD滤波器需要已知目标强度函数的不足。本发明给出了基于提出的I-RANSAC算法的自适应GM-PHD多目标跟踪框架，解决了原始GM-PHD多目标跟踪方法需要已知新生目标强度函数的问题，具有实现简单和鲁棒性强的优点。
本发明是通过以下技术方案实现的：
本发明包括：I-RANSAC新生目标检测模块和PHD滤波模块，其中：I-RANSAC新生目标检测模块与“滑窗”中的量测数据和当前的状态估计相连传输当前滑窗中可能出现的新生目标的位置信息，所谓“滑窗”是指固定帧数的连续的量测集合，该集合随着时间的推移而向前移动，即当新的一帧量测到来后，把该帧量测加入“滑窗”集合的最后面，而舍弃“滑窗”集合中最前面的一帧量测，从而保持“滑窗”集合量测帧数的固定。确切的说，“滑窗”是一个随着时间向前移动的量测队列。PHD滤波模块与I-RANSAC新生目标检测模块和当前量测相连传输全部目标的状态估计随机集和目标数估计随机集。
所述I-RANSAC新生目标检测模块包括：有效量测生成子模块、假设生成子模块、假设检验子模块，其中：有效量测生成子模块分别与当前量测和当前状态估计相连传输当前“滑窗”的有效量测，假设生成子模块与有效量测子模块相连传输从有效量测中通过随机采样而生成的航迹假设，假设检验子模块与假设生成子模块相连传输对当前假设航迹的经验结果，即确认的新生目标的位置信息。
所述PHD滤波模块包括：新生目标强度函数构建子模块、预测子模块、更新子模块、高斯元修剪子模块和状态抽取子模块，其中：新生目标强度函数构建子模块与I-RANSAC新生目标检测模块中的假设检验子模块相连传输新生目标强度函数，预测子模块与新生目标强度函数构建子模块相连传输目标状态随机集PHD的预测高斯元参数，更新子模块分别与预测子模块和量测相连传输目标状态随机集PHD的更新后的高斯元参数，高斯元修剪子模块与更新子模块相连传输对更新后的高斯元合并和删除后的目标状态随机集PHD的高斯元参数，状态抽取子模块与高斯元修剪子模块相连传输目标状态估计随机集和目标数估计随机集。
与现有技术相比，本发明有益效果是：基于I-RANSAC的新生目标检测模块，为PHD滤波器提供了新生目标强度信息；解决了PHD滤波器需要知道新生目标强度函数的问题，保证了PHD滤波器的鲁棒性和可靠性。该方法简单有效、易于实施，鲁棒性好，可广泛应用于军事与民用多目标跟踪领域。
附图说明
图1为本发明实施例中基于I-RANSAC的自适应GM-PHD多目标跟踪系统框图；
图2为本发明实施例中跟踪场景及目标航迹图；
图3为本发明实施例中每一帧x方向和y方向的包括杂波的量测图；
图4为本发明实施例中在x和y方向上对时间的位置估计结果；
图5为本发明实施例给出了由标准GM-PHD滤波器估计同一多目标场景的结果；
图6为本发明实施例中给出100次蒙特卡洛平均的真实目标数、标准GM-PHD滤波器和本发明方法估计的目标数的比较结果，；
图7为本发明实施例中了100次蒙特卡洛平均的标准GM-PHD滤波器和本发明方法估计的OSPA距离的比较图；
具体实施方式
下面结合附图对本发明的实施例作详细说明：本实施例在以本发明技术方案为前提下进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。
如图1所示，本实施例包括：I-RANSAC新生目标检测模块和PHD滤波模块，其中：I-RANSAC新生目标检测模块包括：有效量测生成子模块、假设生成子模块、假设检验子模块，其中：有效量测生成子模块分别与当前量测和当前状态估计相连传输当前“滑窗”的有效量测，假设生成子模块与有效量测子模块相连传输从有效量测中通过随机采样而生成的航迹假设，假设检验子模块与假设生成子模块相连传输对当前假设航迹的经验结果，即确认的新生目标的位置信息。
所述有效量测生成子模块实现用于检测新生目标的有效量测。
所述假设生成子模块实现各种可能的新生目标假设。
所述假设检验子模块实现新生目标假设的检验，以确认可能的新生目标。
PHD滤波模块包括：新生目标强度函数构建子模块、预测子模块、更新子模块、高斯元修剪子模块和状态抽取子模块。其中：新生目标强度函数构建子模块与I-RANSAC新生目标检测模块中的假设检验子模块相连传输新生目标强度函数，预测子模块与新生目标强度函数构建子模块相连传输目标状态随机集PHD的预测高斯元参数，更新子模块分别与预测子模块和量测相连传输目标状态随机集PHD的更新后的高斯元参数，高斯元修剪子模块与更新子模块相连传输对更新后的高斯元合并和删除后的目标状态随机集PHD的高斯元参数，状态抽取子模块与高斯元修剪子模块相连传输目标状态估计随机集和目标数估计随机集。
所述新生目标强度函数子模块实现新生目标的强度函数的构建。
所述预测子模块实现目标状态随机集的PHD的高斯元参数的预测。
所述更新子模块实现更新目标预测状态随机集PHD的高斯元参数。
所述高斯元修剪子模块实现合并距离很近的高斯元，以及去除极小权值的高斯元，用于减少计算量和杂波。
所述状态抽取子模块抽取权值大于阈值(一般取0.5)的高斯元所对应的期望值即为状态随机集，高斯元数为目标数随机集。
本实施例中效量测生成子模块的具体工作过程为：
1)输入新帧量测，放入“滑窗”队列的末尾，去掉“滑窗”队列的首帧量测，实现“滑窗”队列的移动。
2)设当前状态估计为xk，计算“滑窗”中对应帧的量测残差dk(i)＝zk(i)-HkFkxk，其中Hk，Fk分别是系统的量测矩阵和状态转移矩阵；以及残差的范数：式中Sk是残差向量的协方差矩阵。若gk(i)距离小于设定的门限，则对应的量测被认为是与存在航迹关联的，把该量测去掉。如此处理则得到有效的量测集合。
本实施例中假设生成子模块的具体工作过程为：从有效量测集的前面两帧量测中各随机抽取包含一个样本的子集，根据这两个样本点由最小二乘法计算出其模型M；有效量测集中除去抽取的两个样本后的集合叫做原量测集的余集。余集中与模型M的误差小于某一设定阈值的样本叫做内点，它们构成集合叫做一致集。模型M即当前的假设模型。
本实施例中假设检验子模块的具体工作过程为：若一致集的基数大于给定的阈值，则认为得到了正确的模型参数，并利用一致集中的样本由最小二乘法重新计算模型，该新模型即是一个新生的目标模型；确定一个新目标后，由余集替换有效量测集，重复上述假设生成和假设检验过程，即可确认多个新生目标。
本实施例中目标模型为：
1)目标的状态方程为非线性渐近转弯(CT)模型：
xk＝F(ω)xk-1+Gvk-1(3)
其中，ω是转弯速率，当以目标质心坐标及其速度来描述目标时，状态xk可表示为：
xk＝(locx，k，locy，k，velx，k，vely，k，ω)T(4)
F(ω)为状态转移矩阵，
状态噪声vk是以Q为协方差的零均值高斯噪声，σv为系统噪声标准差。
2)系统观测模型为zk＝Hxk+wk，其中，观测矩阵H为：
观测噪声wk是零均值高斯噪声，协方差为σw为观测噪声标准差。
本实施例PHD滤波模块的一些参数设置如下：取每一个目标的生存概率PS＝0.99，目标检测概率PD＝0.99。状态噪声标准差σv＝2m/s2，观测噪声标准差σw＝8。采样时间间隔取为T＝1s。杂波在跟踪区域[0，20000]×[0，20000]内服从均匀分布，而每一时刻杂波的数目服从参数为λc＝10-7m-2的泊松分布(即平均每一帧有40个杂波)。
本实施例中新生目标强度函数子模块的具体工作过程为：
假设k时刻由I-RANSAC新生目标检测模块检测出m个新生目标，则以每一个新生目标位置作为位置均值，以给定值作为速度均值，并以给定的较大的正定数值矩阵为方差阵构造出一个多维高斯分布函数作为每一个新生目标的强度函数，而其权值取为0.03.
本实施例中预测子模块的具体工作过程为：
假设k时刻描述PHD的高斯元参数为其中，Jk为k时刻的高斯元的数目，为k时刻第i个高斯元的均值，是其权值，是相应的协方差。
1)对新产生目标进行预测，设k+1时刻新产生目标个数为Jγ，k+1，则对j＝1，…，Jγ，k+1有其中分别为新生目标高斯元的权值，状态期望(均值)及协方差；
2)对孵化目标预测，设孵化目标个数为Jβ，k+1，则对j＝1，…，Jβ，k+1，l＝1，…，Jk有其中为孵化目标高斯元模型的权值，分别为(5)式及(6)式中的状态转移矩阵及状态噪声协方差。
3)对继续存在目标进行预测计算，设其生存概率为pS，则对j＝1，…，Jk，按照以下公式更新权值，均值及协方差：
本实施例中更新子模块的具体工作过程为：
1)利用运动目标检测模块得到的量测随机集记为Zk+1，及(7)式的观测矩阵H和量测噪声协方差R更新权值，均值及协方差。pD为检测概率如上实施例中所述，则对未检测到的目标用式(8)～(10)更新：对j＝1，…，Jk+1｜k，其中Jk+1｜k＝Jγ，k+1+lJβ，k+1+Jk+1，有：
2)对检测到的目标更新，即利用运动目标检测模块的质心坐标作为量测随机集进行更新计算，对每一个z∈Zk+1，计算：
设泊松分布的杂波RFS的概率为κk(z)，对j＝1，…，Jk+1｜k有：
实施效果
本实施例的跟踪场景为二维多目标场景，监视区域内有8个目标。新产生目标服从泊松分布。为了与原始PHD滤波器进行比较，设置5个目标起始位置已知，新目标的出现服从泊松分布，强度为：其中ωγ＝0.03，Ｐγ＝diag([501050100.1]′)。其余3个目标的起始位置和时刻均未知。
图2给出了仿真场景中的目标航迹。
图3是每一帧x方向和y方向的包括杂波的量测图。
图4是本发明方法在x和y方向上对时间的位置估计结果。从图4可以看出，本发明方法不仅可以准确估计出5个已知起始位置的新生目标，而且可以准确估计出3个未知起始位置的新生目标。
图5给出了由标准GM-PHD滤波器估计同一多目标场景的结果。其中PHD滤波部分的参数设置与本发明方法完全相同。可以看出，GM-PHD滤波器仅能准确估计出5个已知起始位置的新生目标，而不能估计出3个未知起始位置的新生目标。这是由于GM-PHD滤波器不具有起始位置未知新目标能力所决定的。
图6给出100次蒙特卡洛平均的真实目标数、标准GM-PHD滤波器和本发明方法估计的目标数的比较结果，其中红色实线表示真实目标数，蓝色虚线表示本发明方法估计的目标数，黑色带点的点画线表示标准GM-PHD滤波器估计的目标数。可以清除的看出，标准GM-PHD滤波器估计在新生目标出现后的时间段内给出了错误的目标数目估计，而本发明方法则全程均给出正确的目标数目估计。
图7给出了100次蒙特卡洛平均的标准GM-PHD滤波器和本发明方法估计的OSPA距离的比较图。OSPA距离是两个随机集合之间“距离”的一种度量，它可以用来评价跟踪算法的性能，其值越小说明算法的跟踪性能越好。
因此，本发明方法很好地解决了标准GM-PHD滤波器需要已知目标初始强度函数的问题，具有很强的实用性，有效性和鲁棒性。