CN106372646A - 基于srck‑gmcphd滤波的多目标跟踪方法 - Google Patents

Abstract

本发明提供了一种基于SRCK‑GMCPHD滤波的多目标跟踪方法，利用容积数值积分方法来处理非线性变换后随机变量的均值和方差，并同时使用GMCPHD算法对目标状态和数目进行估计；将平方根方法引入到CPHD的预测、更新和高斯元修剪过程中，使得仅有误差方差阵的平方根在整个过程中传递，增强了算法的数值精度和稳定性。本发明的多目标跟踪效果优于传统的EK‑GMCPHD方法，并且具有较广的适用性，能够被应用于雷达、红外目标跟踪、移动机器人定位等领域。

Description

基于SRCK-GMCPHD滤波的多目标跟踪方法

技术领域

本发明涉及目标跟踪领域的多目标跟踪方法，具体地，涉及一种基于平方根容积卡尔曼(square-root cubature Kalman,SRCK)-高斯混合基数概率假设密度(Gaussianmixture cardinalized probability hypothesis density,GMCPHD)滤波的多目标跟踪方法。

背景技术

多目标跟踪是根据传感器探测到的、由多个运动目标和环境噪声产生的观测序列，来估计多个目标运动状态和目标数目的目标跟踪方法。在实际多目标跟踪问题中，传感器的观测具有随机误差，且检测概率小于1，存在着虚警和漏警问题；目标在观测区域随机的出现和消失，目标的数目往往是未知的；存在航迹的交叉和分叉。上面的问题使得多目标跟踪具有一定的挑战性。

传统的多目标跟踪方法(如概率数据关联、联合概率数据关联和多假设跟踪)把关联和估计分成两个独立的部分进行，关联的精度对跟踪的影响比较大，在目标数目时并且未知的情况下难以应用。另外，这些方法的计算量随着目标数目和杂波密度的增加而急剧增长，存在“组合爆炸”问题，这限制了传统的多目标跟踪方法的实际应用。而基于随机有限集的方法在多目标跟踪的应用中有天然的优势，它避免了数据关联问题，可以在目标时变并且未知时应用，轨迹起始、维持与终结都是自然完成的，不需要单独列出。Mahler提出概率假设密度(probability hypothesis density，PHD)滤波器，使得基于随机有限集的多目标跟踪方法得以实用并得到了广泛的应用。实现PHD滤波器的主要方法有序贯Monte Carlo方法和高斯混合(Gaussian mixture probability hypothesis density GM-PHD)方法。序贯Monte Carlo滤波器可以处理非线性、非高斯情况，通用性较强，然而这也带来了计算量的增加，且需要单独采用聚类的方法得到目标的状态。GM-PHD滤波器峰值提取方法简单，在GM-PHD函数中提取出目标的状态较为容易，且可以通过高斯元的合并和剪枝有效的控制计算量。然而，PHD滤波器只递推目标的一阶矩，且假设虚警目服从泊松分布，这些都带来了一定的信息损失。为此，Mahler提出了基数概率假设密度(cardinalized probabilityhypothesis density，CPHD)滤波器，放松了泊松假 设，能够在传递PHD函数的同时传递目标数目分布的概率密度函数。B.T.Vo给出了CPHD的高斯混合实现。

基于RFS的多目标跟踪方法处理的系统往往是非线性的，这就需要采用能够处理非线性问题的跟踪方法。目前存在的高斯非线性滤波方法主要有扩展卡尔曼滤波-基数概率假设密度滤波、无迹卡尔曼-基数概率假设密度滤波等。扩展卡尔曼滤波-基数概率假设密度滤波方法在对系统进行线性化时需要计算雅克比矩阵，且仅具有一阶精度。无迹卡尔曼-基数概率假设密度滤波方法是基于确定性采样的滤波方法，可以得到三阶的精度，但当引入的尺度参数小于零时，可能导致更新后的方差阵为非正定阵，从而影响滤波的稳定性。另外，一般的CPHD方法都没有考虑在实际应用中可能出现的数值稳定性和数值精度问题。

发明内容

针对现有技术中的缺陷，本发明的目的是提供一种基于SRCK-GMCPHD滤波的多目标跟踪方法。

根据本发明提供基于SRCK-GMCPHD滤波的多目标跟踪方法，包括以下步骤：

步骤1：给定初始的基数分布以及高斯元集合，得到初始时刻多目标强度；

步骤2：初始步时，利用步骤1中的初始参数对目标跟踪的每一步进行基数分布预测，非初始步时，利用步骤i得到的参数进行基数分布预测；

步骤3：初始步时，利用步骤1中的初始参数对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，非初始步时，利用步骤i得到的参数对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，所述预测值包括对新生目标、衍生目标、存活目标以及已存在目标的预测；

步骤4：当存在新的传感器观测值时，采用SRCKF的更新方程对步骤3中的高斯元均值和方差的平方根进行更新，得到更新后的概率密度函数，并对基数分布进行更新，得到基数分布的更新值；

步骤5：对更新后的高斯元进行修剪，消除多余的高斯元；

步骤i：对修剪后的高斯元进行多目标状态提取，得到多目标状态。

优选地，所述步骤1中给定初始的基数分布记为：p₀(n)，高斯元集合记为：初始时刻多目标的强度记为：D₀(x)；上标i表示高斯元索引，n表示最大目标数，表示第i个高斯元的初始权重，表示第i个高斯元的初始状态， 表示第i个高斯元的初始状态协方差，J₀表示初始时刻高斯元数目；

通过概率假设密度函数进行预测，得到不同时刻的多目标强度，所述概率假设密度函数的计算公式如下：

其中：式中：表示经过平方根卡尔曼滤波预测得到的预测方差的平方根，D_k|k-1(x)表示多目标预测状态的概率假设密度函数，表示高斯分布随机变量x服从均值为方差为的高斯分布，表示的平方根，表示k-1时刻第j个高斯元的初始权重；下标_k|k-1表示k-1时刻到k时刻的预测；上标j表示第j个高斯元的参数，下标k表示k时刻。

优选地，所述步骤2包括：对目标跟踪的每一步进行基数分布预测，预测公式如下：

式中：p_Γ,k(·)表示k时刻新出现目标的基数分布函数；表示二项式系数，p_k|k-1(n)表示目标数目分布的预测概率密度函数，p_k-1(l)表示目标数目分布的先验概率密度函数，p_s,k表示目标生存概率；l为整数表示索引。

优选地，所述步骤3包括：采用SRCK中的预测步骤对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，所述预测值中包括对新生目标、衍生目标、存活目标以及已存在目标的预测；

所述概率假设密度函数的预测值计算公式如下：

D_k|k-1(x)＝D_S,k|k-1+β_k|k-1(x)+γ_k(x)；

式中：D_S,k|k-1表示k时刻存活目标预测强度，β_k|k-1(x)表示k时刻衍生目标预测强度，γ_k(x)表示k时刻新生目标强度；

新生目标的概率假设密度PHD预测公式如下：

式中：表示第j个新生目标高斯元的权重，和分别表示第j个的高斯元的均值和方差的平方根，J_γ,k表示k时刻新生目标的个数，下标γ表示新生目标；

衍生目标的PHD预测公式如下：

式中：表示第j个高斯元的权重，和分别表示第j个衍生目标高斯元的权重均值和方差的平方根；下标_(j,l)表示第j个高斯元的第l个衍生目标，下标β,k|k-1表示k-1时刻到k时刻衍生目标的预测；

存活目标的PHD预测公式如下：

式中：和表示存活目标的均值和方差；

已存在目标预测公式如下：

其中：m＝2(n_x+n_w+n_v)，n_x，n_w和n_v分别表示状态，状态误差方差和观测噪声方差的维数；[1]_p表示点集[1]的第p列，若[1]∈R²，则[1]表示点集如下：

式中：表示扩维的状态向量，表示状态扩维后方差的平方根，diag(·)表示对角阵运算，表示目标状态误差方差的平方根，S_w表示状态噪声方差的平方根，S_v表示观测噪声方差的平方根，表示k-1时刻扩维后状态sigma点，表示k时刻状态预测sigma点，表示k-1时刻状态sigma点，表示k-1时刻观测噪声sigma点，表示k时刻扩维后的状态预测sigma点，表示状态的第i个预测sigma 点，表示状态预测均值，qr(·)表示QR分解运算，表示过程噪声均值；表示观测噪声均值；J_k-1表示k-1时刻高斯元数目，m表示扩维后状态维数，上标x,(j)表示状态x的第j个sigma点，下标p,k-1表示k-1时刻第p个sigma点，下标S,k|k-1表示k-1时刻到k时刻存活目标的预测，下标m,k|k-1表示k-1时刻到k时刻第m个sigma点的预测，下标p,k|k-1表示k-1时刻到k时刻第p个的预测，qr(·)表示QR分解运算。

优选地，所述步骤4包括：当存在新的传感器观测值时，采用SRCK的更新方程对高斯元的均值和方差的平方根进行更新，得到更新后的概率密度函数，记为D_k(x)；对基数分布进行更新，得到基数分布的更新值，更新后的基数分布记为p_k(n)；

对高斯元的均值和方差的平方根进行更新的公式如下：

式中：表示k时刻观测预测的sigma点，f(·)表示状态方程，表示k时刻状态预测，表示k时刻观测噪声预测，表示新息，表示k时刻观测预测集合，表示第i个观测预测sigma点，表示观测估计误差方差的平方根，下标zz,k|k-1表示k-1时刻到k时刻量测预测协方差，z表示量测，表示k时刻状态的和观测的协方差，下标xz,k|k-1表示k-1时刻到k时刻量测和状态的互协方差；表示k时刻滤波增益，表示k时刻通过观测z得到的状态更新，表示k时刻状态预测，cholupdate{·}运算代表对矩阵进行Cholesky分解，cholupdate{S,U,±1}表示对矩阵进行Cholesky更新，S表示某一平方根矩阵，U表示某一向量或矩阵；即计算chol(SS^T±UU^T)，若U不是向量而是矩阵，那么cholupdate{·}表示用U矩阵的每一个列 向量连续进行更新，则矩阵A的QR分解可以表示为：A^T＝QR，其中，R表示上三角矩阵，而S＝qr(A)，则有S＝R^T；

p_k(n)和D_k(x)的计算公式为：

其中

式中：表示与p_k|k-1的内积；ω_k|k-1表示高斯元权重集合，Z_k表示为k时刻观测的集合，p_k|k-₁(n)表示预测分布基数，z表示观测，|Z|表示Z中元素个数，p_D,k表示检测概率，下标D,k表示k时刻检测概率，J_k|k-₁表示k时刻预测高斯元数目；表示量测z第j个高斯元的权重，表示量测z第j个高斯元的均值，＜1,ω＞^j+u表示ω和1内积的j+u次幂，q_k(z)表示k时刻量测z的似然，表示第J_k|k-1个高斯元的预测权重，表示k时刻量测z关于第j个高斯元的似然，κ_k(z)表示杂波强度函数，表示排列组合系数，<α,β>表示实函数α和β的内积，σ_j(·)表示非空实数集合Z阶数为j的均衡函数。

优选地，所述步骤5包括：采用平方根高斯元修剪方法对更新后的高斯元进行修剪，消除多余的高斯元。

优选地，所述步骤7包括：对修剪后的高斯元进行多目标状态提取，得到多目标状态输出具体地，通过提取权重大于某个阈值ω_Th的高斯元的方法得到多目标状态，计算公式如下：

${\hat{X}}_k=\{m_k^{\left(i\right)}:{\mathrm{\&omega;}}_k^{\left(i\right)}>{\mathrm{\&omega;}}_{Th},i=1,...,J_k\};$

式中：表示k时刻第i高斯元的均值，ω_Th表示高斯元修剪阈值，J_k表示k时刻更新后高斯元个数，表示k时刻高斯元权重。

与现有技术相比，本发明具有如下的有益效果：

1、本发明提供的基于SRCK-GMCPHD滤波的多目标跟踪方法利用容积数值积分方法来处理非线性变换后随机变量的均值和方差，实现较为为简单，使用GMCPHD方法同时对目标状态和数目进行估计，并将平方根方法引入到CPHD的预测、更新和高斯元修剪过程中，使得仅有误差方差阵的平方根在整个过程中传递，增强了方法的数值精度和稳定性。

2、本发明的多目标跟踪效果优于传统的EK-GMCPHD方法，并且具有较广的适用性，能够被应用于雷达、红外目标跟踪、移动机器人定位等领域。

附图说明

通过阅读参照以下附图对非限制性实施例所作的详细描述，本发明的其它特征、目的和优点将会变得更明显：

图1为二维平面四个目标的运动轨迹示意图；

图2为SRCK-GMCPHD对目标的估计值和真值的对比图；

图3为SRCK-GMCPHD方法的与EK-GMCPHD方法的OPSA距离对比图；

图4为SRCK-GMCPHD方法的与EK-GMCPHD方法在50次仿真中对目标数目估计的平均值对比图。

具体实施方式

下面结合具体实施例对本发明进行详细说明。以下实施例将有助于本领域的技术人员进一步理解本发明，但不以任何形式限制本发明。应当指出的是，对本领域的普通技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进。这些都属于本发明的保护范围。

根据本发明提供的基于SRCK-GMCPHD滤波的多目标跟踪方法，包括以下步骤：

步骤1：给定初始的基数分布以及高斯元集合，得到初始时刻多目标强度；

步骤2：初始步时，利用步骤1中的初始参数对目标跟踪的每一步进行基数分布预测，非初始步时，利用步骤i得到的参数进行基数分布预测；

步骤3：初始步时，利用步骤1中的初始参数对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，非初始步时，利用步骤i得到的参数对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，所述预测值包括对新生目标、衍生目标、存活目标以及已存在目标的预测；

步骤4：当存在新的传感器观测值时，采用SRCKF的更新方程对步骤3中的高斯元均值和方差的平方根进行更新，得到更新后的概率密度函数，并对基数分布进行更新，得到基数分布的更新值；

步骤5：对更新后的高斯元进行修剪，消除多余的高斯元；

步骤i：对修剪后的高斯元进行多目标状态提取，得到多目标状态。

所述步骤1中给定初始的基数分布记为：p₀(n)，高斯元集合记为：初始时刻多目标的强度记为：D₀(x)；上标i表示高斯元索引，n表示最大目标数，表示第i个高斯元的初始权重，表示第i个高斯元的初始状态，表示第i个高斯元的初始状态协方差，J₀表示初始时刻高斯元数目；

通过概率假设密度函数进行预测，得到不同时刻的多目标强度，所述概率假设密度函数的计算公式如下：

其中：式中：表示经过平方根卡尔曼滤波预测得到的预测方差的平方根，D_k|k-1(x)表示多目标预测状态的概率假设密度函数，表示高斯分布随机变量x服从均值为方差为的高斯分布，表示的平方根，表示k-1时刻第j个高斯元的初始权重；下标_k|k-1表示k-1时刻到k时刻的预测；上标j表示第j个高斯元的参数，下标k表示k时刻。

所述步骤2包括：对目标跟踪的每一步进行基数分布预测，预测公式如下：

式中：p_Γ,k(·)表示k时刻新出现目标的基数分布函数；表示二项式系数，p_k|k-1(n)表示目标数目分布的预测概率密度函数，p_k-1(l)表示目标数目分布的先验概率密度函数，p_s,k表示目标生存概率；l为整数表示索引。

所述步骤3包括：采用SRCK中的预测步骤对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，所述预测值中包括对新生目标、衍生目标、存活目标以及已存在目标的预测；

所述概率假设密度函数的预测值计算公式如下：

D_k|k-1(x)＝D_S,k|k-1+β_k|k-1(x)+γ_k(x)；

式中：D_S,k|k-1表示k时刻存活目标预测强度，β_k|k-1(x)表示k时刻衍生目标预测强度，γ_k(x)表示k时刻新生目标强度；

新生目标的概率假设密度PHD预测公式如下：

式中：表示第j个新生目标高斯元的权重，和分别表示第j个的高斯元的均值和方差的平方根，J_γ,k表示k时刻新生目标的个数，下标γ表示新生目标；

衍生目标的PHD预测公式如下：

式中：表示第j个高斯元的权重，和分别表示第j个衍生目标高斯元的权重均值和方差的平方根；下标(j,l)表示第j个高斯元的第l个衍生目标，下标β,k|k-1表示k-1时刻到k时刻衍生目标的预测；

存活目标的PHD预测公式如下：

式中：和表示存活目标的均值和方差；

已存在目标预测公式如下：

其中：m＝2(n_x+n_w+n_v)，n_x，n_w和n_v分别表示状态，状态误差方差和观测噪声方差的维数；[1]_p表示点集[1]的第p列，若[1]∈R²，则[1]表示点集如下：

式中：表示扩维的状态向量，表示状态扩维后方差的平方根，diag(·)表示对角阵运算，表示目标状态误差方差的平方根，S_w表示状态噪声方差的平方根，S_v表示观测噪声方差的平方根，表示k-1时刻扩维后状态sigma点，表示k时刻状态预测sigma点，表示k-1时刻状态sigma点，表示k-1时刻观测噪声sigma点，表示k时刻扩维后的状态预测sigma点，表示状态的第i个预测sigma点，表示状态预测均值，qr(·)表示QR分解运算，表示过程噪声均值；表示观测噪声均值；J_k-1表示k-1时刻高斯元数目，m表示扩维后状态维数，上标_x,(j)表示状态x的第j个sigma点，下标p,k-1表示k-1时刻第p个sigma点，下标S,k|k-1表示k-1时刻到k时刻存活目标的预测，下标m,k|k-1表示k-1时刻到k时刻第m个sigma点的预测，下标p,k|k-1表示k-1时刻到k时刻第p个的预测，qr(·)表示QR分解运算。

所述步骤4包括：当存在新的传感器观测值时，采用SRCK的更新方程对高斯元的均值和方差的平方根进行更新，得到更新后的概率密度函数，记为D_k(x)；对基数分布进行更新，得到基数分布的更新值，更新后的基数分布记为p_k(n)；

对高斯元的均值和方差的平方根进行更新的公式如下：

式中：表示k时刻观测预测的sigma点，f(·)表示状态方程，表示k时刻状态预测，表示k时刻观测噪声预测，表示新息，表示k时刻观测预测集合，表示第i个观测预测sigma点，表示观测估计误差方差的平方根，下标zz,k|k-1表示k-1时刻到k时刻量测预测协方差，z表示量测，表示k时刻状态的和观测的协方差，下标xz,k|k-1表示k-1时刻到k时刻量测和状态的互协方差；表示k时刻滤波增益，表示k时刻通过观测z得到的状态更新，表示k时刻状态预测，cholupdate{·}运算代表对矩阵进行Cholesky分解，cholupdate{S,U,±1}表示对矩阵进行Cholesky更新，S表示某一平方根矩阵，U表示某一向量或矩阵；即计算chol(SS^T±UU^T)，若U不是向量而是矩阵，那么cholupdate{·}表示用U矩阵的每一个列向量连续进行更新，则矩阵A的QR分解可以表示为：A^T＝QR，其中，R表示上三角矩阵，而S＝qr(A)，则有S＝R^T；

p_k(n)和D_k(x)的计算公式为：

其中

式中：表示与p_k|k-1的内积；ω_k|k-1表示高斯元权重集合，Z_k表示为k时刻观测的集合，p_k|k-1(n)表示预测分布基数，z表示观测，|Z|表示Z中元素个数，p_D,k表示检测概率，下标D,k表示k时刻检测概率，J_k|k-1表示k时刻预测高斯元数目；表示量测z第j个高斯元的权重，表示量测z第j个高斯元的均值，＜1,ω＞^j+u表示ω和1内积的j+u次幂，q_k(z)表示k时刻量测z的似然，表示第J_k|k-1个高斯元的预测权重，表示k时刻量测z关于第j个高斯元的似然，κ_k(z)表示杂波强度函数，表示排列组合系数，<α,β>表示实函数α和β的内积，σ_j(·)表示非空实数集合Z阶数为j的均衡函数。

所述步骤5包括：采用平方根高斯元修剪方法对更新后的高斯元进行修剪，消除多余的高斯元。

所述步骤7包括：对修剪后的高斯元进行多目标状态提取，得到多目标状态输出，具体地，通过提取权重大于某个阈值ω_Th的高斯元的方法得到多目标状态，计算公式如下：

式中：表示k时刻第i高斯元的均值，ω_Th表示高斯元修剪阈值，J_k表示k时刻更新后高斯元个数，表示k时刻高斯元权重。

下面结合具体实施例对本发明做更加详细的说明。

考虑一个二维平面的跟踪问题，假设目标的状态为其中ω_k为转弯速率；该向量包含了目标的位置、速度和加速度。场景中共有四个目标，其运动如图1所示：目标1和2从t＝0时刻开始存在；目标3为新生目标在t＝80s时开始出现；目标4为衍生目标，在t＝200s时出现，到300s时消失。假设雷达的观测量为斜距，径向速度和方位角。雷达距离观测的噪声标准差100m，径向速度观测标准差为10m/s，角度标准差0.2°，观测采样周期4s。目标存活概率p_S＝0.99，检测概率p_D＝0.99。观测中的杂波随机有限集服从泊松分布，其概率密度为κ_k(z)＝λ_cVu(z)。杂波密度为λ_c＝1.9×10^-10m^-2,观测区域面积V为＝2.6×10¹⁰m²。高斯元剪枝的阈值 T＝10^-4，合并阈值U＝5，状态提取的阈值ω_Th＝0.5，最大高斯元数目J_max＝200。

给定上面的初值和仿真参数后，每个仿真周期内，具体步骤如下所述：

步骤S1：基数预测

步骤S2：新生目标预测和更新元素构造；

步骤S3：衍生目标预测和更新元素构造；

步骤S4：存活目标预测和更新元素构造

步骤S5：更新元素构造和更新更新

步骤S6：平方根高斯元修剪；

步骤S7：多目标状态提取。

本实施例使用Matlab语言对所提出的方法进行了测试，并和传统的扩展卡尔曼-高斯混合基数概率假设密度滤波器(EK-GMCPHD)进行了对比。分别对EK-GMCPHD和SRCK-GMCPHD进行50次蒙特卡洛仿真，结果如图2-图4和表1所示。

图2给出了SRCK-GMCPHD对目标的估计值和真值，由图可以看出SRCK-GMCPHD方法对存活目标、新生目标和衍生目标三类目标均能进行较为有效跟踪，“错跟”和“漏跟”的次数都很少，这是因为方法中对这三类目标均进行了“显式”的建模。

图3给出了两种方法的OPSA距离，可以看出SRCK-GMCPHD方法的OSPA距离总体来说要小于EK-GMCPHD，这说明其对目标的跟踪精度更高。

图4给出了两种方法50次仿真中对目标数目估计的平均值，由图可以看出EK-GMCPHD方法和SRCK-GMCPHD方法均可以正确的估计目标数目，SRCK-GMCPHD方法要略好于EK-GMCPHD方法，对目标数目的准确估计主要是因为采用了CPHD方法，在对目标状态进行递推的同时还对目标数目的分布进行递推。

表1给出了两种方法OSPA距离的平均值和对目标数目估计的均方根误差，可以看出SRCK-GMCPHD方法对目标状态和数目的估计均优于EK-GMCPHD方法。另外，在仿真中还发现EK-GMCPHD会出现病态矩阵的情况，而SRCK-GMCPHD则一直具有良好的数值稳定性。

表1两种方法的比较

综上可以看出，本发明提出的方法能有效的实现对目标状态和目标数目的估计，其效果优于EK-GMCPHD方法。

以上对本发明的具体实施例进行了描述。需要理解的是，本发明并不局限于上述特定实施方式，本领域技术人员可以在权利要求的范围内做出各种变形或修改，这并不影响本发明的实质内容。

Claims (7)

1.一种基于SRCK-GMCPHD滤波的多目标跟踪方法，其特征在于，包括以下步骤：

步骤1：给定初始的基数分布以及高斯元集合，得到初始时刻多目标强度；

步骤2：初始步时，利用步骤1中的初始参数对目标跟踪的每一步进行基数分布预测，非初始步时，利用步骤i得到的参数进行基数分布预测；

步骤3：初始步时，利用步骤1中的初始参数对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，非初始步时，利用步骤i得到的参数对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，所述预测值包括对新生目标、衍生目标、存活目标以及已存在目标的预测；

步骤4：当存在新的传感器观测值时，采用SRCKF的更新方程对步骤3中的高斯元均值和方差的平方根进行更新，得到更新后的概率密度函数，并对基数分布进行更新，得到基数分布的更新值；

步骤5：对更新后的高斯元进行修剪，消除多余的高斯元；

步骤i：对修剪后的高斯元进行多目标状态提取，得到多目标状态。

2.根据权利要求1所述的基于SRCK-GMCPHD滤波的多目标跟踪方法，其特征在于，所述步骤1中给定初始的基数分布记为：p₀(n)，高斯元集合记为：初始时刻多目标的强度记为：D₀(x)；上标i表示高斯元索引，n表示最大目标数，表示第i个高斯元的初始权重，表示第i个高斯元的初始状态，表示第i个高斯元的初始状态协方差，J₀表示初始时刻高斯元数目；

通过概率假设密度函数进行预测，得到不同时刻的多目标强度，所述概率假设密度函数的计算公式如下：

$D_{k|k-1}\left(x\right)=\underset{i=1}{\overset{J_{k|k-1}}{\mathrm{\&Sigma;}}}{\mathrm{\&omega;}}_{k|k-1}^{\left(i\right)}N(x;m_{k|k-1}^{\left(i\right)},S_{k|k-1}^{\left(i\right)});$

其中：式中：表示经过平方根卡尔曼滤波预测得到的预测方差的平方根，D_k|k-1(x)表示多目标预测状态的概率假设密度函数，表示高斯分布随机变量x服从均值为方差为的高斯分布，表示的平方根，表示k-1时刻第j个高斯元的初始权重；下标k|k-1表示k-1时刻到k时刻的预测；上标j表示第j个高斯元的参数，下标k表示k时刻。

3.根据权利要求2所述的基于SRCK-GMCPHD滤波的多目标跟踪方法，其特征在于，所述步骤2包括：对目标跟踪的每一步进行基数分布预测，预测公式如下：

$p_{k|k-1}\left(n\right)=\underset{i=0}{\overset{n}{\&Sigma;}}{p}_{\mathrm{\&Gamma;},k}\left(n-i\right)\underset{l=i}{\overset{\mathrm{\&infin;}}{\&Sigma;}}{C}_i^l{p}_{k-1}\left(l\right)p_{s,k}^i{(1-p_{s,k})}^{l-i};$

式中：p_Γ,k(·)表示k时刻新出现目标的基数分布函数；表示二项式系数，p_k|k-1(n)表示目标数目分布的预测概率密度函数，p_k-1(l)表示目标数目分布的先验概率密度函数，p_s,k表示目标生存概率；l为整数表示索引。

4.根据权利要求3所述的基于SRCK-GMCPHD滤波的多目标跟踪方法，其特征在于，所述步骤3包括：采用SRCK中的预测步骤对高斯元集合的均值和方差的平方根进行预测后得到概率假设密度函数的预测值，所述预测值中包括对新生目标、衍生目标、存活目标以及已存在目标的预测；

所述概率假设密度函数的预测值计算公式如下：

D_k|k-1(x)＝D_S,k|k-1+β_k|k-1(x)+γ_k(x)；

式中：D_S,k|k-1表示k时刻存活目标预测强度，β_k|k-1(x)表示k时刻衍生目标预测强度，γ_k(x)表示k时刻新生目标强度；

新生目标的概率假设密度PHD预测公式如下：

${\mathrm{\&gamma;}}_k\left(x\right)=\underset{j=1}{\overset{J_{\mathrm{\&gamma;},k}}{\&Sigma;}}{\mathrm{\&omega;}}_{\mathrm{\&gamma;},k}^{\left(j\right)}N(x;m_{\mathrm{\&gamma;},k}^{\left(j\right)},S_{\mathrm{\&gamma;},k}^{\left(j\right)});$

式中：表示第j个新生目标高斯元的权重，和分别表示第j个的高斯元的均值和方差的平方根，J_γ,k表示k时刻新生目标的个数，下标γ表示新生目标；

衍生目标的PHD预测公式如下：

${\mathrm{\&beta;}}_{k|k-1}\left(x\right)=\underset{j=1}{\overset{J_{k-1}}{\&Sigma;}}\underset{l=1}{\overset{J_{\mathrm{\&beta;},k}}{\&Sigma;}}{\mathrm{\&omega;}}_{k-1}^{\left(j\right)}{\mathrm{\&omega;}}_{\mathrm{\&beta;},k}^{\left(l\right)}N(x;m_{\mathrm{\&beta;},k|k-1}^{(j,l)},S_{\mathrm{\&beta;},k|k-1}^{(j,l)}{\left(S_{\mathrm{\&beta;},k|k-1}^{(j,l)}\right)}^T);$

式中：表示第j个高斯元的权重，和分别表示第j个衍生目标高斯元的权重均值和方差的平方根；下标(j,l)表示第j个高斯元的第l个衍生目标，下标β,k|k-1表示k-1时刻到k时刻衍生目标的预测；

存活目标的PHD预测公式如下：

$D_{S,k|k-1}\left(x\right)=p_{S,k}\underset{j=1}{\overset{J_{k-1}}{\mathrm{\&Sigma;}}}{\mathrm{\&omega;}}_{k-1}^{\left(j\right)}N(x;m_{S,k|k-1}^{\left(j\right)},S_{S,k|k-1}^{\left(j\right)}{\left(S_{S,k|k-1}^{\left(j\right)}\right)}^T);$

式中：和表示存活目标的均值和方差；

已存在目标预测公式如下：

${\hat{x}}_{k-1}^{\left(j\right)}=\left[\begin{array}{ccc}{m}_{k-1}^{\left(j\right)}& \stackrel{\&OverBar;}{w}& \stackrel{\&OverBar;}{v}\end{array}\right];$

$S_{k-1}^{a,\left(j\right)}=diag(S_{k-1}^{\left(j\right)},S_w,S_v);$

${\mathrm{\&chi;}}_{p,k-1}^{\left(j\right)}=S_{k-1}^{a,\left(j\right)}{\mathrm{\&xi;}}_p+{\hat{x}}_{k-1}^{\left(j\right)};$

${\mathrm{\&chi;}}_{p,k|k-1}^{x,\left(j\right)}=f\left({\mathrm{\&chi;}}_{p,k-1}^{x,\left(j\right)},{\mathrm{\&chi;}}_{p,k-1}^{w,\left(j\right)}\right);$

$m_{S,k|k-1}^{\left(j\right)}=\frac{1}{m}\underset{p=1}{\overset{m}{\&Sigma;}}{\mathrm{\&chi;}}_{p,k|k-1}^{x,\left(j\right)};$

${\mathrm{\&chi;}}_{k|k-1}^{\left(j\right)}=\frac{1}{\sqrt{m}}\&lsqb;{\mathrm{\&chi;}}_{1,k|k-1}^{x,\left(j\right)}-m_{S,k|k-1}^{\left(j\right)},...,{\mathrm{\&chi;}}_{m,k|k-1}^{x,\left(j\right)}-m_{S,k|k-1}^{\left(j\right)}\&rsqb;;$

$S_{S,k|k-1}^{\left(j\right)}=qr\left({\mathrm{\&chi;}}_{k|k-1}^{\left(j\right)}\right);$

其中：n_x，n_w和n_v分别表示状态，状态误差方差和观测噪声方差的维数；[1]_p表示点集[1]的第p列，若[1]∈R²，则[1]表示点集如下：

$\left\{\left(\begin{array}{c}1\\ 0\end{array}\right),\left(\begin{array}{c}0\\ 1\end{array}\right),\left(\begin{array}{c}-1\\ 0\end{array}\right),\left(\begin{array}{c}0\\ -1\end{array}\right)\right\};$

式中：表示扩维的状态向量，表示状态扩维后方差的平方根，diag(·)表示对角阵运算，表示目标状态误差方差的平方根，S_w表示状态噪声方差的平方根，S_v表示观测噪声方差的平方根，表示k-1时刻扩维后状态sigma点，表示k时刻状态预测sigma点，表示k-1时刻状态sigma点，表示k-1时刻观测噪声sigma点，表示k时刻扩维后的状态预测sigma点，表示状态的第i个预测sigma点，表示状态预测均值，qr(·)表示QR分解运算，表示过程噪声均值；表示观测噪声均值；J_k-1表示k-1时刻高斯元数目，m表示扩维后状态维数，上标x,(j)表示状态x的第j个sigma点，下标p,k-1表示k-1时刻第p个sigma点，下标S,k|k-1表示k-1时刻到k时刻存活目标的预测，下标m,k|k-1表示k-1时刻到k时刻第m个sigma点的预测，下标p,k|k-1表示k-1时刻到k时刻第p个的预测，qr(·)表示QR分解运算。

5.根据权利要求4所述的基于SRCK-GMCPHD滤波的多目标跟踪方法，其特征在于，所述步骤4包括：当存在新的传感器观测值时，采用SRCK的更新方程对高斯元的均值和方差的平方根进行更新，得到更新后的概率密度函数，记为D_k(x)；对基数分布进行更新，得到基数分布的更新值，更新后的基数分布记为p_k(n)；

对高斯元的均值和方差的平方根进行更新的公式如下：

$Z_{p,k|k-1}^{\left(j\right)}=f\left({\mathrm{\&chi;}}_{p,k|k-1}^{x,\left(j\right)},{\mathrm{\&chi;}}_{p,k|k-1}^{v,\left(j\right)}\right)$

${\mathrm{\&eta;}}_{k|k-1}^{\left(j\right)}=\frac{1}{m}\underset{p=1}{\overset{m}{\mathrm{\&Sigma;}}}{Z}_{p,k|k-1}^{\left(j\right)}$

$Z_{k|k-1}^{\left(j\right)}=\frac{1}{\sqrt{m}}\&lsqb;Z_{1,k|k-1}^{\left(j\right)}-{\mathrm{\&eta;}}_{k|k-1}^{\left(j\right)},...,Z_{m,k|k-1}^{\left(j\right)}-{\mathrm{\&eta;}}_{k|k-1}^{\left(j\right)}\&rsqb;$

$S_{zz,k|k-1}^{\left(j\right)}=qr\left(Z_{k|k-1}^{\left(j\right)}\right)$

$P_{xz,k|k-1}^{\left(j\right)}={\mathrm{\&chi;}}_{k|k-1}^{\left(j\right)}{\left(Z_{k|k-1}^{\left(j\right)}\right)}^T$

$K_k^{\left(j\right)}=(P_{xz,k|k-1}^{\left(j\right)}/{\left(S_{zz,k|k-1}^{\left(j\right)}\right)}^T)/S_{zz,k|k-1}^{\left(j\right)}$

$m_k^{\left(j\right)}\left(z\right)=m_{k|k-1}^{\left(j\right)}+K_k^{\left(j\right)}(z-{\mathrm{\&eta;}}_{k|k-1}^{\left(j\right)})$

$S_k^{\left(j\right)}=cholupdate\{S_{k|k-1}^{\left(j\right)},K_k^{\left(j\right)}{S}_{zz,k|k-1}^{\left(j\right)},-1\}$

式中：表示k时刻观测预测的sigma点，f(·)表示状态方程，表示k时刻状态预测，表示k时刻观测噪声预测，表示新息，表示k时刻观测预测集合，表示第i个观测预测sigma点，表示观测估计误差方差的平方根，下标zz,k|k-1表示k-1时刻到k时刻量测预测协方差，z表示量测，表示k时刻状态的和观测的协方差，下标xz,k|k-1表示k-1时刻到k时刻量测和状态的互协方差；表示k时刻滤波增益，表示k时刻通过观测z得到的状态更新，表示k时刻状态预测，cholupdate{·}运算代表对矩阵进行Cholesky分解，cholupdate{S,U,±1}表示对矩阵进行Cholesky更新，S表示某一平方根矩阵，U表示某一向量或矩阵；即计算chol(SS^T±UU^T)，若U不是向量而是矩阵，那么cholupdate{·}表示用U矩阵的每一个列向量连续进行更新，则矩阵A的QR分解可以表示为：A^T＝QR，其中，R表示上三角矩阵，而S＝qr(A)，则有S＝R^T；

p_k(n)和D_k(x)的计算公式为：

$p_k\left(n\right)=\frac{{\mathrm{\&Psi;}}_k^0\&lsqb;{\mathrm{\&omega;}}_{k|k-1},Z_k\&rsqb;\left(n\right)p_{k|k-1}\left(n\right)}{<{\mathrm{\&Psi;}}_k^0\&lsqb;{\mathrm{\&omega;}}_{k|k-1},Z_k\&rsqb;,p_{k|k-1}>};$

$D_k\left(x\right)=\frac{<{\mathrm{\&Psi;}}_k^1\&lsqb;{\mathrm{\&omega;}}_{k|k-1},Z_k\&rsqb;,p_{k|k-1}>}{<{\mathrm{\&Psi;}}_k^0\&lsqb;{\mathrm{\&omega;}}_{k|k-1},Z_k\&rsqb;,p_{k|k-1}>}(1-p_{D,k})D_{k|k-1}\left(x\right)+\underset{z\&Element;Z_k}{\&Sigma;}\underset{j=1}{\overset{J_{k|k-1}}{\mathrm{\&Sigma;}}}{\mathrm{\&omega;}}_k^{\left(j\right)}\left(z\right)N(x;m_k^{\left(j\right)}(z),S_k^{\left(j\right)}{\left(S_k^{\left(j\right)}\right)}^T);$

其中

${\mathrm{\&Psi;}}_k^u\&lsqb;\mathrm{\&omega;},Z\&rsqb;\left(n\right)=\underset{j=0}{\overset{\min \left(\left|Z\right|,n\right)}{\mathrm{\&Sigma;}}}\left(\left|Z\right|-j\right)p_{K,k}\left(\left|Z\right|-j\right)P_{j+u}^n\&times;\frac{{\left(1-p_{D,k}\right)}^{n-\left(j+u\right)}}{<1,\mathrm{\&omega;}{>}^{j+u}}{\mathrm{\&sigma;}}_j\left({\mathrm{\&Lambda;}}_k\left(\mathrm{\&omega;},Z\right)\right);$

${\mathrm{\&Lambda;}}_k\left(\mathrm{\&omega;},Z\right)=\left\{\frac{<1,{\mathrm{\&kappa;}}_k>}{{\mathrm{\&kappa;}}_k\left(z\right)}{p}_{D,k}{\mathrm{\&omega;}}^T{q}_k\left(z\right):z\&Element;Z\right\};$

${\mathrm{\&omega;}}_{k|k-1}={\&lsqb;{\mathrm{\&omega;}}_{k|k-1}^{\left(1\right)},...,{\mathrm{\&omega;}}_{k|k-1}^{\left(J_{k|k-1}\right)}\&rsqb;}^T;$

$q_k\left(z\right)={\&lsqb;q_k^{\left(1\right)}\left(z\right),...,q_k^{\left(J_{k|k-1}\right)}\left(z\right)\&rsqb;}^T;$

$q_k^{\left(j\right)}\left(z\right)=N\left(z,{\mathrm{\&eta;}}_{k|k-1}^{\left(j\right)},S_{zz,k|k-1}^{\left(j\right)}\right);$

$P_j^n=\frac{n!}{(n-j)!}$

式中：表示与p_k|k-1的内积；ω_k|k-1表示高斯元权重集合，Z_k表示为k时刻观测的集合，p_k|k-1(n)表示预测分布基数，z表示观测，|Z|表示Z中元素个数，p_D,k表示检测概率，下标D,k表示k时刻检测概率，J_k|k-1表示k时刻预测高斯元数目；表示量测z第j个高斯元的权重，表示量测z第j个高斯元的均值，＜1,ω＞^j+u表示ω和1内积的j+u次幂，q_k(z)表示k时刻量测z的似然，表示第J_k|k-1个高斯元的预测权重，表示k时刻量测z关于第j个高斯元的似然，κ_k(z)表示杂波强度函数，表示排列组合系数，＜α,β＞表示实函数α和β的内积，σ_j(·)表示非空实数集合Z阶数为j的均衡函数。

6.根据权利要求5所述的基于SRCK-GMCPHD滤波的多目标跟踪方法，其特征在于，所述步骤5包括：采用平方根高斯元修剪方法对更新后的高斯元进行修剪，消除多余的高斯元。

7.根据权利要求6所述的基于SRCK-GMCPHD滤波的多目标跟踪方法，其特征在于，所述步骤7包括：对修剪后的高斯元进行多目标状态提取，得到多目标状态输出具体地，通过提取权重大于某个阈值ω_Th的高斯元的方法得到多目标状态，计算公式如下：

${\hat{X}}_k=\{m_k^{\left(i\right)}:{\mathrm{\&omega;}}_k^{\left(i\right)}>{\mathrm{\&omega;}}_{Th},i=1,...,J_k\};$

式中：表示k时刻第i高斯元的均值，ω_Th表示高斯元修剪阈值，J_k表示k时刻更新后高斯元个数，表示k时刻高斯元权重。