具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。

需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。

具体实施方式一：下面结合图1说明本实施方式，本实施方式所述一种被动声呐多目标方位检测跟踪方法，包括：

步骤一、采用波束形成技术对声呐观测信号进行处理，获得多目标的方位历程图，利用所述多目标的方位历程图，将方位历程图中每个采样时刻的方位谱的极大值点作为检测单元，并获取每个检测单元的自适应噪声背景；

步骤二、利用N-P准则检测器，利用每个检测单元的自适应噪声背景，对检测单元进行检测，获得量测点迹，利用初始两个时刻的量测点迹，获取多个目标初始航迹；

步骤三、采用跟踪波门对多个目标初始航迹进行筛选，获得多个真实目标航迹，根据每个真实目标跟踪波门内当前采样点时刻落入量测点迹的个数确定目标跟踪的方法；

若真实目标跟踪波门内仅落入一个量测点迹，则采用MPUKF滤波器对跟踪波门内量测点迹滤波；实现目标方位检测跟踪；

若真实目标踪波门内落入多个点迹，判断跟踪波门是否相交，若是，则真实目标航迹相交，否则真实目标航迹不相交；

若真实目标航迹相交，采用MHT算法对跟踪波门内量测形成假设事件，并计算假设事件的概率，获得目标的状态更新方程；执行步骤四；

若真实目标航迹不相交，采用PDA算法分别对每个目标的量测点迹进行处理，利用波门内所有量测点迹与目标的关联概率对量测点迹进行加权，以加权和作为等效结果，获得目标的状态更新方程；执行步骤四；

步骤四、根据目标航迹是否相交，利用对应的目标状态更新方程与MPUKF滤波技术结合对目标方位进行跟踪；

步骤五、跟踪过程中，跟踪波门在A个时刻内均未出现量测点迹时，则判定目标航迹终结，完成对目标方位检测跟踪，其中，A为大于3的正整数。本实施方式中在跟踪的过程中，当波门内未出现量测点迹时，利用滤波器的一步预测结果对目标状态进行估计以保持对目标的跟踪。若一段时间后目标波门内仍然无量测点迹出现，则说明目标航迹确实终结，结束对该目标的跟踪。

判断目标航迹是否终结具体过程为：

设连续无检测点迹时刻的门限个数A，连续无检测点迹时刻个数为i，当i＜A时进行以下循环：

步骤五一、设当前时刻k无量测点迹，i＝1，判断k+i时刻波门内是否有量测点迹；

若是，则利用量测点迹进行滤波跟踪，保持目标航迹，航迹终结判决结束；

否则，利用一步预测值对目标状态进行估计，将一步预测结果作为目标跟踪滤波结果，同时连续无检测点迹时刻个数为i+1。

步骤五二、判断i＜A是否成立；

若是，返回步骤四一继续递推循环；否则判定目标航迹终结。

进一步地，本实施方式中，步骤二中，利用初始两个时刻的量测点迹，获取多个目标初始航迹的具体方法为:

步骤A1、采用N-P准则检测器，利用每个检测单元的自适应噪声背景，对每个采样时刻的检测单元进行检测，判断每个采样时刻的检测单元是否存在量测点迹，若存在量测点迹，执行步骤A2；

步骤A2、根据初始两时刻检测出的量测点迹，获得目标初始航迹。进一步地，本实施方式中，

步骤A1中判断检测单元是否为量测点迹的具体方法为：

步骤A11、计算声呐测量环境背景噪声功率的水平估计值

步骤A12、计算检测单元的功率

步骤A13、利用所述估计值和已知的虚警概率P_fa，获取检测门限D_T；

步骤A14、将检测单元的功率与检测门限D_T比较，当时，存在量测点迹。

进一步地，本实施方式中，

步骤A2中获得目标初始航迹的具体方法为：

步骤A21、将初始时刻检测到的所有量测点迹作为航迹起始点，以每一个量测点迹为中心建立一个初始方位波门；

步骤A22、通过公式：

|z(2)-z(1)|≤2θ_BW

判断每个初始方位波门内第二时刻是否有量测点迹落入，若有量测点迹落入，则根据初始方位波门内点迹的个数和距离获取一条初始航迹，若无有点迹落入，则删除对应作航迹起始点，其中，θ_BW＝2arcsin(λ/md)为波束主瓣的半功率点波束宽度，λ为信号波长，m为阵元个数，d为阵元间距，z(k)(k＝1,2)是初始时刻和第二时刻检测到的量测点迹。

若第二时刻落入初始方位波门内一个点迹，则两个时刻点迹形成一条初始航迹；

若第二时刻落入初始方位波门内多个点迹，则取距离最近的点迹与起始点形成初始航迹。

本实施方式中，利用波束形成技术可得到多目标的方位历程图，方位历程图是对所有采样时刻目标方位谱的记录和显示，将每个采样时刻得到的方位谱的极大值点作为检测单元进行检测(检测的过程即为将所述极大值点的功率与门限值比较)。由于波束形成的波束存在主瓣宽度，取主瓣的半功率点波束宽度θ_BW为检测时保护单元的长度，即左右保护单元的长度均为θ_BW/2，并通过保护单元左右长度均为θ_BW的参考单元求取背景噪声功率的水平估计值

θ_BW由式：

θ_BW＝2arcsin(λ/md)

计算获得，式中，λ为接收信号的波长，m为阵元个数，d为阵元间距，对于均匀直线阵d＝λ/2。

由于受到环境噪声、水声信道以及波束旁瓣的影响，在选取参考单元背景时剔除掉最大的r₁个和最小的r₂个噪声背景，防止不均匀的背景使门限过高或过低。对参考单元内剩余的噪声背景求平均，即可得到噪声背景的平均功率

检测单元的功率即为功率谱的峰值功率信噪比即为：

再应用N-P准则，通过与给定的虚警概率P_fa可由式：

计算求出检测门限D_T，其中，Q(·)为标准正态分布的右尾函数。

将与D_T比较，若满足式：

则判为有目标，该点能量不变，否则判为无目标，能量置为0。

根据信噪比SNR可通过式：

求出检测概率P_D，式中，H₀代表目标不存在的情况，H₁代表目标存在的情况。

本实施方式中，经过检测后得到的每个点迹都对应一个量测方位值，表示k时刻测得的量测点迹的集合，其中，z_i(k)为k时刻第i个点迹，n(k)为k时刻全部点迹的个数。将z(1)中的所有点迹看作航迹起始点，以每一个点迹为中心建立初始方位波门，根据z(2)中的点迹是否满足：

|z(2)-z(1)|≤2θ_BW

判断第二时刻是有点迹落入初始波门内，若第二时刻初始波门内落入一个点迹，则两次点迹形成一条初始航迹；若初始波门内落入多个点迹，则取距离最近的点迹与起始点形成初始航迹；若无点迹落入初始波门内，则认为该轨迹起始点为虚假量测，将该起始点删除。

进一步地，本实施方式中，

步骤三中获得真实目标航迹的具体方法为：

步骤B1、利用目标初始航迹对修正极坐标系无迹卡尔曼滤波器进行初始化，利用初始化后的修正极坐标系无迹卡尔曼滤波器对目标初始航迹作预测，获得下一采样点时刻的预测值，以所述预测值为中心设置跟踪波门；

步骤B2、判断第三至第五采样点时刻是否有量测点迹落入对应的跟踪波门内，若是至少一个采样点时刻有量测点迹落入对应的跟踪波门内，则所述目标初始航迹为真实目标航迹。

本实施方式中通过形成初始轨迹的两个点迹对卡尔曼滤波器进行初始化，预测出下一时刻目标的一步预测值，以预测值为中心设置跟踪波门，即目标的量测z(k)是否满足：

[z(k)-z(k|k-1)]′·S^-1(k)·[z(k)-z(k|k-1)]≤γ

式中，γ为波门参数；P(k|k-1)为预测值；S(k)为k时刻新息协方差：

S(k)＝H(k)·P(k|k-1)·H′(k)+R(k)

由时间方位历程图得到目标方位量测信息随时间的变化结果，某时刻检测得到的量测点迹会出现在一条线形区域上，取跟踪波门以预测值为中心，长度为且量测值z(k)为n_z＝1维，由下表给出量测维数n_z为1时，不同参数λ对应的门概率P_G。

表1n_z＝1时落入波门内的概率P_G

由于航迹起始的时刻数需由目标数量及其相对位置、检测概率、虚警概率而定，具体实施时可选取起始步数N＝5。

若下一时刻有点迹落入波门内，则可认为初始航迹是真实存在的目标航迹；若没有点迹落入波门内，则以一步预测值近似代替该时刻的量测值再进行一次预测，重复判断二或三次。若后续两时刻内，有点迹落入波门的情况，也将其确认为目标航迹；若后续两时刻波门内仍无点迹落入，则认为该航迹是虚假航迹予以取消。

进一步地，本实施方式中，

步骤三中若真实目标跟踪波门内落入多个点迹，判断跟踪波门是否相交的具体方法为：

判断目标t在第k个采样点时刻的量测预测值z^t(k|k-1)(t＝1,2)是否满足：

若是，则波门相交，否则波门不相交，其中，z¹(k|k-1)为目标1在k采样点时刻的量测预测值，z²(k|k-1)为目标2在k采样点时刻的量测预测值,γ为波门参数。

进一步地，本实施方式中，

步骤二中，若真实目标航迹相交，获得目标的状态更新方程的具体方法为：

步骤C1、根据上一采样点时刻假设事件和当前时刻的量测集z(k)，产生当前采样点时刻假设事件；所述假设事件为：假设跟踪波门内量测点迹属于真实目标航迹的个数为：

波门内属于虚假航迹的量测点迹的个数为：

φ＝m_k-τ

其中，τ为跟踪波门内属于真实目标航迹的量测点迹的个数，τ_i为第i个属于真实目标航迹的量测点迹，φ为波门内属于虚假航迹的量测点迹的个数，m_k为跟踪波门内量测点迹的总个数；

步骤C2、计算步骤C1中假设事件成立的概率，根据假设事件成立的概率，获得跟踪波门内每个量测点迹与交叉的两个目标关联的概率；

步骤C3、对跟踪波门内每个量测点迹与交叉的两个目标关联的概率求和，获得每个量测点迹与交叉的两个目标的关联概率；

步骤C4、根据每个量测点迹与交叉的两个目标的关联概率，获得跟踪波门对应目标的状态更新方程：

式中：V^t(k)是互联于目标t的组合新息；是k时刻第i个量测点迹与目标t的组合新息；X^t(k|k)是k时刻目标t的状态更新值；为k时刻的第i个量测点迹z_i(k)与目标t互联的概率；X^t(k|k-1)是目标t在状态k采样点时刻的预测值；K^t(k)是目标t的滤波增益。

在本发明中，只用MHT算处理两个目标航迹交叉时的情况，所以点迹z_i(k)的假设是目标1、目标2航迹的延续或是虚假量测，而不考虑有新目标的产生，每个目标只能和一个量测互联或无量测互联。

若k-1时刻结束时产生m种假设事件Θ^k-1,s(s＝1,…,m)，在时刻k生成n种可能的假设，则在k时刻结束时一共会产生mn种假设事件Θ^k,l(l＝1,…,mn)。

设关于当前量测的事件Θ(k)中包含：τ个源于已建立航迹的量测；φ个虚假量测。对于i＝1，…，m_k定义与事件Θ(k)有关的以下标记变量

则源于已建立航迹的量测数：

虚假量测数：

φ＝m_k-τ

计算假设概率：

对任一假设事件Θ^k,l可以导出其概率与假设事件Θ^k-1,s之间的递归关系，则其概率计算公式为：

式中：c＝P{z(k)|z^k-1}为归一化常数因子；

为两航迹目标的检测概率；N_t[z_i(k)]表示与目标t关联的量测服从高斯分布：

N_t[z_i(k)]＝N[z_i(k)；z_t(k|k-1)S_t(k)]

其中，z_t(k|k-1)是目标t的预测量测；S_t(k)是其相应的新息协方差。

量测点迹z_i(k)与目标t(t＝1,2)互联的概率是k时刻结束时所有假设事件Θ^k,l中第i个量测z_i(k)源于目标t的假设事件概率之和。

状态更新与协方差更新：

MHT算法的状态更新和协方差更新，针对单个目标t，其状态更新方程的表达式为：

式中：V^t(k)是互联于目标t的组合新息；X^t(k|k)是k时刻目标t的状态估计更新值；X^t(k|k-1)是目标t的状态预测值；K^t(k)是目标t的滤波增益。

目标t的协方差更新方程为：

式中：p^t(k|k)为k时刻目标t的协方差更新值；

p^t(k|k-1)为目标t的协方差一步预测值；

为无量测与目标t关联的假设事件的概率；

应用MHT算法进行滤波更新时，采用上述目标t的协方差更新方程与态更新方程进行计算，

采用宽容法进行假设删除，来解决MHT假设事件组合爆炸的问题，即保留P(P＞0)个概率较大的假设，这里取P＝5。

进一步地，本实施方式中，

步骤三中，若真实目标航迹不相交，获得目标的状态更新方程的具体方法为：

步骤D1、计算跟踪波门内所有量测点迹属于真实目标航迹的概率；

k采样点时刻的第i个量测点迹z_i(k)属于目标的概率β_i(k)为：

其中，μ为虚假量测的空间密度，即单位体积的虚假量测数；V_i(k)是与k采样点时刻第i个量测点迹对应的新息，V_i′(k)为V_i(k)的转置，S(k)为k采样点时刻新息协方差，S^-1(k)为S(k)的逆，P_D为真实目标航迹不交叉时目标航迹的检出概率，P_G为真实目标航迹不交叉时量测点迹落入波门内的概率；

步骤D2、利用跟踪波门内所有量测点迹属于真实目标航迹的概率，获取目标的状态更新方程：

其中，为组合新息，V_i(k)为k采样点时刻第i个量测点迹的新息；

X_i(k|k)＝X(k|k-1)+K(k)·V_i(k)

X_i(k|k)是k采样点时刻以事件θ_i(k)为条件的目标的状态估计；X(k|k-1)是k-1采样点时刻通过状态方程得到的对k采样点时刻目标状态的预测，当i＝0时，用预测值作为估计值，即X₀(k|k)＝X(k|k-1)；K(k)是k采样点时刻的滤波增益。

本实施方式，当目标波门不相交时应用PDA算法对目标状态进行滤波更新的过程为：

首先计算波门内所有量测点迹来自于目标的概率：

用表示k时刻落入波门内点迹的集合，m(k)为落入波门的点迹数目；z^k表示直到k时刻的确认点迹的累计合集。

定义事件：

θ_i(k)表示z_i(k)是来自目标的正确量测的事件；

θ₀(k)表示k时刻没有源于目标的量测。

k时刻的第i个点迹z_i(k)与目标互联的概率为：

β_i(k)＝P{θ_i(k)|z^k}

在实际情况中，当波门内存在较多虚假量测或虚警率较高时，落入波门内的虚假量测个数近似服从参数为μV_k的泊松分布，此时互联概率为：

式中，其中，μ为虚假量测的空间密度，即单位体积的虚假量测数，以直到k时刻落入波门内的点迹数进行平均可对μ进行估计：

其中，V_k为波门体积，当n_z＝1时，

状态更新与协方差更新：

考虑波门内所有量测点迹，对每一个量测点迹赋以一定的关联概率，并进行加权求和，得到最终的目标在k时刻的状态：

式中，X_i(k|k)是以事件θ_i(k)为条件的状态估计，即：

X_i(k|k)＝X(k|k-1)+K(k)·V_i(k)

式中，V_i(k)是与该量测值对应的新息

若没有量测与目标关联，即i＝0时

X₀(k|k)＝X(k|k-1)

可得目标状态更新方程：

式中，为组合新息。

对应的协方差更新为：

式中：

P_c(k|k)＝P(k|k-1)-P(k|k-1)·K(k)·H(k)

进一步地，本实施方式中，

步骤四中，利用目标的状态更新方程与MPUKF滤波技术结合对目标方位进行跟踪的具体方法为：

步骤E1、建立目标在极坐标系下的状态方程和量测方程；

步骤E2、利用状态方程和量测方程得到k时刻目标状态的预测值X(k|k-1)和量测的预测值z(k|k-1)；

步骤E3、利用k采样点时刻的量测预测值z(k|k-1)与k采样点时刻观测得到的量测值z(k)，获得k采样点时刻的新息及其协方差；

步骤E4、利用k采样点时刻的目标状态预测值X(k|k-1)和新息，通过目标状态更新方程，得到k采样点时刻目标状态的估计值X(k|k)，实现对目标方位的跟踪。

进一步地，本实施方式中，

步骤E4中，获取k采样点时刻目标状态的估计值X(k|k)，实现对目标方位进行跟踪的具体方法为：

建立目标在极坐标系下的状态方程：

其中，为k时刻的状态向量；X(k-1)为k-1时刻的状态向量，f[·]为非线性状态函数； 为k采样点时刻的方位角，为k-1采样点时刻的方位角，为k采样点时刻的方位角变化率，采样点时刻的方位角变化率，r(k)为k采样点时刻的目标距离观测站的距离，为k采样点时刻的距离变化率；r(k-1)为k-1采样点时刻的目标距离观测站的距离，为k-1采样点时刻的距离变化率；T为量测时间间隔；

利用目标在极坐标系下的状态方程，建立量测方程：

z(k)＝H(k)X(k)+v(k)

其中，H(k)为k采样点时刻的量测矩阵，所述H(k)＝[0 0 1]；v(k)为k采样点时刻的量测噪声，是0均值高斯白噪声；

由状态方程得到k采样点时刻对目标状态的预测值X(k|k-1)及预测值X(k|k-1)和k采样点时刻的状态向量X(k)之间的协方差P(k|k-1)：

其中，ΔX_i(k|k-1)＝ξ_i(k|k-1)-X(k|k-1)，ΔX′_i(k|k-1)是ΔX_i(k|k-1)的转置；Q(k)为过程噪声的协方差矩阵；W_i为采样点估计值ξ_i(k-1|k-1)对应的权值，n_x是状态向量的维度；

ξ_i(k|k-1)＝f[ξ_i(k-1|k-1)]

是由k-1采样点时刻的采样点ξ_i(k-1|k-1)在k时刻的预测值，ξ_i(k-1|k-1)由k-1采样点时刻估计出的状态向量X(k-1|k-1)和协方差P(k-1|k-1)计算得到：

每个采样点估计值ξ_i(k-1|k-1)对应的权值W_i为：

式中：κ是尺度参数，满足n_x+κ≠0；是(n_x+κ)P(k-1|k-1)矩阵的第i行或列，n_x是状态向量的维度，P(k-1|k-1)为k-1时刻更新的目标状态的估计值和k采样点时刻的状态向量X(k)之间的协方差；

利用k采样点时刻对目标状态的预测值X(k|k-1)及预测值X(k|k-1)和k采样点时刻的状态向量X(k)之间的协方差P(k|k-1)获取k采样点时刻目标状态的估计值X(k|k)：

X(k|k)＝X(k|k-1)+K(k)·V(k)

式中：

V(k)＝z(k)-z(k|k-1)

z(k|k-1)＝H(k)·X(k|k-1)

K(k)＝P(k|k-1)·H′(k)·S^-1(k)

S(k)＝H(k)·P(k|k-1)·H′(k)+R(k)

P(k|k)＝P(k|k-1)-K(k)·S(k)·K'(k)

其中，K(k)为k采样点时刻的增益，K′(k)为K(k)的转置；S(k)为新息的协方差；z(k|k-1)是k采样点时刻的量测预测值；H′(k)为量测矩阵H(k)的转置；R(k)为k采样点时刻的量测噪声的协方差矩阵；P(k|k)为k采样点时刻更新的目标状态的估计值X(k|k)和真实值X(k)之间的协方差。

具体实施例：海试实验：结合图2至图4进行说明；

从32元均匀直线阵采集到的数据中选取25分钟时域信号，对目标跟踪750步，每步时长2秒。利用波束形成技术对时域信号进行处理，得到如图2所示的方位历程图。

然后利用检测器对目标进行检测，设置虚警概率0.1％，左右保护单元与参考单元长度均为1°。检测后判决为有目标处能量不变，判决为无目标处能量置为0。将能量值取分贝(dB)并归一化处理，得到如图3的检测结果。

将实验数据中信号的信噪比、阵元个数、阵元间距、快拍数以及目标大致方位等信息代入克拉美罗界，对方位估计克拉美罗界进行计算，以此为依据调整滤波器中量测噪声标准差，门概率P_G＝1，航迹终结流程取A＝5，跟踪结果如图4所示。

从图中结果可以看出，利用多目标跟踪方法可对上述环境中的多目标进行有效纯方位跟踪。对于能量相当的多目标航迹交叉、能量相差明显的多目标航迹交叉、目标航迹间断以及目标航迹终结等情况均可进行较好的处理。

虽然在本文中参照了特定的实施方式来描述本发明，但是应该理解的是，这些实施例仅仅是本发明的原理和应用的示例。因此应该理解的是，可以对示例性的实施例进行许多修改，并且可以设计出其他的布置，只要不偏离所附权利要求所限定的本发明的精神和范围。应该理解的是，可以通过不同于原始权利要求所描述的方式来结合不同的从属权利要求和本文中所述的特征。还可以理解的是，结合单独实施例所描述的特征可以使用在其他所述实施例中。