Search

Article

x

留言板

尊敬的读者、作者、审稿人, 关于本刊的投稿、审稿、编辑和出版的任何问题, 您可以本页添加留言。我们将尽快给您答复。谢谢您的支持!

姓名
邮箱
手机号码
标题
留言内容
验证码

Shallow sea matching field continuous tracking method based on trajectory Poisson multi-Bernoulli hybrid filter

Zhou Yu-Yuan Sun Chao Xie Lei

Zhou Yu-Yuan, Sun Chao, Xie Lei. Shallow sea matching field continuous tracking method based on trajectory Poisson multi-Bernoulli hybrid filter. Acta Phys. Sin., 2023, 72(18): 184301. doi: 10.7498/aps.72.20230124
Citation:

Shallow sea matching field continuous tracking method based on trajectory Poisson multi-Bernoulli hybrid filter

Zhou Yu-Yuan, Sun Chao, Xie Lei
PDF
HTML
Get Citation
  • In the shallow water waveguide, matched field tracking methods use the continuity of the peak position of the moving source and the disorder of pseudo-peaks on the sequential ambiguity surfaces to track the underwater source trajectory. However, owing to the dual influence of the space-time fluctuating shallow water waveguide environment and the complex sources motion scene, the existing matching field tracking methods are prone to track interruption, switches and false track phenomena, leading to discontinuous tracking results. Using the consistency between the peak position distance likelihood and the peak amplitude likelihood of sequential ambiguity surfaces, a continuous matched field tracking method is proposed based on the trajectory Poisson multi-Bernoulli mixture filter in this paper. The proposed method is applied to SWellEx-96 experimental data, and the tracking performance is measured by the linear programming metric. The results show that compared with the existing matching field tracking method and multi-target tracking method via random finite set, the proposed method achieves continuous tracking and accurate quantity estimation of moving sources trajectory. Among them, the prediction step and updating step in the trajectory space can avoid the phenomenon of trajectory interruption and switches in unvoiced periods.
      PACS:
      43.30.Wi(Passive sonar systems and algorithms, matched field processing in underwater acoustics)
      43.60.Uv(Model-based signal processing)
      43.60.-c(Acoustic signal processing)
      Corresponding author: Sun Chao, csun@nwpu.edu.cn
    • Funds: Project supported by the National Natural Science Foundation of China (Grant Nos. 12274348, 11534009).

    实际浅海波导呈时空非均匀分布, 且海底边界对声传播影响较大. 由浅海水下声源运动引起的声场扰动和沿声源轨迹的水深变化等因素易造成接收声场与拷贝场的相关性下降, 此时匹配场处理(matched field processing, MFP)[1]输出的模糊度函数出现较多伪峰, 导致运动声源位置不确定. 于是业内学者在MFP基础上陆续提出了多种匹配场跟踪(matched field tracking, MFT)方法[27], 利用模糊度函数时间序列中声源所在位置移动的连续性和伪峰所在位置的无序性, 降低声源位置的不确定性, 实现浅海水下运动声源轨迹跟踪.

    最初Bucker[2]和Fialkowski等[3]假设声源轨迹完全由起始位置确定, 通过对模糊度函数时间序列分别采取求和与平均的方式, 跟踪直线运动声源的轨迹; Maranda和Fawcett[4,5]又以声源直线运动作为约束条件, 将模糊度函数峰值数量最多的搜索路径当成声源轨迹的跟踪结果; Zala等 [6]鉴于声源运动会改变测量协方差矩阵的秩, 将各拷贝轨迹矩阵的最大特征向量与时间平均测量协方差矩阵进行匹配, 从而实现直线运动声源的轨迹和速度估计. 由于这些MFT方法均将声源直线运动当作先验假设, 以致其跟踪性能往往在实际场景中受限. 后来Tantum和Nolte [7]提出了最优不确定场跟踪, 引入马尔可夫模型以捕获声源运动的随机性, 是一种未假设声源直线运动的MFT方法, 但仅限于单声源跟踪情形. 由于目前仍缺少有关声源数量未知的MFT方法研究, 本文将针对运动形式复杂且数量未知的水下声源进行匹配场跟踪, 除模糊度函数峰值位置与声源位置不一致的影响, 还面临声源未发声和轨迹交叠时间段, 跟踪过程易出现轨迹中断、交叉混叠和虚假轨迹现象, 导致不连续的轨迹跟踪结果, 这是具有挑战性的问题.

    基于随机有限集(random finite set, RFS)[8]的多目标贝叶斯滤波方法[9]是实现未知数量目标跟踪的有效途径, 包含多种滤波器形式. RFS滤波器将多个状态变量和测量向量建模为元素数量和状态随机的RFS形式, 可有效表征多目标跟踪过程中的不确定性, 其中高斯混合概率假设密度 (Gaussian mixture-probability hypothesis density, GMPHD)滤波器[10]在水下被动跟踪场景的应用范围最广[1115]. 但由于GMPHD将目标状态RFS中所有元素视为同分布, 以致无法区分多个紧邻目标. 与GMPHD相比, 泊松多伯努利混合(Poisson multi-Bernoulli mixture, PMBM)滤波器[16]考虑了多个目标状态的差异性, 在轨迹交叠情形中也能获得较好的跟踪结果. 在PMBM基础上又发展了轨迹泊松多伯努利混合(trajectory Poisson multi-Bernoulli mixture, TPMBM)滤波器[17], TPMBM先将所有状态变量建模为轨迹状态RFS而非目标状态RFS, 再在轨迹空间内进行预测、更新和数据关联等滤波步骤, 每个时间步的滤波输出为轨迹状态后验. 相比其他RFS滤波器需将各目标状态后验添加标签才能获得各轨迹估计结果, TPMBM给出了不同的轨迹估计方式, 无需标签也能保证任意步长之间目标状态的连续性[18]. 由于水下目标辨识是一项研究难题, 难以添加声学特征标签, 因此TPMBM更适合本文考虑的水下被动跟踪场景. 但环境失配和水下声源数量未知等客观因素会导致模糊度函数中存在较高的旁瓣和较多的伪峰, TPMBM基于距离似然的数据关联方式易产生错误关联, 导致虚假轨迹出现.

    针对声源运动形式复杂且数量未知的浅海水下被动跟踪场景, 本文基于TPMBM滤波器, 利用模糊度函数峰值位置距离似然和峰值幅度似然的一致性, 提出一种匹配场连续跟踪 (trajectory Poi sson multi-Bernoulli mixture filter-amplitude-matched field tracking, TPMBM-A-MFT)方法. 该方法首先建立匹配场轨迹状态空间模型, 将模糊度函数变换为delta函数之和以便测量RFS表示, 同时将声源轨迹状态集合分成未检测和已检测两部分并由泊松和多伯努利混合参数集表示, 以准确描述水下声源数量、运动状态及测量关联的不确定性; 然后在每个时间步中于轨迹空间内依次进行预测、更新和数据关联步骤, 实现声源数量估计及其轨迹连续跟踪. 由于本文所提方法考虑了各声源轨迹状态的差异性并依据TPMBM滤波方式在轨迹空间内进行预测和更新, 可减少轨迹混叠现象出现; 每个时间步中由当前测量RFS对起始至当前时间步的轨迹状态先验更新可以避免轨迹中断现象; 此外, 本文所提方法改进了TPMBM原有的数据关联步骤, 通过结合模糊度函数中峰值位置和幅度信息, 将距离似然和幅度似然联合表示关联概率, 从而降低错误关联概率以及抑制虚假轨迹.

    其他部分安排如下: 第2节介绍TPMBM滤波器; 以该滤波器为基础, 在第3节首先建立匹配场轨迹状态空间模型, 然后推导出结合峰值位置和峰值幅度的数据关联公式, 从而给出匹配场连续跟踪方法原理和所提方法的跟踪流程; 第4节将现有匹配场跟踪方法和本文所提方法应用于SWellEx-96实验数据, 并与其他3种RFS多目标跟踪方法对比, 再由线性规划准则评估跟踪性能, 从而验证本文所提方法的优越性能; 第5节给出相关结论.

    浅海水下运动目标辐射声波经远距离传播后到达垂直阵, 将各阵元接收声压时序数据划分K$ {{\boldsymbol{y}}_{1:K}} $, 再与声场计算的拷贝声压$ {{\boldsymbol{w}}_{{\text{look}}}}(r, d) $进行匹配, 得匹配输出功率[1]:

    $$ ^{ } {B_k}(r,d) = {\left| {\frac{{{\boldsymbol{w}}_{{\text{look}}}^{\rm{H}}(r,d){{\boldsymbol{y}}_k}}}{{\left| {{{\boldsymbol{w}}_{{\text{look}}}}(r,d){{\boldsymbol{y}}_k}} \right|}}} \right|^2} = b_k^2\left( {r,d} \right),k \in [1,K] {, } $$ (1)

    式中, 上标${\rm{H}}$为共轭转置, $ b( \cdot ) $为匹配输出幅值, 距离$r$与深度$d$二维搜索空间${\mathcal{Z}}( {{{[ r\;\;d ]}^{\rm{T}}} \in {\mathcal{Z}}} )$中匹配输出功率构成模糊度表面函数(下文统称为模糊度函数), 理想条件下模糊度函数的全局极大峰值点对应目标位置, 非理想条件下模糊度函数出现较多的伪峰与高旁瓣, 导致目标位置不确定.

    MFT方法在模糊度函数时间序列上沿搜索路径累加幅值, 伪峰随步数增加而减少, 可实现运动声源轨迹跟踪. 但是, 实际浅海被动跟踪场景中目标数量通常未知, 现有MFT方法采取的轨迹搜索求和方式在实际应用中受限. 而RFS多目标滤波器在跟踪时无需假设目标数量已知, 将目标状态变量与测量信息建模为数量随机的特殊集合形式:

    $$ {{\bf{X}}_k} = \{ {\boldsymbol{x}}_k^1,{\boldsymbol{x}}_k^2, \cdots ,{\boldsymbol{x}}_k^{{N_k}}\} , $$ (2)
    $$ {{\bf{Z}}_k} = \{ {\boldsymbol{z}}_k^1,{\boldsymbol{z}}_k^2, \cdots ,{\boldsymbol{z}}_k^{{M_k}}\} , $$ (3)

    式中, $ {\boldsymbol{x}}_k^1, {\boldsymbol{x}}_k^2, \cdots , {\boldsymbol{x}}_k^{{N_k}} $是第k时间步${N_k}$个目标状态变量, $ {\boldsymbol{z}}_k^1, {\boldsymbol{z}}_k^2, \cdots , {\boldsymbol{z}}_k^{{M_k}} $${M_k}$个测量向量.

    TPMBM滤波器独特之处在于, ${{\bf{X}}_k}$${N_k}$个轨迹状态变量组成, $ {{\bf{X}}_k} = \{ {{X}}_k^1, {{X}}_k^2, \cdots , {{X}}_k^{{N_k}}\} $, 任意单个轨迹状态变量$ {X_k} $表示为

    $$ {{{X}}_k} = \left( {\begin{array}{*{20}{c}} {\beta ,}&{\varepsilon ,}&{{{\boldsymbol{x}}_{\beta :\varepsilon }}} \end{array}} \right) , $$ (4)
    $$ {{\boldsymbol{x}}_{\beta :\varepsilon }} = {[{{\boldsymbol{x}}_\beta };{{\boldsymbol{x}}_{\beta + 1}};\cdots;{{\boldsymbol{x}}_{k'}};\cdots;{{\boldsymbol{x}}_\varepsilon }]^{\rm{T}}} \in {{\mathbb{R}}^{{n_x}L}} , $$ (5)

    式中, $\beta $$\varepsilon (\varepsilon \leqslant k)$分别为初始与最新目标状态对应时间步, 上标T表示转置, $ {{\boldsymbol{x}}_{\beta :\varepsilon }} $$L = \varepsilon - \beta + 1$个时间步的目标状态变量${{\boldsymbol{x}}_{k'}} \in {\mathbb{R}^{{n_x}}}$$(k' \in [\beta , \varepsilon ])$组成的目标状态序列, ${n_x}$是目标状态$ {{\boldsymbol{x}}_{k'}} $维数. TPMBM的轨迹状态RFS建模方式可视为“隐式标签”, 在匹配场跟踪过程中确保多个目标状态在任意步长之间关联连续性.

    另外, TPMBM将$ {{\bf{Z}}_k} = {{\bf{O}}_k} \lus {{\bf{C}}_k} $分为目标产生的测量RFS $ {{\bf{O}}_k} $与杂波产生的测量RFS $ {{\bf{C}}_k} $, $ \lus $为并集, $ {{\bf{X}}_k} = {\bf{X}}_k^{\rm{d}} \lus {\bf{X}}_k^{\rm{u}} $分为已检测的状态RFS $ {\bf{X}}_k^{\rm{d}} $与未检测的状态RFS $ {\bf{X}}_k^{\rm{u}} $, 能考虑到匹配场跟踪过程中各测量信息、各运动声源状态的差异性, 如非理想条件下模糊度函数中多数伪峰是与声源位置无关的杂波测量, 声源间歇发声会降低检测率并增加声源位置不确定性. TPMBM将$ {{\bf{X}}_k} $的概率密度函数(probability density function, PDF) $ p({{\bf{X}}_k}) $分别由泊松多伯努利混合与泊松参数集表示以度量$ {\bf{X}}_k^{\rm{d}} $$ {\bf{X}}_k^{\rm{u}} $不确定性:

    $$ p({{\bf{X}}_k}) \triangleq \left\{ w_k^h, \left\{ r_k^{i,h},p({\bf{X}}_k^{{\rm{d}},i,h}) \right\} _{i = 1}^{N_k^h}\right\} _{h = 1}^{{{\mathcal{H}}_k}} \lus \lambda ({\bf{X}}_k^{\rm{u}}), $$ (6)

    式中, $ {\bf{X}}_k^{\rm{d}} $$ {\bf{X}}_k^{\rm{u}} $存在不确定性分别由伯努利存在概率$ r_k^{i, h} $和泊松强度$ \lambda ({\bf{X}}_k^{\rm{u}}) $描述; 轨迹状态的不确定性由$ p({\bf{X}}) $描述; 共${{\mathcal{H}}_k}$种全局关联假设描述$ {\bf{X}}_k^{\rm{d}} $$ {{\bf{Z}}_k} $之间关联不确定性, 第h个全局关联假设成立概率为$ w_k^h $, $ N_k^h $$ {\bf{X}}_k^{\rm{d}} $中轨迹状态的预期数量.

    RFS多目标滤波器均是基于贝叶斯滤波理论传递$ p({{\bf{X}}_k}) $参数集, 实现状态变量及其数量估计. 而TPMBM滤波过程的独特之处在于, 参数集的预测、更新和数据关联步骤是在轨迹空间内执行, 并非状态空间.

    2.2.1   $ p({{\bf{X}}_k}) $参数集预测

    k – 1时刻$ p({{\bf{X}}_{k - 1}}) $参数集已知, 单条轨迹状态$ {X_{k - 1}} = \left( {{\beta _0}, k - 1, {{\boldsymbol{x}}_{{\beta _0}:k - 1}}} \right) $ PDF为

    $$ \begin{split}p({{\bf{X}}_{k - 1}}) =\;& \delta (\beta - {\beta _0})\delta (\varepsilon - k + 1){\mathcal{N}}({{\boldsymbol{x}}_{\beta :k - 1}};\\ &{{\boldsymbol{\mu}} _{\beta :k - 1}},{{\boldsymbol{P}}_{\beta :k - 1}}) ,\end{split} $$ (7)
    $$ {{\boldsymbol{\mu}} _{\beta :k - 1}} = {[{{\boldsymbol{\mu}} _\beta };{{\boldsymbol{\mu}} _{\beta + 1}}; \cdots {{\boldsymbol{\mu}} _{k'}}; \cdots ;{{\boldsymbol{\mu}} _{k - 1}}]^{\rm{T}}} \in {{\mathbb{R}}^{{n_x}L}} , $$ (8)
    $$\begin{split} {{\boldsymbol{P}}_{\beta :k - 1}} =\;& {\rm{diag}}([{{\boldsymbol{P}}_\beta },{{\boldsymbol{P}}_{\beta + 1}}, \cdots {{\boldsymbol{P}}_{k'}}, \cdots ,{{\boldsymbol{P}}_{k - 1}}]) \\ &\in {{\mathbb{R}}^{{n_x}L \times {n_x}L}},\end{split} $$ (9)

    式中, $ \delta $为delta函数, $ {\mathcal{N}} $表示高斯分布, 轨迹状态的均值$ {{\boldsymbol{\mu}} _{\beta :k - 1}} $与协方差$ {{\boldsymbol{P}}_{\beta :k - 1}} $分别由$ L = k - {\beta _0} $个时间步上目标状态的均值${{\boldsymbol{\mu}} _{k'}} \in {{\mathbb{R}}^{{n_x}}}$与协方差$ {{\boldsymbol{P}}_{k'}} \in {{\mathbb{R}}^{{n_x} \times {n_x}}} $组成.

    TPMBM在轨迹空间$ {{\mathbb{R}}^{{n_x}L}} $$ {{\mathbb{R}}^{{n_x}\left( {L + 1} \right)}} $中进行预测以获得$ p({{\bf{X}}_k}) $先验参数集, 其中, 第k时间步的任意轨迹状态先验均值$ {{\boldsymbol{\mu}} _{\beta :k|\beta :k - 1}} $和协方差$ {{\boldsymbol{P}}_{\beta :k|\beta :k - 1}} $分别为

    $$ {{\boldsymbol{\mu}} _{\beta :k|\beta :k - 1}} = {\left[ {{{\boldsymbol{\mu}} _{\beta :k - 1}},{{\boldsymbol{\mu }}_{k|k - 1}}} \right]^{\rm{T}}} \in {{\mathbb{R}}^{{n_x}(L + 1)}}, $$ (10)
    $$ \begin{split} {{\boldsymbol{P}}_{\beta :k|\beta :k - 1}} =& \left[ {\begin{array}{*{20}{c}} {{{\boldsymbol{P}}_{\beta :k - 1}}}&{{{\left( {{{\boldsymbol{F}}^{\text{t}}}{{\boldsymbol{P}}_{\beta :k - 1}}} \right)}^{\rm{T}}}} \\ {{{\boldsymbol{F}}^{\rm{t}} }{{\boldsymbol{P}}_{\beta :k - 1}}}&{{{\boldsymbol{P}}_{k|k - 1}}} \end{array}} \right] \\ &\in {\mathbb{R}^{{n_x}(L + 1) \times {n_x}(L + 1)}},\end{split} $$ (11)

    式中, $ {{\boldsymbol{\mu}} _{k|k - 1}} $$ {{\boldsymbol{P}}_{k|k - 1}} $为目标状态先验均值与协方差, $ {{\boldsymbol{F}}^{\text{t}}} = [{\mathbb{O}^{{n_x} \times {n_x}(L - 1)}}, {\boldsymbol{F}}] $是轨迹状态转移矩阵, $ \mathbb{O} $为零矩阵, $ {\boldsymbol{F}} \in {\mathbb{R}^{{n_x} \times {n_x}}} $是目标状态转移矩阵, 状态噪声协方差为$ {\boldsymbol{Q}} $,

    $$ {{\boldsymbol{\mu}} _{k|k - 1}} = {\boldsymbol{F}}{{\boldsymbol{\mu}} _{k - 1}} \in {\mathbb{R}^{{n_x}}}, $$ (12)
    $$ {{\boldsymbol{P}}_{k|k - 1}} = {\boldsymbol{F}}{{\boldsymbol{P}}_{k - 1}}{{\boldsymbol{F}}^{\rm{T}}} + {\boldsymbol{Q}} \in {\mathbb{R}^{{n_x} \times {n_x}}}. $$ (13)

    (10), (11)式与(12), (13)式分别为轨迹空间与目标状态空间预测示例, 对比可知, (10)式与(11)式不仅将k –1时刻目标状态$ {{\boldsymbol{\mu}} _{k - 1}} $$ {{\boldsymbol{P}}_{k - 1}} $作为先验以预测第k时间步目标状态$ {{\boldsymbol{\mu}} _{k|k - 1}} $$ {{\boldsymbol{P}}_{k|k - 1}} $, 还在轨迹状态先验$ {{\boldsymbol{\mu}} _{\beta :k|\beta :k - 1}} $中保留以前时刻的目标状态, 可确保任意步长间各目标状态连续性, 从而减少匹配场跟踪过程中出现轨迹混叠次数.

    由此, 任意第h个全局关联假设中第i条已检测轨迹状态先验PDF可表示为

    $$ p({\bf{X}}_{k|k - 1}^{{\rm{d}} ,i,h}) = \delta (\beta - {\beta _0}){\mathcal{N}}({{\boldsymbol{x}}_{\beta :k}};{\boldsymbol{\mu}} _{\beta :k|\beta :k - 1}^{i,h},{\boldsymbol{P}}_{\beta :k|\beta :k - 1}^{i,h}) , $$ (14)

    存在概率$ r_{k|k - 1}^{i, h} = r_{k - 1 |k - 1}^{i, h}{P^{\text{S}}} $, 存活概率$ {P^{\text{S}}} $; 其余参数除$ \lambda ({\bf X}_{k|k - 1}^{\rm{u}} ) $外均与k –1时刻一致,

    $$ \begin{split} &\lambda ({\bf{X}}_{k|k - 1}^{\rm{u}} ) = w_k^{\text{b}}{\mathcal{N}}({{\boldsymbol{x}}_k};{\boldsymbol{\mu}} _k^{\text{b}},{\boldsymbol{P}}_k^{\text{b}}) \\ & ~~~~ + \sum\limits_{i' = 1}^{N_{k - 1|k - 1}^{\rm{u}} } w_{k - 1}^{{\rm{u}} ,i'} \delta (\beta - {\beta _0}) \\ & ~~~~ \times {\mathcal{N}}({{\boldsymbol{x}}_{\beta :k}};{\boldsymbol{\mu}} _{\beta :k|\beta :k - 1}^{{\rm{u}} ,i'},{\boldsymbol{P}}_{\beta :k|\beta :k - 1}^{{\rm{u}} ,i'}){P^{\text{S}}} ,\end{split} $$ (15)

    式中, $ w_{k - 1}^{{\rm{u}} , i'} $, $ {\boldsymbol{\mu}} _k^{\text{b}} $$ {\boldsymbol{P}}_k^{\text{b}} $k时刻$ {\bf{X}}_k^{\rm{u}} $中新生轨迹的出现概率、状态均值与协方差.

    2.2.2   $ p({{\bf{X}}_k}) $参数集更新

    TPMBM在轨迹空间$ {\mathbb{R}^{{n_x}(L + 1)}} $中由测量RFS对$ p({{\bf{X}}_k}) $先验参数集更新, 获得$ p({{\bf{X}}_k}) $后验参数集. 其中, 任意轨迹状态后验均值$ {{\boldsymbol{\mu}} _{\beta :k|\beta :k}} $和协方差$ {{\boldsymbol{P}}_{\beta :k|\beta :k}} $分别为

    $$ {{\boldsymbol{\mu}} _{\beta :k|\beta :k}} = {{\boldsymbol{\mu}} _{\beta :k|\beta :k - 1}} + {{\boldsymbol{K}}_{\beta :k}}({{\boldsymbol{z}}_k} - {\hat {\boldsymbol{z}}_k}) \in {\mathbb{R}^{{n_x}(L + 1)}}, $$ (16)
    $$ \begin{split}{{\boldsymbol{P}}_{\beta :k|\beta :k}} =\;& {{\boldsymbol{P}}_{\beta :k|\beta :k - 1}} - {{\boldsymbol{K}}_{\beta :k}}{{\boldsymbol{H}}^{\text{t}}}{P{\boldsymbol{}}_{\beta :k|\beta :k - 1}} \\ &\in {\mathbb{R}^{{n_x}(L + 1) \times {n_x}(L + 1)}},\end{split} $$ (17)

    式中$ {{\boldsymbol{K}}_{\beta :k}} = {{\boldsymbol{P}}_{\beta :k|\beta :k - 1}}{({{\boldsymbol{H}}^{\text{t}}})^{\rm{T}}}{{\boldsymbol{S}}_k}^{ - 1} $, $ {\hat {\boldsymbol{z}}_k} = {\boldsymbol{H}}{{\boldsymbol{\mu}} _{k - 1 |k}} $, $ {{\boldsymbol{S}}_k} = {{\boldsymbol{H}}^{\text{t}}}{{\boldsymbol{P}}_{\beta :k|\beta :k - 1}}{({{\boldsymbol{H}}^{\text{t}}})^{\rm{T}}} + {\boldsymbol{R}} $, 轨迹测量矩阵$ {{\boldsymbol{H}}^{\text{t}}} = [{\mathbb{O}^{{n_z} \times {n_x}L}},{\boldsymbol{ H}}] $, $ {\boldsymbol{H}} \in {\mathbb{R}^{{n_z} \times {n_x}}} $是目标测量矩阵, $ {\boldsymbol{R}} $是测量噪声协方差. 从(16)式乘积项$ {{\boldsymbol{K}}_{\beta :k}}({{\boldsymbol{z}}_k} - {\hat {\boldsymbol{z}}_k}) \in {\mathbb{R}^{{n_x}(L + 1)}} $维数可知, 当前测量会对初始至当前时刻的目标状态先验序列矫正, 可减少匹配场跟踪过程中间歇发声引起的轨迹中断现象.

    由此, 若第h个先验全局关联假设中第i条已检测轨迹状态先验$ {\bf{X}}_{k|k - 1}^{{\rm{d}}, i, h} $与第j个测量向量${\boldsymbol{z}}_k^j$关联, 其存在概率$ r_{k|k}^{j \ne 0, i, h} = 1 $, 轨迹状态后验PDF可表示为

    $$ {p^{j \ne 0}}({\bf{X}}_{k|k}^{{\rm{d}},i,h}) = \delta (\beta - {\beta _0}){\mathcal{N}}({{\boldsymbol{x}}_{\beta :k}};{\boldsymbol{\mu}} _{\beta :k|\beta :k}^{i,h},{\boldsymbol{P}}_{\beta :k|\beta :k}^{i,h}) . $$ (18)

    $ {\bf{X}}_{k|k - 1}^{{\rm{d}}, i, h} $与任意${z_k}$均无关联, 存在概率$ r_{k|k}^{j = 0, i, h} = r_{k|k - 1}^{i, h}( {1 - {P^{\text{D}}}})/( {1 - r_{k|k - 1}^{i, h}{P^{\text{D}}}} ) $, $ {P^{\text{D}}} $是检测概率, PDF$ {p^{j = 0}}(X_{k|k}^{{\rm{d}}, i, h}) = p(X_{k|k - 1}^{{\rm{d}}, i, h}) $; 若未检测轨迹$ {\bf{X}}_{k|k - 1}^{\rm{u}} $${\boldsymbol{z}}_k^j$关联, 将在$ {\bf{X}}_k^{\rm{d}} $中生成新检测轨迹状态后验$ {\bf{X}}_{k|k}^{i = 0, h} $, 其存在概率和PDF分别为

    $$ r_{k|k}^{j,i = 0} = \frac{{{P^{\text{D}}}\sum\limits_{i' = 1}^{N_{k|k - 1}^{\rm{u}}} {w_{k|k - 1}^{{\rm{u}},i'}{\mathcal{N}}({\boldsymbol{z}}_k^j;{\boldsymbol{H}}{\boldsymbol{\mu}} _{k|k - 1}^{{\rm{u}},i'},{\boldsymbol{S}}_k^{u,i'})} }}{{{\lambda ^{\text{c}}} + {P^{\text{D}}}\sum\limits_{i' = 1}^{N_{k|k - 1}^{\rm{u}}} {w_{k|k - 1}^{{\rm{u}},i'}{\mathcal{N}}({\boldsymbol{z}}_k^j;{\boldsymbol{H}}{\boldsymbol{\mu }}_{k|k - 1}^{{\rm{u}},i'},{\boldsymbol{S}}_k^{{\rm{u}},i'})} }} , $$ (19)
    $$ {p^j}({\bf{X}}_{k|k}^{i = 0,h}) \approx \delta (\beta - \beta _0^{{i_{\max }}}){\mathcal{N}}({{\boldsymbol{x}}_{\beta :k}};{\boldsymbol{\mu}} _{\beta :k|\beta :k}^{{i_{\max }}},{\boldsymbol{P}}_{\beta :k|\beta :k}^{{i_{\max }}}) , $$ (20)

    式中, $ {i_{\max }} = \arg \max \{ w_{k|k - 1}^{{\rm{u}}, i'}{\mathcal{N}}({\boldsymbol{z}}_k^j;{\boldsymbol{H\mu}} _{k|k - 1}^{{\rm{u}}, i'}, {\boldsymbol{S}}_k^{{\rm{u}}, i'})\} $, $ {\lambda ^{\text{c}}} $是杂波密度; $ {\bf{X}}_{k|k - 1}^{\rm{u}} $与任意${{\boldsymbol{z}}_k}$均无关联情形由$ \lambda ({\bf{X}}_{k|k}^{\rm{u}}) = (1 - {P^D})\lambda (X_{k|k - 1}^{\rm{u}}) $表示.

    2.2.3   数据关联

    h个先验全局关联假设中, $ {\bf{X}}_{k|k - 1}^{\rm{d}} $, $ {\bf{X}}_{k|k - 1}^{\rm{u}} $$ {\boldsymbol{z}}_k^j $关联的概率分别为

    $$ \begin{split}\tilde w_k^{j,i \ne 0} = \;&r_{k|k - 1}^{i,h}{P^{\text{D}}} \dfrac{{\mathcal{N}}({\boldsymbol{z}}_k^j;{\boldsymbol{H\mu}} _{k|k - 1}^{i,h},{\boldsymbol{S}}_k^{i,h})}{ \left( {1 - r_{k|k - 1}^{i,h}{P^{\text{D}}}} \right)},\end{split} $$ (21)
    $$ \begin{split} \tilde w_k^{j,i = 0} =\;& {P^{\text{D}}}\sum\limits_{i' = 1}^{N_{k|k - 1}^{\rm{u}}} {w_{k|k - 1}^{{\rm{u}},i'}{\mathcal{N}}({\boldsymbol{z}}_k^j;{\boldsymbol{H\mu}} _{k|k - 1}^{{\rm{u}},i'},{\boldsymbol{S}}_k^{{\rm{u}},i'})} \\ &+ {\lambda ^{\text{c}}},\\[-15pt]\end{split} $$ (22)

    式中, $ {\mathcal{N}}({\boldsymbol{z}}_k^j;{\boldsymbol{H}}{{\boldsymbol{\mu}} _{k|k - 1}}, {{\boldsymbol{S}}_k}) $是状态先验$ {{\boldsymbol{\mu}} _{k|k - 1}} $与测量$ {\boldsymbol{z}}_k^j $的距离似然值. 距离似然矩阵为

    $$ {\boldsymbol{W}}_{{\rm{dis}}}^h = \left[ {\begin{array}{*{20}{c}} {{{\tilde w}^{1,1}}}&{{{\tilde w}^{1,2}}}& \cdots &{{{\tilde w}^{1,N_{k|k - 1}^h}}}&{{{\tilde w}^{1,0}}}& \cdots &0 \\ \vdots & \vdots & \ddots & \vdots & \vdots & \ddots & \vdots \\ {{{\tilde w}^{{M_k},1}}}&{{{\tilde w}^{{M_k},2}}}& \cdots &{{{\tilde w}^{{M_k},N_{k|k - 1}^h}}}&0& \cdots &{{{\tilde w}^{{M_k},0}}} \end{array}} \right] . $$ (23)

    各元素由(21)式与(22)式组成, 描述各状态先验与各测量关联概率. 在任意测量最多与一个状态先验关联的约束条件下, TPMBM依据$ {\boldsymbol{W}}_{{\rm{dis}}}^h $以求解最优后验关联假设:

    $$ \begin{split} &{\boldsymbol{A}}_k^{\dagger ,h} = \mathop {\arg \min }\limits_{A_k^h} {\text{tr}}\left[ {{{ - {\rm{lg}}(}}{\boldsymbol{W}}_{{\rm{dis}}}^h{\text{)}} \times {{({\boldsymbol{A}}_k^h)}^{\rm{T}}}} \right]{,} \\ &{\text{s}}{\text{.t }}\left\{ \begin{aligned} &a_k^{j,i} \in \{ 0,1\} ,~~\forall~ i,j, \\ &\sum\limits_i {a_k^{j,i} = 1,~~\forall ~j = 1,2, \cdots ,{M_k},} \\ &\sum\limits_j {a_k^{j,i} \in \{ 0,1\} ,~~\forall~ i = 1,2, \cdots ,N_{k|k - 1}^h + {M_k}. } \\ \end{aligned} \right. \end{split} $$ (24)

    式中$ {\boldsymbol{A}}_k^h $元素$ a_k^{j, i} $非0即1, ${{\rm{t}}{\rm{r}}} [ \cdot ]$是矩阵的迹, $ w_{k|k}^h = {\text{tr}}[{\boldsymbol{W}}_{{\rm{dis}}}^h{({\boldsymbol{A}}_k^{\dagger , h})^{\rm{T}}}] $归一化后得后验全局关联假设概率$ w_k^h, h = 1, 2, \cdots , {{\mathcal{H}}_k} $. 若$ {h^\dagger } = \arg \max \{ w_k^h\} $, k时刻轨迹状态估计结果:

    $$ {\hat{\bf{ X}}_k} = \{{\boldsymbol{ \mu}} _{\beta :k|\beta :k}^{i,{h^\dagger }}:r_{k|k}^{i,{h^\dagger }} \geqslant \varGamma \} _{i = 1}^{{{\hat N}^{{h^\dagger }}}} \lus \{ {\boldsymbol{\mu}} _{\beta :k|\beta :\varepsilon <k}\} ,$$ (25)

    式中, ${\boldsymbol{\mu }}_{\beta :k|\beta :k}^i$是存在概率$ r_{k|k}^{i, {h^\dagger }} \geqslant \varGamma $阈值的轨迹状态均值, 不等式成立个数${\hat N^{{h^\dagger }}}$k时刻仍存活的轨迹数量估计值, $\{ {{\boldsymbol{\mu}} _{\beta :\varepsilon \lt k|\beta :k}}\} $k时刻前已结束的轨迹状态均值集合.

    3.1.1   模糊度函数变换

    本文所提方法在各模糊度函数$ {B_{1:K}}(r, d) $上依次保留${M_{1:K}}$个最高局部极大峰值作为测量信息, 其中, 峰值个数${M_{1:K}}$服从泊松分布$ {\rm{P}}\left( {{M_{1:K}};{\lambda ^z}} \right) $, 任意保留的峰值位置${{\boldsymbol{z}}_k} = {\left[ {{z_{k, r}}}\;\;{{z_{k, r}}} \right]^{\rm{T}}} \in {\mathcal{Z}}$, 且搜索区间${\mathcal{Z}}$中除$z{\boldsymbol{}}_k^m\left( {m \in [1, {M_k}], k \in [1, K]} \right)$外其余位置处幅值均设为零. 变换后模糊度函数$ \tilde B({{\boldsymbol{z}}_k}) $可表示以目标测量向量${\boldsymbol{o}}_k^j$和杂波测量向量${\boldsymbol{c}}_k^{j'}$为自变量的delta函数之和形式:

    $$ \tilde B({{\boldsymbol{z}}_k}) = \sum\limits_{j = 1}^{{N_k}} {{b^2}({{\boldsymbol{z}}_k})\delta ({{\boldsymbol{z}}_k}{{ - }}{\boldsymbol{o}}_k^j)} + \sum\limits_{j' = 1}^{N_k^c} {{b^2}({{\boldsymbol{z}}_k})\delta ({z_k} - {\boldsymbol{c}}_k^{j'})} , $$ (26)

    式中, $ {N_k} + N_k^c = {M_k} $, $ {\boldsymbol{o}}_k^j \in {{\bf{O}}_k} $与目标位置有关, 其映射关系由测量方程后续给出, 而$ {\boldsymbol{c}}_k^{j'} \in {{\bf{C}}_k} $与目标位置无关. 对比现有MFT在模糊度函数时间序列上搜索轨迹求和方式, (26)式变换成测量RFS形式, 既考虑杂波存在可能性, 又去除了模糊度函数中冗余信息, 避免在多目标跟踪中引入更大测量关联不确定性, 减少计算量.

    3.1.2   状态方程和测量方程

    依据2.1节状态变量建模方式, 各水下运动声源轨迹状态$ {{\bf{X}}_k} $由PMBM RFS表示. 其中, 目标状态${{\boldsymbol{x}}_k}$由声源距离$ {x_{k, r}} $、距离维速度$ {\dot x_{k, r}} $、声源深度$ {x_{k, d}} $和深度维速度$ {\dot x_{k, d}} $组成, 维数${n_x}$=4, ${{\boldsymbol{x}}_k}$在相邻步长间的转移关系由状态方程给出:

    $$ \begin{split} &{{\boldsymbol{x}}_k} = \left[ {\begin{array}{*{20}{c}} {{x_{k,r}}} \\ {{{\dot x}_{k,r}}} \\ {{x_{k,d}}} \\ {{{\dot x}_{k,d}}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} 1&T&0&0 \\ 0&1&0&0 \\ 0&0&1&T \\ 0&0&0&1 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{x_{k - 1,r}}} \\ {{{\dot x}_{k - 1,r}}} \\ {{x_{k - 1,d}}} \\ {{{\dot x}_{k - 1,d}}} \end{array}} \right]\\ &~~ + \left[ {\begin{array}{*{20}{c}} {{T^2}/2}&0 \\ T&0 \\ 0&{{T^2}/2} \\ 0&T \end{array}} \right] \left[ {\begin{array}{*{20}{c}} {{n_{k,r}}} \\ {{n_{k,d}}} \end{array}} \right] = {\boldsymbol{F}}{{\boldsymbol{x}}_{k - 1}} + {\boldsymbol{G}}{{\boldsymbol{n}}_k} , \end{split}$$ (27)

    式中T是时间间隔; $ {\boldsymbol{G}} $是状态噪声驱动矩阵; 状态噪声${{\boldsymbol{n}}_k}$服从均值${[0, 0]^{\rm{T}}}$、协方差${\boldsymbol{Q}} = {{\boldsymbol{G}}^{\rm{T}}}{{\rm{diag}}} (q_r^2, q_d^2){\boldsymbol{G}}$的高斯分布. ${{\boldsymbol{x}}_k}$${{\boldsymbol{o}}_k}$的映射关系由测量方程给出:

    $$ \begin{split}{{\boldsymbol{o}}_k} = \;&\left[ {\begin{array}{*{20}{c}} {{o_{k,r}}} \\ {{o_{k,d}}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} 1&0&0&0 \\ 0&0&1&0 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{x_{k,r}}} \\ {{{\dot x}_{k,r}}} \\ {{x_{k,d}}} \\ {{{\dot x}_{k,d}}} \end{array}} \right] \\ &+ \left[ {\begin{array}{*{20}{c}} {{v_{k,r}}} \\ {{v_{k,d}}} \end{array}} \right] = {\boldsymbol{H}}{{\boldsymbol{x}}_k} + {{\boldsymbol{v}}_k} ,\end{split} $$ (28)

    式中, 测量噪声$ {\boldsymbol{v_k}} $服从均值$[0, 0]^{\rm{T}}$、协方差$ {\boldsymbol{R}} = {{\rm{diag}}} (\sigma _r^2, \sigma _d^2) $的高斯分布. 杂波测量${{\boldsymbol{c}}_k}$${{\boldsymbol{x}}_k}$无关, 在${\mathcal{Z}}$空间均匀分布, 杂波密度$ {\lambda ^{\text{c}}} = S({\mathcal{Z}})/N_k^{\text{c}} $, $S( \cdot )$是空间面积.

    本节基于TPMBM建模方式, 建立匹配场轨迹状态空间模型. 将模糊度函数时间序列变换为delta函数之和以便测量RFS表示, 水下声源运动状态建模为PMBM RFS $ {{\bf{X}}_k} $以便描述各轨迹状态及其数量的不确定性, 依据(27)式与(28)式所示状态方程和测量方程, 在轨迹空间执行TPMBM滤波过程, 可实现匹配场跟踪. 但需注意, TPMBM的数据关联步骤以状态先验与测量的距离似然为关联依据, 在匹配场跟踪过程中易将模糊度函数中伪峰与声源位置关联而出现虚假轨迹. 为抑制杂波的影响, 本文将在3.2节提出结合峰值位置与峰值幅度的数据关联步骤.

    3.2.1   模糊度函数峰值幅度模型

    归一化模糊度函数幅度PDF $ p(b({{\boldsymbol{z}}_k})) $可表示为方差$ \sigma _b^2 = 1 $的瑞利分布[19]:

    $$ p(b({{\boldsymbol{z}}_k})) = \frac{{b({{\boldsymbol{z}}_k})}}{{\sigma _b^2}}\exp \left( { - \frac{{{b^2}({{\boldsymbol{z}}_k})}}{{2\sigma _b^2}}} \right) . $$ (29)

    当模糊度函数变换为$ \tilde B({{\boldsymbol{z}}_k}) $后, 非零峰值幅度$ b({\boldsymbol{z}}_k^j) $需高于门限幅度获得, 由Shalom的瑞利目标振幅波动模型[20]可知, $ b({\boldsymbol{z}}_k^j) $是目标幅度和杂波幅度的概率可表示为

    $$ {p}^{\circ}(b({{\boldsymbol{z}}}_{k}^{j})|S)=\frac{b({{\boldsymbol{z}}}_{k}^{j})}{1+S}{\rm{exp}}\bigg[\frac{{B}^{\text{c}}-{b}^{2}({{\boldsymbol{z}}}_{k}^{j})}{2(1+S)}\bigg], $$ (30)
    $$ {p}^{\text{c}}\left(b({{\boldsymbol{z}}}_{k}^{j})\right)=b({{\boldsymbol{z}}}_{k}^{j}){\rm{exp}}\bigg(\frac{{B}^{\text{c}}-{b}^{2}({{\boldsymbol{z}}}_{k}^{j})}{2}\bigg), $$ (31)

    式中, 参数S是模糊度函数中所有声源位置处幅度均方值与模糊度函数平均功率之比, $ {B^{\text{c}}} = \mathop {\min }\limits_{\forall k, j} ( {\tilde B({\boldsymbol{z}}_k^j)} t) $是模糊度函数时间序列变换时最低门限, (30)式和(31)式分别为目标与杂波的幅度似然函数. 由于声源位置是未知变量, 参数S及其相关的目标幅度似然函数无法计算, 只能通过估计$S \in \left[ {{d_1}, {d_2}} \right]$取值范围来近似, 其中,

    $$ {d_1} = \mathop {\min }\limits_k \left( {\mathop {\min }\limits_{\forall j \in {M_k}} \left( {\tilde B({\boldsymbol{z}}_k^j)} \right)} \right)/{E_N} , $$ (32)
    $$ {d_2} = \mathop {\max }\limits_k \left( {\mathop {{\text{max}}}\limits_{\forall j \in {M_k}} \left( {\tilde B({\boldsymbol{z}}_k^j)} \right)} \right)/{E_N} . $$ (33)

    这里模糊度函数平均功率${E_N} = \mathop {{\text{mean}}}\limits_{k, \forall {\boldsymbol{z}} \in {\mathcal{Z}}} ( {B({\boldsymbol{z}}_k^j)} )$. 目标幅度似然函数可近似表示为[21]

    $$ \begin{split} &{p^{\circ}}\left( {b({\boldsymbol{z}}_k^j)} \right) \\ \approx\;& \frac{{2\left\{ {\exp \left[ {\dfrac{{ - {b^2}({\boldsymbol{z}}_k^j)}}{{2(1 + {d_2})}}} \right] - \exp \left[ {\dfrac{{ - {b^2}({\boldsymbol{z}}_k^j)}}{{2(1 + {d_1})}}} \right]} \right\}}}{{b({\boldsymbol{z}}_k^j)\left[ {\ln (1 + {d_2}) - \ln (1 + {d_1})} \right]}}.\end{split} $$ (34)

    (31)式和(34)式分别描述了第k时间步测量集合$ {{\bf{Z}}_k} $中任意元素$ {\boldsymbol{z}}_k^j $是杂波测量或是目标测量的可能性, 可作为匹配场跟踪过程中区分杂波和目标测量的重要依据.

    3.2.2   距离似然和幅度似然实现数据关联

    假定目标幅度与其运动状态相互独立, 首先将3.2.1节(34)式与(31)式所示目标幅度似然与杂波幅度似然的比值$ {\overset{\lower0.5 em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w}{}^j} = {p^{\circ}}( {b({\boldsymbol{z}}_k^j)} )/{p^{\text{c}}}( {b({\boldsymbol{z}}_k^j)} ) $构成幅度似然矩阵:

    $$ {{\boldsymbol{W}}_{{\rm{amp}}}} = \left[ {\begin{array}{*{20}{c}} {{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w} }{}^1}}&{{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w}}{}^1}}& \cdots &{{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w} }{}^1}}&{{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w} }{}^1}}& \cdots &0 \\ \vdots & \vdots & \ddots & \vdots & \vdots & \ddots & \vdots \\ {{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w} }{}^{{M_k}}}}&{{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w} }{}^{{M_k}}}}& \cdots &{{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w} }{}^{{M_k}}}}&0& \cdots &{{{\overset{\lower0.5em\hbox{$\smash{\scriptscriptstyle\frown}$}}{w} }{}^{{M_k}}}} \end{array}} \right] , $$ (35)

    再由$ {{\boldsymbol{W}}_{{\rm{amp}}}} $与(23)式所示距离似然矩阵$ {\boldsymbol{W}}_{{\rm{dis}}}^h $联合表示各状态先验与测量之间的关联概率, 以(24)式约束条件不变为前提, 将原有最优关联假设求解方式更改为

    $$ {\boldsymbol{A}}_k^{\dagger ,h} = \mathop {\arg \min }\limits_{{\boldsymbol{A}}_k^h} {\text{tr}}\left[ {{{ - [{\rm{lg}}(}}{\boldsymbol{W}}_{{\rm{dis}}}^h) +{\text{ {\rm{lg}}}}({{\boldsymbol{W}}_{{\rm{amp}}}}){\text{}}]{\boldsymbol{A}}_k^{\rm{T}}} \right] , $$ (36)

    式中, $ {\boldsymbol{W}}_{{\rm{dis}}}^h $各元素大小与模糊度函数峰值位置$ {\boldsymbol{z}}_k^j $相关, 取决于各峰值位置$ {\boldsymbol{z}}_k^j $与各状态先验$ {\boldsymbol{\mu}} _{k|k - 1}^i $的距离似然值$ {\mathcal{N}}({\boldsymbol{z}}_k^j;{\boldsymbol{H\mu}} _{k|k - 1}^i, {{\boldsymbol{S}}_k}) = ( {\boldsymbol{z}}_k^j - {\boldsymbol{H\mu}} _{k|k - 1}^i )^{\rm{T}} \cdot {{\boldsymbol{S}}_k} ( {{\boldsymbol{z}}_k^j -{\boldsymbol{ H\mu }}_{k|k - 1}^i} ) $; $ {{\boldsymbol{W}}_{{\rm{amp}}}} $各元素大小与模糊度函数各峰值幅度$ b({\boldsymbol{z}}_k^j) $相关, 由各峰值幅度似然值决定.

    在声源数量未知的匹配场跟踪过程中, 模糊度函数各峰值位置$ {\boldsymbol{z}}_k^j $与各轨迹状态先验存在多种关联假设并且产生多个轨迹状态后验. 由于模糊度函数在非理想条件中不可避免地存在伪峰并且伪峰出现位置具有无序性, TPMBM原有关联步骤仅依据峰值位置距离似然求解最优关联, 易将伪峰与轨迹状态先验错误关联, 导致虚假轨迹出现; (36)式中加入峰值幅度信息, 即使最高峰值未对应声源位置, 以峰值位置距离似然和峰值幅度似然一致最大为准则, 仍能更加准确地反映各测量与各轨迹状态之间的关联程度, 从而抑制虚假轨迹以及取得更好的跟踪性能.

    作为总结, 表1列出了本文所提TPMBM-A-MFT方法在每个时间步的算法流程, 基于匹配场轨迹状态空间模型, 在轨迹空间依次进行预测、更新以及结合峰值位置和峰值幅度的数据关联步骤, 从而实现水下声源轨迹连续跟踪和数量估计.

    表 1  TPMBM-A-MFT算法流程
    Table 1.  Steps of TPMBM-A-MFT algorithm.
     步骤1 建立匹配场轨迹状态空间模型     变换模糊度函数: (26)式;      给出状态方程和观测方程: (27)式和(28)式;
     步骤2 根据(27)式对$ p({{\bf X}_k}) $参数集预测: 2.2.1节;
     步骤3 根据(28)式对$ p({{\bf X}_k}) $参数集更新: 2.2.2节;
     步骤4 结合峰值位置与幅度实现数据关联: (36)式;
     步骤5 估计轨迹状态和数量: (25)式.
    注: 初始时间步中$ {\bf{X} }_{k = 1}^{\rm{d}} = \varnothing $且泊松强度$ \lambda _{k = 1}^{\rm{u}} = \lambda _{k = 1}^{\rm{b}} $.
    下载: 导出CSV 
    | 显示表格

    本节利用SWellEx-96实验数据构建出双目标运动场景, 先展示现有MFT方法跟踪结果及局限性, 再给出本文所提TPMBM-A-MFT方法的跟踪结果, 并与3种RFS滤波器TPMBM, PMBM和GMPHD进行对比, 最后采用线性规划准则度量跟踪性能.

    SWellEX-96海试实验[22]是在加州圣地亚哥附近开展的浅海实验, 环境信息和实验数据集可在线查阅. 如图1(a)所示, 测线S5中发射船沿蓝色轨迹拖曳了9 m和54 m声源并以5节速度向垂直阵方向行进. 本节选取54 m声源发出的166 Hz信号进行数据分析, 首先将垂直阵(采样频率1500 Hz)记录时长75 min的声压数据分成156段(每段时长约28.4 s), 每段分成21个快拍, 重叠50%, 将第1—142段与第15—156段声压数据相加, 获得142段同频双目标声压数据, 时长缩短为68 min; 然后使用图1(b)所示波导环境参数和声速剖面, 由Kraken模型[23]计算搜索空间${\mathcal{Z}}{=}\left[ {0, 10{\kern 1 pt} {\kern 1 pt} {\kern 1 pt} {\text{km}}} \right] \times \left[ {0, 215{\kern 1 pt} {\kern 1 pt} {\kern 1 pt} {\kern 1 pt} {\text{m}}} \right]$内各拷贝向量(网格精度为0.05 km和2 m), 并与每段双目标声压数据进行匹配从而获得模糊度函数时间序列.

    图 1 SWellEx-96测线S5实验 (a)发射船轨迹; (b)波导环境示意图\r\nFig. 1. The SWellEx-96 event S5: (a) The launch ship track; (b) SWellEx-96 waveguide.
    图 1  SWellEx-96测线S5实验 (a)发射船轨迹; (b)波导环境示意图
    Fig. 1.  The SWellEx-96 event S5: (a) The launch ship track; (b) SWellEx-96 waveguide.

    为与本文所提MFT方法进行对比, 首先将Fialkowski等[3]所提MFT方法应用于上述模糊度函数时间序列, 该方法假设目标在一段时间保持直线运动且数量已知, 令搜索轨迹由初始与终点位置决定来缩减轨迹搜索维数, 以模糊度函数时间序列上沿声源轨迹累加幅值最大为原则给出跟踪结果, 如图2(a)中短折线所示, 圆圈与三角标记分别为每5 min内轨迹初始与终点位置估计结果. 图2(b)为双目标S1和S2距离维轨迹真值, 如点划线和点线所示, 目标S1与S2的最近距离点大约分别出现在第50 min和第60 min, 轨迹交叠发生在第55 min, 双目标运动过程中166 Hz信号被中断多次, 如图2(b)星号标记所示. 由于现有MFT方法受直线运动假设条件约束, 图2(a)所示轨迹跟踪结果需以合适时间区间分段给出. 对比图2(a)图2(b)可知, 图2(a)所示S1与S2轨迹跟踪结果均出现轨迹中断现象, 在第20—50 min内有多段虚假轨迹. 此外双目标运动过程还存在轨迹交叠, 伴随轨迹中断与虚假轨迹现象出现, 因此图2(a)中无法判断第55 min前后各条轨迹跟踪结果归属目标S1或S2, 易导致轨迹混叠.

    图 2 双目标距离维轨迹 (a) 现有匹配场跟踪方法轨迹估计结果; (b) 轨迹真值\r\nFig. 2. Range dimension trajectory of two targets: (a) MFT result; (b) the truth.
    图 2  双目标距离维轨迹 (a) 现有匹配场跟踪方法轨迹估计结果; (b) 轨迹真值
    Fig. 2.  Range dimension trajectory of two targets: (a) MFT result; (b) the truth.

    图3(a)(c)给出了影响匹配场跟踪方法性能的客观因素, 其中, 图3(a)为沿轨迹的水深变化, 可以看到, 前10 min出现较大地形差异, 易引起环境失配. 图3(b)图3(c)是平均阵元接收信噪比[24]和相关系数的变化曲线, 其中相关系数[25]为数据快拍与最佳拷贝向量的归一化内积, 失配程度越大则相关系数越小. 如图3(b)图3(c)所示, 在前20 min内, 受地形变化和目标距离相对较远的影响, 两目标的平均阵元接收信噪比和数据相关系数普遍较低; 在第20—45 min内, 目标S1信噪比和相关系数均低于目标S2; 另外在未发声时间段, 两目标的信噪比与相关系数极低. 低信噪比易引起模糊度函数中非目标位置出现伪峰, 低相关系数带来模糊度函数峰值偏移, 因此真实目标位置与峰值位置存在较大误差. 结合图2(b)图3(a)(c)可知, 本节构建的双目标匹配场跟踪数据集中不仅存在声源运动轨迹交叠和未持续发声情形, 还包含环境失配、低信噪比等客观因素影响, 这是图2(a)所示MFT方法出现轨迹中断、混叠和虚假轨迹现象, 导致不连续轨迹跟踪结果的原因.

    图 3 匹配场跟踪性能影响因素 (a)沿轨迹的水深变化; (b)平均阵元接收信噪比; (c)相关系数\r\nFig. 3. Influencing factors of matching field tracking performance: (a) The bottom depth along the trajectory; (b) average element level SNR; (c) correlation coefficient.
    图 3  匹配场跟踪性能影响因素 (a)沿轨迹的水深变化; (b)平均阵元接收信噪比; (c)相关系数
    Fig. 3.  Influencing factors of matching field tracking performance: (a) The bottom depth along the trajectory; (b) average element level SNR; (c) correlation coefficient.

    实际浅海被动跟踪场景中, 水下目标数量和运动状态通常难以确知. 现有MFT方法需将目标数量已知与直线运动作为先验条件, 其应用场景存在较大限制. 而本文所提方法假设目标数量未知, 对各模糊度函数提取${M_{1:K}}$个最高局部峰值, 峰值数量${M_{1:K}}$服从均值为5的泊松分布, 获得$ {{\bf{Z}}_{k = 1}}, {{\bf{Z}}_{k = 2}}, \cdots , {{\bf{Z}}_{k = 142}} $测量集合. 图4(a)图4(b)分别展示某次随机提取的各时间步测量集合的距离维空间状态分布及其元素数量, 从图4(a)可知, 前20 min内杂波测量的覆盖区域较大, 几乎无法分辨目标测量与杂波测量以及判断轨迹起始位置. 20 min后杂波测量的覆盖面积减少, 但与垂直阵相距较远的目标S1在第20—60 min内出现多次目标测量连续缺失, 受不发声时间段的影响较大. 此外距离8 km区域附近, 杂波测量密度较大.

    图 4 所有时间步测量集合 (a) 距离维空间分布; (b) 元素数量\r\nFig. 4. Measurements sets of all steps: (a) Spatial distribution of measurement elements in range dimension; (b) the number of elements in the measurement set.
    图 4  所有时间步测量集合 (a) 距离维空间分布; (b) 元素数量
    Fig. 4.  Measurements sets of all steps: (a) Spatial distribution of measurement elements in range dimension; (b) the number of elements in the measurement set.

    本节构建的双目标匹配场跟踪数据集中假设声源数量未知, 不仅存在声源轨迹交叉和未持续发声等情形, 而且环境失配等客观因素导致模糊度函数测量空间中的杂波密度和覆盖区域均较大. 综合这些因素以及现有MFT方法跟踪结果考虑, 由该双目标数据集实现轨迹连续跟踪具有较大挑战性.

    本节将TPMBM-A-MFT与TPMBM, PMBM和GMPHD三种RFS滤波器一同应用于4.1节建立的双目标跟踪数据集, 基于匹配场轨迹状态空间模型执行各自的滤波过程, 比较四种方法的匹配场跟踪结果. 其中PMBM和GMPHD只估计当前时间步的目标状态, 需以各时间步目标状态估计之间的最小距离为标签, 形成轨迹估计结果.

    表2列出了匹配场跟踪过程中的滤波参数设置, 其中, 检测概率${P^{\text{D}}}$依据图3(b)中信噪比曲线设置, 前20 min的$ {P^{\text{D}}} = 0.2 $, 第20—68 min内$ {P^{\text{D}}} = 0.5 $; 测量噪声标准差${\sigma _r}$${\sigma _d}$${\mathcal{Z}}$空间网格精度保持一致; 各时间步$k = 1, 2, \cdots , 142$出生的未检测轨迹状态先验均值设置为${\boldsymbol{\mu}} _k^b = {[z_{k, r}^ * , 0, z_{k, d}^ * , 0]^{\rm{T}}}$, 与模糊度函数最大峰值位置$ [z_{k, r}^ * , {\boldsymbol{z}}_{k, d}^ * ] = \mathop {\arg \max }\limits_{{z_k} \in {\mathcal{Z}}} B({{\boldsymbol{z}}_k}) $对应, 其他滤波参数如表2所列.

    表 2  匹配场跟踪过程的滤波参数设置
    Table 2.  Filter parameters setup for matched field tracking process.
    ${q_r}$${q_d}$${\sigma _r}$${\sigma _d}$$T$${w^{\text{b}}}$${P^{\text{S}}}$${P^{\text{D}}}$$\varGamma $
    $6 \times {10^{ - 3}}$$6 \times {10^{ - 3}}$$0.05$$2$$1$$5 \times {10^{ - 3}}$$0.99$$\{ 0.2, 0.5\} $$0.1$
    下载: 导出CSV 
    | 显示表格

    当4种方法应用于图4测量集合时, 图5(a)(d)依次给出了距离维轨迹跟踪结果. 图5(a)(d)中点划线和点线分别表示目标S1和S2的轨迹真值, 实线、虚线和加号标记分别表示S1, S2的估计轨迹和虚假轨迹.

    图 5 四种方法的距离维轨迹估计结果 (a) TPMBM-A-MFT; (b) TPMBM-MFT; (c) PMBM-MFT; (d) GMPHD-MFT\r\nFig. 5. Range dimension trajectory estimation results for four methods: (a) TPMBM-A-MFT; (b) TPMBM-MFT; (c) PMBM-MFT; (d) GMPHD-MFT.
    图 5  四种方法的距离维轨迹估计结果 (a) TPMBM-A-MFT; (b) TPMBM-MFT; (c) PMBM-MFT; (d) GMPHD-MFT
    Fig. 5.  Range dimension trajectory estimation results for four methods: (a) TPMBM-A-MFT; (b) TPMBM-MFT; (c) PMBM-MFT; (d) GMPHD-MFT.

    首先可以看出, 图5(a)图5(b)中目标S1和S2的估计轨迹完整, 而图5(c)图5(d)中S1和S2的估计轨迹均有中断现象出现, 原因在于漏检时间段中目标测量缺失, 导致目标状态后验仅依赖于目标状态先验, 随着缺失步长数增加, 估计轨迹的存在概率逐渐降低至小于阈值$\varGamma $时会出现轨迹中断现象. 对比TPMBM与PMBM, GMPHD滤波过程中的更新步骤, 每个时间步中PMBM和GMPHD只对当前时间步的目标状态进行更新, 而TPMBM会对估计轨迹中起始至当前时间步的目标状态序列进行更新. 当目标再次发声时, TPMBM和PMBM, GMPHD中估计轨迹的存在概率均上升至大于阈值, 但只有TPMBM中漏检时间段的后验目标状态序列会被后续时间步的目标测量重新修正, 轨迹中断现象消失.

    对比图5(d)图5(a)(c)可以发现, 图5(d)中不仅估计轨迹的破碎程度最大, 而且目标S1的估计轨迹分别在约第10 min, 20—35 min, 48 min和62 min共出现7次交叉混叠现象. 结合两目标的运动过程和RFS建模方式进一步分析原因, 从图2(b)知目标S1和S2在前40 min内保持约1 km间距平行运动且交替地出现未发声情形, 后28 min内存在轨迹交叠, 运动过程较为复杂. 对比TPMBM, PMBM与GMPHD的RFS建模方式, TPMBM和PMBM将轨迹(目标)状态RFS分成未检测和已检测两部分并分别由泊松和多伯努利混合参数集表示, GMPHD未考虑目标差异性并由泊松强度参数统一表示目标状态RFS, 其状态先验存在更大不确定性, 导致GMPHD的滤波过程更依赖于测量集合. 当S1和S2交替未发声时, GMPHD的目标状态后验估计易被测量集合迁移到另一条轨道上, 从而出现轨迹交叉混叠现象(图5(d)). 而在TPMBM和PMBM的滤波过程中, 鉴于多伯努利混合参数集通过多个参数准确描述了各已检测轨迹(目标)状态先验的确定程度, 已检测轨迹(目标)状态后验不轻易跟随测量位置跳变, 减少了估计轨迹交叉混叠的次数. 可以看到, 图5(c)中估计轨迹仅在第10 min出现一次交叉混叠, 图5(a)图5(b)中未有轨迹交叉混叠现象, 说明TPMBM将各轨迹状态建模为PMBM RFS并在轨迹空间滤波的方式, 保证了任意时间步之间目标S1和S2状态估计的连续性.

    图5(a)以外, 图5(b)(d)中的估计轨迹数量均大于真实值, 出现虚假轨迹. 即只有图5(a)中未出现轨迹中断、交叉混叠和虚假轨迹等现象, 实现了连续的轨迹跟踪结果. 图5(b)(d)中出现虚假轨迹的时间和空间集中在前10 min以及8—10 km区域, 从图4(a)可知前10 min内杂波测量覆盖区域较大且目标测量连续缺失, 可能导致轨迹状态先验与杂波关联, 使得轨迹状态后验逐渐偏离真值轨迹. 此外8 km区域附近的杂波测量密度较大同时存在时间和空间连续, 当未检测轨迹状态先验与杂波测量的距离似然值较小且被杂波测量更新时, 将驱动新的轨迹生成, 其存在概率随时间累积逐渐增大, 易形成虚假轨迹. 图5(a)(d)中仅图5(a)没有虚假轨迹出现, 对比图5(a)图5(b), 印证了TPMBM-A-MFT方法在数据关联步骤中通过增加峰值幅度似然信息来判决关联概率, 可以有效抑制杂波与轨迹状态先验关联以及虚假轨迹形成的可能性.

    为对TPMBM-A-MFT的跟踪过程进一步分析, 图6(a)(i)给出9个不同时间步的轨迹估计结果. 按照时间顺序进行说明, 如图6(a)所示, 在跟踪开始的前6 min内未出现轨迹估计, 原因是S1和S2距垂直阵较远且S1在第1—4 min未发声, 各轨迹状态后验的存在概率未超出阈值$\varGamma $, 轨迹均值未被提取; 从图6(b)知, 第6.77 min开始出现较近目标S2的轨迹估计; 如图6(c)图6(d)所示, 大约第11 min两目标轨迹估计由初始位置向正确位置移动, 由于前10 min杂波密度较大且目标测量连续缺失, 图6(c)出现一条虚假轨迹, 说明此时杂波与轨迹状态先验有错误关联, 但随着步长数增加以及目标测量出现, 图6(d)中轨迹状态先验与杂波的关联概率逐渐降低并不再是最优后验关联假设, 虚假轨迹已消失.

    图 6 TPMBM-A-MFT的跟踪过程 (a)第6.29 min; (b)第6.77 min; (c)第11.13 min; (d)第11.61 min; (e)第21.77 min; (f)第22.26 min; (g)第42.10 min; (h)第54.19 min; (i)第60.00 min\r\nFig. 6. Tracking process of TPMBM-A-MFT method: (a) 6.29 minutes; (b) 6.77 minutes; (c) 11.13 minutes; (d) 11.61 minutes; (e) 21.77 minutes; (f) 22.26 minutes; (g) 42.10 minutes; (h) 54.19 minutes; (i) 60.00 minutes.
    图 6  TPMBM-A-MFT的跟踪过程 (a)第6.29 min; (b)第6.77 min; (c)第11.13 min; (d)第11.61 min; (e)第21.77 min; (f)第22.26 min; (g)第42.10 min; (h)第54.19 min; (i)第60.00 min
    Fig. 6.  Tracking process of TPMBM-A-MFT method: (a) 6.29 minutes; (b) 6.77 minutes; (c) 11.13 minutes; (d) 11.61 minutes; (e) 21.77 minutes; (f) 22.26 minutes; (g) 42.10 minutes; (h) 54.19 minutes; (i) 60.00 minutes.

    由于目标S1和S2分别存在多个漏检时间段, 在图6(e)(i)所示后续时间步中估计轨迹出现了多次中断, 但在每段漏检发生的后续检测时间步中, 可观察到估计轨迹的中断部分又被重新连接. 以图6(e)图6(f)举例说明, 对比图6(e)中的点划线与实线所示S1真值和估计轨迹可知, S1估计轨迹的中断时间步约为第15—21.77 min, 当跟踪持续到第22.26 min时, 如图6(f)所示, S1估计轨迹与真值已基本重合, S1估计轨迹状态中第15至当前时间步的目标状态序列已被重新估计.

    以上现象归因于TPMBM-A-MFT方法在轨迹空间内进行滤波的方式, 其更新步骤利用当前时间步测量集合对整条轨迹状态修正, 轨迹状态包含了之前时间步的目标状态序列, 即使在目标漏检时间段内的估计轨迹中断, 一旦目标被重新检测到, 会将轨迹中断时间步的目标状态重新估计, 最终跟踪结果如图5(a)所示, 整个跟踪过程实现了连续轨迹估计. 可见, 相比传统MFT方法跟踪结果(图2(a)所示), 本文所提方法能在声源数量未知且间歇发声的运动场景中实现匹配场连续跟踪, 突破了传统MFT方法局限性.

    线性规划准则[26](linear programming metric, LP)通过匹配位置误差、轨迹缺失代价、虚假轨迹代价和轨迹混叠代价4个度量指标(文献[26]中(24)式有详细说明)对跟踪性能进行量化, 是衡量多目标跟踪算法综合性能的准则. 令匹配截止参数c和惩罚参数$\gamma $均为2、敏感参数p为1, 此时匹配位置误差是两条真值轨迹与各条估计轨迹进行位置匹配的最小绝对误差, 轨迹缺失、虚假轨迹和轨迹混叠代价分别为两条真值轨迹与其匹配的估计轨迹长度之差、两条真值轨迹与各估计轨迹数量之差和估计轨迹的匹配跳变次数, LP度量是4种度量指标之和, LP值越小则综合性能越优.

    4.1节中142段模糊度函数时间序列进行1000次独立随机的测量集合获取, 其中峰值数量${M_{1:K}}$服从泊松分布${\rm{ P}}\left( {{M_{1:K}};{\lambda ^z}} \right) $, 均值$ {\lambda ^z} $$[3,8]$区间任意选取. 表2跟踪参数不变, 在跟踪过程的每个时间步中对当前估计轨迹应用LP准则并按步长数进行归一化, 从而给出LP度量值以及四种度量指标随时间的变化情况用于跟踪性能评估, 分别如图7图8(a)(d)所示, 其中实线、虚线、点划线和点线分别为TPMBM-A-MFT, TPMBM-MFT, PMBM-MFT和GMPHD-MFT这4种跟踪方法的度量均值, 相应半透明填充区间表示4种跟踪方法度量标准差.

    图 7 四种跟踪方法在每个时间步的LP度量结果\r\nFig. 7. LP metrics at each time step for four methods.
    图 7  四种跟踪方法在每个时间步的LP度量结果
    Fig. 7.  LP metrics at each time step for four methods.
    图 8 LP的4种分解度量 (a)匹配位置误差; (b)轨迹缺失代价; (c)虚假轨迹代价; (d)轨迹混叠代价\r\nFig. 8. Four decomposition metrics of LP: (a) Matched localization error; (b) missed trajectory cost; (c) false trajectory cost; (d) switching cost.
    图 8  LP的4种分解度量 (a)匹配位置误差; (b)轨迹缺失代价; (c)虚假轨迹代价; (d)轨迹混叠代价
    Fig. 8.  Four decomposition metrics of LP: (a) Matched localization error; (b) missed trajectory cost; (c) false trajectory cost; (d) switching cost.

    图7可知, 跟踪前10 min内4种方法均保持较大LP均值, 在后续时间段4种跟踪方法的LP均值以不同速率减小, 说明跟踪过程在逐渐收敛. 最后时间步的LP值为最终轨迹估计的性能评估结果, 对比最后时间步的4个LP均值, TPMBM-A-MFT最小, TPMBM次之, 其他两种方法较大. 相比其他3种方法, 跟踪过程中TPMBM-A-MFT的LP标准差最低, 说明对测量集合均值${M_k}$的选取具有鲁棒性. 综合考虑, 本文所提TPMBM-A-MFT方法的跟踪性能最佳.

    图8(a)(d)展示了4种不同度量指标随时间的变化, 跟踪过程包含142个时间步(约68 min). 从最后时间步的度量结果来看, 其中, 图8(a)所示平均每个时间步TPMBM-A-MFT, TPMBM-MFT, PMBM-MFT和GMPHD-MFT对目标S1和S2位置估计的绝对误差之和分别约为0.14 km, 0.14 km, 0.28 km和0.49 km; 图8(b)所示TPMBM-A-MFT和TPMBM-MFT最终估计轨迹的缺失次数约为0次, 平均每个时间步PMBM-MFT和GMPHD-MFT的轨迹缺失次数约为0.5次和1.6次; 图8(c)中TPMBM-A-MFT的最终估计结果未包含虚假轨迹, TPMBM-MFT和PMBM-MFT约有$ 5.68 ( 0.04 \times 142) $条虚假轨迹, GMPHD-MFT最终包含的虚假轨迹条数最多, 是有$ 19.8\;(0.14 \times 142) $条; 图8(d)中TPMBM-A-MFT和TPMBM-MFT最终估计轨迹的混叠次数约为0次, PMBM-MFT和GMPHD-MFT最终估计轨迹混叠总次数约为$ 1.13\;(0.008 \times 142)$次和$ 5.68\;(0.04 \times 142) $次; 图8(a)(d)通过1000次蒙特卡罗实验表明本文所提方法具有低匹配位置误差和轨迹数量估计误差, 对两个存在轨迹重叠和未持续发声的目标能实现长时间稳健连续跟踪.

    对比图8(a)(d)中TPMBM-A-MFT和TPMBM-MFT的4种度量曲线, 其中两方法的匹配位置误差、轨迹缺失代价和轨迹混叠代价随时间的变化趋势大致相当, 而虚假轨迹代价除外. 两方法的唯一差别是, TPMBM-MFT只依赖于测量位置与轨迹先验位置距离似然值进行数据关联, 而TPMBM-A-MFT在数据关联步骤结合峰值位置和幅度测量信息, 将模糊度函数峰值幅度似然比作为区分目标测量与杂波的判别依据, 改变TPMBM原有数据关联方式. 对比图8(c)中两种方法的虚假轨迹代价曲线可知, TPMBM-A-MFT有效提高了关联决策准确性, 说明在数据关联中将距离似然与幅度似然结合的方式能简单有效约束虚假轨迹数量.

    现有匹配场跟踪方法在水下运动声源未持续发声、数量未知且轨迹交叠的跟踪场景中失效, 针对此问题, 本文提出一种基于轨迹泊松多伯努利混合滤波器的浅海匹配场连续跟踪方法(TPMBM-A-MFT). 该方法建立匹配场轨迹状态空间模型, 利用模糊度函数时间序列中峰值位置距离似然和峰值幅度似然的一致性, 在轨迹空间内进行预测、更新和数据关联步骤, 可以避免跟踪过程中出现轨迹中断、交叉混叠和虚假轨迹等现象, 最终实现轨迹连续跟踪及声源数量估计.

    SWellEx-96实验数据处理结果表明: 1)相比现有MFT方法需以声源数量已知为先验条件, 并在模糊度函数时间序列上沿搜索轨迹求和以实现声源跟踪, 所提方法将模糊度函数变换为delta函数之和并表示成测量RFS形式, 声源运动状态由轨迹状态RFS表示, 无需假设声源数量已知; 2)所提方法考虑了各声源轨迹状态的差异性并在轨迹空间内执行预测与更新步骤, 每个时间步测量RFS会对起始至当前时刻的目标状态先验序列矫正, 可减少声源间歇发声和轨迹交叉引起的轨迹中断与混叠现象; 3)所提方法更改了TPMBM原有数据关联步骤, 将模糊度函数峰值位置距离似然和峰值幅度似然联合表示关联概率, 可抑制虚假轨迹; 4)由线性规划准则中匹配位置误差、轨迹缺失、虚假轨迹和轨迹混叠代价4种度量指标, 验证了本文所提方法综合跟踪性能优势.

    [1]

    Bucker H 1976 J. Acoust. Soc. Am. 59 368Google Scholar

    [2]

    Bucker H 1994 J. Acoust. Soc. Am. 96 3809Google Scholar

    [3]

    Fialkowski L T, Perkins J S, Collins M D 2001 J. Acoust. Soc. Am. 110 739Google Scholar

    [4]

    Maranda B H, Fawcett J A 1991 IEEE J. Ocean. Eng. 16 189Google Scholar

    [5]

    Fawcett J A, Maranda B H 1993 J. Acoust. Soc. Am. 94 1363Google Scholar

    [6]

    Zala C A, Ozard J M, Wilmut M J 1998 J. Acoust. Soc. Am. 103 374Google Scholar

    [7]

    Tantum S L, Nolte L W 2002 J. Acoust. Soc. Am. 112 119Google Scholar

    [8]

    Mahler R P 2014 Advances in Statistical Multisource-Multitarget Information Fusion (Boston, London: Artech house) p83

    [9]

    Yardim C, Michalopoulou Z, Gerstoft P 2011 IEEE J. Ocean. Eng. 36 71Google Scholar

    [10]

    Vo B N, Ma W K 2006 IEEE Trans. Signal. Proces. 54 4091Google Scholar

    [11]

    Gruden P, White P R 2016 J. Acoust. Soc. Am. 140 1981Google Scholar

    [12]

    Gruden P, White P R 2016 J. Acoust. Soc. Am. 148 3014Google Scholar

    [13]

    Gruden P, Nosal E M 2021 J. Acoust. Soc. Am. 150 3399Google Scholar

    [14]

    Kupilik M J, Petersen T 2014 J. Acoust. Soc. Am. 136 1736Google Scholar

    [15]

    Georgescu R, Willett P 2012 IEEE J. Ocean. Eng. 37 220Google Scholar

    [16]

    García-Fernández Á F, Williams J L 2018 IEEE Trans. Signal. Proces. 54 1883Google Scholar

    [17]

    García-Fernández Á F, Svensson L 2020 IEEE Trans. Signal. Proces. 68 4933Google Scholar

    [18]

    Vo B T, Vo B N 2013 IEEE Trans. Signal. Proces. 61 3460Google Scholar

    [19]

    Tracey B H 2005 J. Acoust. Soc. Am. 118 1372Google Scholar

    [20]

    Lerro D, Bar-Shalom Y 1993 IEEE Trans. Aerosp. Electron. Syst. 29 404Google Scholar

    [21]

    Bendat J S, Piersol A G 2010 Random Data: Analysis and Measurement Procedures (4th Ed.) (Hoboken, NJ: Wiley) p604

    [22]

    Murray J, Ensberg D The SWellEx-96 experiment http://swellex96.ucsd.edu/ (Last viewed July 2021

    [23]

    Porter M B 1991 The KRAKEN Normal Mode Program (La Spezia: SACLANT Undersea Research Centre

    [24]

    Gemba K L, Nannuru S 2017 J. Acoust. Soc. Am. 141 3411Google Scholar

    [25]

    Booth N O, Baxley P A 1996 IEEE J. Ocean. Eng. 21 402Google Scholar

    [26]

    García-Fernández Á F, Rahmathullah A S 2020 IEEE Trans. Signal. Proces. 68 3917Google Scholar

  • 图 1  SWellEx-96测线S5实验 (a)发射船轨迹; (b)波导环境示意图

    Figure 1.  The SWellEx-96 event S5: (a) The launch ship track; (b) SWellEx-96 waveguide.

    图 2  双目标距离维轨迹 (a) 现有匹配场跟踪方法轨迹估计结果; (b) 轨迹真值

    Figure 2.  Range dimension trajectory of two targets: (a) MFT result; (b) the truth.

    图 3  匹配场跟踪性能影响因素 (a)沿轨迹的水深变化; (b)平均阵元接收信噪比; (c)相关系数

    Figure 3.  Influencing factors of matching field tracking performance: (a) The bottom depth along the trajectory; (b) average element level SNR; (c) correlation coefficient.

    图 4  所有时间步测量集合 (a) 距离维空间分布; (b) 元素数量

    Figure 4.  Measurements sets of all steps: (a) Spatial distribution of measurement elements in range dimension; (b) the number of elements in the measurement set.

    图 5  四种方法的距离维轨迹估计结果 (a) TPMBM-A-MFT; (b) TPMBM-MFT; (c) PMBM-MFT; (d) GMPHD-MFT

    Figure 5.  Range dimension trajectory estimation results for four methods: (a) TPMBM-A-MFT; (b) TPMBM-MFT; (c) PMBM-MFT; (d) GMPHD-MFT.

    图 6  TPMBM-A-MFT的跟踪过程 (a)第6.29 min; (b)第6.77 min; (c)第11.13 min; (d)第11.61 min; (e)第21.77 min; (f)第22.26 min; (g)第42.10 min; (h)第54.19 min; (i)第60.00 min

    Figure 6.  Tracking process of TPMBM-A-MFT method: (a) 6.29 minutes; (b) 6.77 minutes; (c) 11.13 minutes; (d) 11.61 minutes; (e) 21.77 minutes; (f) 22.26 minutes; (g) 42.10 minutes; (h) 54.19 minutes; (i) 60.00 minutes.

    图 7  四种跟踪方法在每个时间步的LP度量结果

    Figure 7.  LP metrics at each time step for four methods.

    图 8  LP的4种分解度量 (a)匹配位置误差; (b)轨迹缺失代价; (c)虚假轨迹代价; (d)轨迹混叠代价

    Figure 8.  Four decomposition metrics of LP: (a) Matched localization error; (b) missed trajectory cost; (c) false trajectory cost; (d) switching cost.

    表 1  TPMBM-A-MFT算法流程

    Table 1.  Steps of TPMBM-A-MFT algorithm.

     步骤1 建立匹配场轨迹状态空间模型     变换模糊度函数: (26)式;      给出状态方程和观测方程: (27)式和(28)式;
     步骤2 根据(27)式对$ p({{\bf X}_k}) $参数集预测: 2.2.1节;
     步骤3 根据(28)式对$ p({{\bf X}_k}) $参数集更新: 2.2.2节;
     步骤4 结合峰值位置与幅度实现数据关联: (36)式;
     步骤5 估计轨迹状态和数量: (25)式.
    注: 初始时间步中$ {\bf{X} }_{k = 1}^{\rm{d}} = \varnothing $且泊松强度$ \lambda _{k = 1}^{\rm{u}} = \lambda _{k = 1}^{\rm{b}} $.
    DownLoad: CSV

    表 2  匹配场跟踪过程的滤波参数设置

    Table 2.  Filter parameters setup for matched field tracking process.

    ${q_r}$${q_d}$${\sigma _r}$${\sigma _d}$$T$${w^{\text{b}}}$${P^{\text{S}}}$${P^{\text{D}}}$$\varGamma $
    $6 \times {10^{ - 3}}$$6 \times {10^{ - 3}}$$0.05$$2$$1$$5 \times {10^{ - 3}}$$0.99$$\{ 0.2, 0.5\} $$0.1$
    DownLoad: CSV
  • [1]

    Bucker H 1976 J. Acoust. Soc. Am. 59 368Google Scholar

    [2]

    Bucker H 1994 J. Acoust. Soc. Am. 96 3809Google Scholar

    [3]

    Fialkowski L T, Perkins J S, Collins M D 2001 J. Acoust. Soc. Am. 110 739Google Scholar

    [4]

    Maranda B H, Fawcett J A 1991 IEEE J. Ocean. Eng. 16 189Google Scholar

    [5]

    Fawcett J A, Maranda B H 1993 J. Acoust. Soc. Am. 94 1363Google Scholar

    [6]

    Zala C A, Ozard J M, Wilmut M J 1998 J. Acoust. Soc. Am. 103 374Google Scholar

    [7]

    Tantum S L, Nolte L W 2002 J. Acoust. Soc. Am. 112 119Google Scholar

    [8]

    Mahler R P 2014 Advances in Statistical Multisource-Multitarget Information Fusion (Boston, London: Artech house) p83

    [9]

    Yardim C, Michalopoulou Z, Gerstoft P 2011 IEEE J. Ocean. Eng. 36 71Google Scholar

    [10]

    Vo B N, Ma W K 2006 IEEE Trans. Signal. Proces. 54 4091Google Scholar

    [11]

    Gruden P, White P R 2016 J. Acoust. Soc. Am. 140 1981Google Scholar

    [12]

    Gruden P, White P R 2016 J. Acoust. Soc. Am. 148 3014Google Scholar

    [13]

    Gruden P, Nosal E M 2021 J. Acoust. Soc. Am. 150 3399Google Scholar

    [14]

    Kupilik M J, Petersen T 2014 J. Acoust. Soc. Am. 136 1736Google Scholar

    [15]

    Georgescu R, Willett P 2012 IEEE J. Ocean. Eng. 37 220Google Scholar

    [16]

    García-Fernández Á F, Williams J L 2018 IEEE Trans. Signal. Proces. 54 1883Google Scholar

    [17]

    García-Fernández Á F, Svensson L 2020 IEEE Trans. Signal. Proces. 68 4933Google Scholar

    [18]

    Vo B T, Vo B N 2013 IEEE Trans. Signal. Proces. 61 3460Google Scholar

    [19]

    Tracey B H 2005 J. Acoust. Soc. Am. 118 1372Google Scholar

    [20]

    Lerro D, Bar-Shalom Y 1993 IEEE Trans. Aerosp. Electron. Syst. 29 404Google Scholar

    [21]

    Bendat J S, Piersol A G 2010 Random Data: Analysis and Measurement Procedures (4th Ed.) (Hoboken, NJ: Wiley) p604

    [22]

    Murray J, Ensberg D The SWellEx-96 experiment http://swellex96.ucsd.edu/ (Last viewed July 2021

    [23]

    Porter M B 1991 The KRAKEN Normal Mode Program (La Spezia: SACLANT Undersea Research Centre

    [24]

    Gemba K L, Nannuru S 2017 J. Acoust. Soc. Am. 141 3411Google Scholar

    [25]

    Booth N O, Baxley P A 1996 IEEE J. Ocean. Eng. 21 402Google Scholar

    [26]

    García-Fernández Á F, Rahmathullah A S 2020 IEEE Trans. Signal. Proces. 68 3917Google Scholar

  • [1] Wang Lei, Huang Yi-Wang, Guo Lin, Ren Chao. Acoustic scattering modeling and sound field characteristics of rough seafloor in shallow sea. Acta Physica Sinica, 2024, 73(3): 034301. doi: 10.7498/aps.73.20231472
    [2] Li Yong-Fei, Guo Rui-Ming, Zhao Hang-Fang. Sparse reconstruction of acoustic interference fringes in shallow water and internal wave environment. Acta Physica Sinica, 2023, 72(7): 074301. doi: 10.7498/aps.72.20221932
    [3] Zhou Yu-Yuan, Sun Chao, Xie Lei, Liu Zong-Wei. A method of estimating depth of moving sound source in shallow sea based on incoherently matched beam-wavenumber. Acta Physica Sinica, 2023, 72(8): 084302. doi: 10.7498/aps.72.20222361
    [4] Tao Jian-Fei, Xia Qin-Zhi, Liao Lin-Gu, Liu Jie, Liu Xiao-Jing. Theory and application of photoelectron trajectory interference holography for atomic ionization in intense laser field. Acta Physica Sinica, 2022, 71(23): 233206. doi: 10.7498/aps.71.20221296
    [5] Wang Xuan, Sun Chao, Li Ming-Yang, Zhang Shao-Dong. Detection by angle-domain subspace with horizontal array in uncertain shallow-water environment. Acta Physica Sinica, 2022, 71(8): 084304. doi: 10.7498/aps.71.20211742
    [6] Zhang Shao-Dong, Sun Chao, Xie Lei, Liu Xiong-Hou, Wang Xuan. Influence of environmental uncertainty on source power estimation in shallow water waveguide. Acta Physica Sinica, 2021, 70(24): 244301. doi: 10.7498/aps.70.20210852
    [7] Zhou Yao-Zhi, Li Chun, Li Chen-Yang, Li Qing-Lian. Prediction of liquid jet trajectory in supersonic crossflow and continuous liquid column model. Acta Physica Sinica, 2020, 69(23): 234702. doi: 10.7498/aps.69.20200903
    [8] Kong De-Zhi, Sun Chao, Li Ming-Yang, Zhuo Jie, Liu Xiong-Hou. Dimension-reduced generalized likelihood ratio detection based on sampling of normal modes in deep ocean. Acta Physica Sinica, 2019, 68(17): 174301. doi: 10.7498/aps.68.20190700
    [9] Qian Zhi-Wen, Shang De-Jiang, Sun Qi-Hang, He Yuan-An, Zhai Jing-Sheng. Acoustic radiation from a cylinder in shallow water by finite element-parabolic equation method. Acta Physica Sinica, 2019, 68(2): 024301. doi: 10.7498/aps.68.20181452
    [10] Meng Rui-Jie, Zhou Shi-Hong, Li Feng-Hua, Qi Yu-Bo. Identification of interference normal mode pairs of low frequency sound in shallow water. Acta Physica Sinica, 2019, 68(13): 134304. doi: 10.7498/aps.68.20190221
    [11] Li Xiao-Man, Zhang Ming-Hui, Zhang Hai-Gang, Piao Sheng-Chun, Liu Ya-Qin, Zhou Jian-Bo. A passive range method of broadband impulse source based on matched-mode processing. Acta Physica Sinica, 2017, 66(9): 094302. doi: 10.7498/aps.66.094302
    [12] Lin Cheng, Zhang Hua-Tang, Sheng Zhi-Hao, Yu Xian-Huan, Liu Peng, Xu Jing-Wen, Song Xiao-Hong, Hu Shi-Lin, Chen Jing, Yang Wei-Feng. Strong field photoelectron holography studied by a generalized quantum-trajectory Monte Carlo method. Acta Physica Sinica, 2016, 65(22): 223207. doi: 10.7498/aps.65.223207
    [13] Qi Yu-Bo, Zhou Shi-Hong, Zhang Ren-He. Warping transform of the refractive normal mode in a shallow water waveguide. Acta Physica Sinica, 2016, 65(13): 134301. doi: 10.7498/aps.65.134301
    [14] Zhang Tong-Wei, Yang Kun-De. A virtual time reversal method for passive source localization in a range-dependent waveguide. Acta Physica Sinica, 2014, 63(21): 214303. doi: 10.7498/aps.63.214303
    [15] Zhang Tong-Wei, Yang Kun-De, Ma Yuan-Liang, Li Xue-Gang. The performance of matched-field localization with a horizontal line array at different depths in shallow water. Acta Physica Sinica, 2010, 59(5): 3294-3301. doi: 10.7498/aps.59.3294
    [16] Yu Yun, Hui Jun-Ying, Chen Yang, Sun Guo-Cang, Teng Chao. Research on target depth classification in low frequency acoustic field of shallow water. Acta Physica Sinica, 2009, 58(9): 6335-6343. doi: 10.7498/aps.58.6335
    [17] Zheng Chun-Lan, Li Tong-Bao, Ma Yan, Ma Shan-Shan, Zhang Bao-Wu. Analysis of Cr atom trajectory and focusing deposition in the standing wave field. Acta Physica Sinica, 2006, 55(9): 4528-4534. doi: 10.7498/aps.55.4528
    [18] ZHANG JIA-SHU, XIAO XIAN-CI. NONLINEAR ADAPTIVE PREDICTIVE TARGETING CONTROL OF THE CONTINUOUS CHAOTIC SYSTEMS. Acta Physica Sinica, 2001, 50(11): 2092-2096. doi: 10.7498/aps.50.2092
    [19] CHENG REN-SHU, XIMEN JI-YE. ION OPTICAL THEORY FOR THE THIRD ORDER TRAJECTORY IN A CROSSED TOROIDAL ELECTRIC AND INHOMOGENEOUS MAGNETIC FIELD (I)——CALCULATION OF TRAJECTORY AND MATRIX FORMALISM. Acta Physica Sinica, 1982, 31(6): 722-737. doi: 10.7498/aps.31.722
    [20] . Acta Physica Sinica, 1975, 24(3): 200-209. doi: 10.7498/aps.24.200
Metrics
  • Abstract views:  3813
  • PDF Downloads:  63
Publishing process
  • Received Date:  31 January 2023
  • Accepted Date:  06 June 2023
  • Available Online:  26 July 2023
  • Published Online:  20 September 2023

/

返回文章
返回