Multi-maneuvering Acoustic Targets Tracking Algorithm Based on Virtual Extension of Single Acoustic Vector Sensor
-
摘要: 为解决单声矢量传声器(Acoustic vector sensor, AVS)可跟踪声目标数目少、跟踪性能差的问题, 提出了基于AVS虚拟扩展的多机动声目标跟踪算法. 首先, 引入高阶累积量预处理过程并建立高阶似然函数, 不仅能够抑制高斯噪声、提高估计精度, 还可通过AVS的虚拟扩展增加可跟踪目标数目. 然后, 在边缘化
$\delta$ 广义标签多伯努利(Marginalized$\delta$ -generalized label multi-bernoulli, M$\delta$ -GLMB)滤波框架下, 提出了基于累积量的增广运动模型状态的M$\delta$ -GLMB (Cumulants-based augumented motion model state M$\delta$ -GLMB, Cum-AMMS-GLMB)算法. 算法引入多种运动模型, 并将表征不同模型的索引标号作为目标状态的增广参数, 通过各模型间的加权混合获取优于单一运动模型的跟踪性能. 除此之外, 算法的序贯蒙特卡洛(Sequential Monte Carlo, SMC)实现过程中, 依据高阶预处理获得的归一化空间谱拟合检测概率函数, 抑制了杂波向可用粒子扩展, 进一步增强了高似然区域的粒子. 最后, 推导了AVS目标跟踪的后验克拉美罗下界(Posterior cram$\acute{e}$ r-rao lower bound, PCRLB), 并通过仿真实验验证了算法的量测噪声抑制能力和声目标跟踪性能.-
关键词:
- 声矢量传声器 /
- 高阶累积量 /
- 虚拟扩展 /
- 广义标签多伯努利滤波 /
- 多目标跟踪
Abstract: To solve the problem of poor performance and fewer trackable targets in the multi-targets tracking of acoustic vector sensor (AVS) in multi-target tracking, a multi-maneuvering acoustic targets tracking algorithm based on virtual extension of single AVS is proposed. First, the higher-order cumulants processing method is introduced to establish a higher-order likelihood function, which can not only improve the estimation accuracy by suppressing the Gaussian noise, but also increase the number of estimable targets by virtually extending the AVS. Then, under the marginalized$\delta$ -generalized label multi-bernoulli (M$\delta$ -GLMB) framework, a cumulants-based augmented motion model state M$\delta$ -GLMB (Cum-AMMS-GLMB) algorithm is proposed. The algorithm introduces multiple models, and uses the model index labels that distinguish different motion models as an augmented parameter for the target state, and obtains a better tracking performance than a single motion model through weighted mixing of the updated states of each model. In addition, in the sequential Monte Carlo (SMC) implementation of the algorithm, the detection probability function is fitted based on the normalized spatial spectrum obtained by higher-order cumulants preprocessing can suppress the diffusion of clutter to the available particles, and further enhance the particles in the high-likelihood region. Finally, the posterior cram$\acute{e}$ r-rao lower bound (PCRLB) for targets tracking of single AVS is derived, and the performance of measurement noise suppression and acoustic targets tracking is verified by simulation experiments. -
声矢量传声器(Acoustic vector sensor, AVS)可同步测量空间一点处声压和笛卡尔坐标系三个正交方向上的声质点振速信息[1-4]. 单个AVS利用其对频率不敏感的“8”字形指向性, 即可实现需由多个声压传声器组成阵列才可实现的声源波达方向(Direction-of-arrival, DOA) 估计. 鉴于AVS的小尺寸便携、声场信息获取等方面的优势, 因此为机器人听觉感知系统的发展提供新思路.
近年来, 最大似然[5]、波束形成[6]、基于子空间分解的MUSIC[7-8]、ESPRIT[9-10]等算法已广泛用于基于AVS的DOA估计问题. 尽管小尺寸AVS可实现DOA估计, 但其孔径一定程度上限制了角分辨率. 高阶累积量处理手段的引入[11], 可通过AVS虚拟扩展和量测噪声抑制达到提高DOA估计精度的目的. 实际场景中多为运动声目标, 且伴随目标出生、死亡, 目标数目是时变的. 然而, 上述算法主要针对固定数目的静止声目标DOA估计, 且未考虑相邻时间步目标状态的高度相关性.
基于随机有限集(Random finite set, RFS)理论的目标跟踪方法[12] 可解决时变多目标跟踪问题. 该类方法将多目标状态和量测建模为RFS, 并通过多目标贝叶斯滤波传递多目标后验概率, 从而达到目标数目与状态参数估计的目的. 目前常用的有概率假设密度(Probability hypothesis density, PHD) [13-15]、势概率假设密度(Cardinalized PHD, CPHD)[16]、多目标多伯努利(Multi-target multi-Bernoulli, MeMBer)[17-20]以及广义标签多伯努利(Generalized label multi-Bernoulli, GLMB)[21-26] 等滤波器. 这些算法能够避免如多目标假设[27]、联合概率数据关联[28]等算法的数据关联, 极大程度上降低计算复杂度, 提高目标跟踪精度与实时性. 在基于AVS的目标跟踪实现中, 文献[29-32]采用RFS来描述状态过程的随机性, 并通过粒子滤波实现时变二维DOA估计. 文献[33]将MeMBer滤波器引入单个AVS的目标跟踪中, 实现了两个不同频率窄带目标检测与跟踪. 文献[34]则进一步引入交互式多模型(Interacting multiple model, IMM)和指数加权的MUSIC伪谱似然函数以提高跟踪精度. 然而, 上述基于单个AVS的目标跟踪算法最多实现双源跟踪, 且无法分辨各目标的航迹信息.
基于标签RFS理论的GLMB、
$ \delta $ -GLMB、M$ \delta $ -GLMB[35, 36]等滤波器通过引入航迹标签信息, 可解决其他非标签RFS滤波器无法区分目标航迹的问题. M$ \delta $ -GLMB滤波器作为一种近似的$ \delta $ -GLMB方法, 通过对量测−航迹的关联映射执行边缘化处理, 能够减少航迹假设数目、降低计算成本. 鉴于此, 本文在M$ \delta $ -GLMB滤波框架下设计AVS虚拟扩展的多机动声目标跟踪算法, 即Cum-AMMS-GLMB算法, 以解决AVS可跟踪声目标数目少、低信噪比下跟踪性能差的问题. 该算法主要贡献为:$ 1) $ 高阶累积量预处理的提出, 能够抑制高斯量测噪声, 提高目标跟踪精度. 除此之外, 高阶似然的设计可实现AVS虚拟扩展. 与现有基于空间谱伪似然的AVS跟踪算法相比, 可增加可跟踪目标数目.$ 2) $ 实际场景中, 真实声目标的运动轨迹无法由单一运动模型描述. Cum-AMMS-GLMB算法在M$ \delta $ -GLMB框架下嵌入IMM思想, 将表征不同运动模型的索引作为扩展状态, 通过各模型间的加权更新提高跟踪性能. 与现有IMM-GLMB算法[37-39]不同, 该算法量测为声目标混合信号, 无需额外引入航迹−量测的关联映射, 能够降低滤波主体的存储容量和计算成本, 抵消预处理过程的部分计算负担.$ 3) $ 实际情况下新生目标先验未知、且检测概率可时变. 算法实现过程中设计了新生目标提取、以及检测概率拟合方法. 其中检测概率函数通过预处理获得的归一化高阶空间谱拟合, 可达到抑制杂波向可用粒子扩散、增强高似然区域粒子的目的.除以上主要贡献外, 本文还推导了多模型交互的AVS目标跟踪的PCRLB, 并通过实验验证了Cum-AMMS-GLMB算法的可行性和有效性.
1. 声目标系统建模
1.1 目标运动模型
在声目标跟踪问题中,
$ k $ 时刻一个单声目标状态为$ {\boldsymbol{\theta}}_k=\{\alpha_k,{\dot\alpha}_k,{\ddot \alpha}_k,\beta_k,{\dot\beta}_k,{\ddot \beta}_k\} $ , 其中$ \alpha_k $ 和$ \beta_k $ 分别表示声目标在$ k $ 时刻的方位角和俯仰角,$ \dot{(\cdot)} $ 和$ \ddot{(\cdot)} $ 分别表示速度和加速度. 考虑到声目标真实运动模型未知, 因此可引入多个运动模型, 并通过一定的概率逼近真实声目标的运动状态. 则$ k -1$ 时刻的单声目标状态向$ k $ 时刻转换的运动方程为:$$ \begin{equation} {\boldsymbol{\theta}}_k={{\boldsymbol{ F}}}^{(o)}{\boldsymbol{\theta}}_{k-1}+{{\boldsymbol{ G}}}^{(o)}{{\boldsymbol{ w}}}_k,o=1,2,3,\cdots \end{equation} $$ (1) 其中
$ {\boldsymbol{\theta}}_{k-1} $ 为声目标在$ k-1 $ 时刻的状态.$ {\boldsymbol{F}}^{\left(o\right)} $ 和$ {\boldsymbol{G}}^{\left(o\right)} $ 分别对应第$ o $ 个运动模型的状态转移矩阵和噪声驱动矩阵,$ o\in {O} $ 为运动模型索引,$O=\{1,2,3,\cdots\}$ 为模型索引空间. 过程噪声$ {\boldsymbol{w}}_k=[w_k^{(\alpha)},w_k^{(\beta)}]^{\mathrm T} $ 服从零均值, 协方差为$ {\boldsymbol{Q}}_k=\sigma_w^2{{\boldsymbol{ I}}}_2 $ 的高斯白噪声, 其中$ \sigma_w $ 表示加速度扰动标准差. I2表示2 × 2维的单位矩阵实际应用场景中, 声目标可能突然出现、消失, 目标状态和数目均随机变化. 每个时刻的声目标状态和数目的变化可由标签RFS描述, 即:
$$ \begin{split} {{\boldsymbol{\Theta }}_k} =\;& \{ ({{\boldsymbol{\theta }}_{k,1}},{l_{k,1}}),\cdots,({{\boldsymbol{\theta }}_{k,M(k)}},{l_{k,M(k)}})\} \in\\ &F({\boldsymbol{\Theta }} \times L) \end{split} $$ 其中
$ M(k) $ 为$ k $ 时刻的声目标数目.$ {L} $ 为离散标签空间.$ l_{k,m} $ 表示$ k $ 时刻第$ m $ 个声目标对应的标签, 每个标签存储二维变量$ (t,i) $ , 其中$ t $ 记录目标出生时刻,$ i $ 用于区别同时出生的目标. 上述标签RFS对每个声目标状态$ {\boldsymbol{\theta}}\in {\boldsymbol{\Theta}} $ 赋予独一无二的标签$ l\in {L} $ , 能够将目标状态扩展为$ \{({\boldsymbol{\theta}},l)\in {{\boldsymbol{ \Theta}}}\times {L}\} $ , 从而解决非标签RFS无法区分声目标航迹的问题.1.2 AVS量测模型
单个AVS可同时测量声场中声压和三个正交方向上声质点振速, 并且不敏感于声目标的频域特征, 直接对时域数据进行处理即可实现声目标DOA估计. 假设
$ k $ 时刻$ M(k) $ 个声目标以$\{\alpha_{k,m},\beta_{k,m}\},m= 1,\cdots,M(k)$ 的来向入射到AVS, 该时刻传声器共接收到$ L $ 个采样快拍, 则$ 4\times L $ 维接收数据模型为:$$ \begin{equation} {\boldsymbol{Z}}_k=\sum\limits_{m=1}^{M(k)}{\boldsymbol{a}(\alpha_{k,m},\beta_{k,m}){\boldsymbol{s}}_{k,m}}+{\boldsymbol{V}}_k ={\boldsymbol{A}}_k{\boldsymbol{s}}_k+{\boldsymbol{V}}_k \end{equation} $$ (2) 其中
$ {\boldsymbol{A}}_k=[{\boldsymbol{a}}(\alpha_{k,1},\beta_{k,1}),\cdots,{\boldsymbol{a}}(\alpha_{k,M(k)},\beta_{k,M(k)})] $ 表示$ 4\times M\left(k\right) $ 维阵列流型,${\boldsymbol{a}}({\alpha _{k,\;m}},\;{\beta _{k,\;m}}) \;= \left[ {\cos {\alpha _{k,m}}\cos {\beta _{k,m}}},\,{\sin {\alpha _{k,m}}\cos {\beta _{k,m}}},\,{\sin {\beta _{k,m}}},1\right]^{\mathrm{T} }$ 表示第$ m $ 个声目标关于AVS的方向向量.${\boldsymbol{s}}_k = [{\boldsymbol{s}}_{k,1}^{\mathrm T},\cdots, {\boldsymbol{s}}_{k,M(k)}^{\mathrm T}]^{\mathrm{T} }$ ,$ 1\times L $ 维$ {\boldsymbol{s}}_{k,m} $ 表示第$ m $ 个目标在$ k $ 时刻的包络.$ 4\times L $ 维$ {\boldsymbol{V}}_k $ 表示均值为0, 协方差为${{\boldsymbol{ \varGamma}}}_0= {\mathrm {diag}}\{\sigma_v^2,\sigma_v^2,\sigma_v^2,\sigma_p^2\}$ 的量测高斯噪声, 其中$ \sigma_v^2 $ 和$ \sigma_p^2 $ 分别表示AVS中声质点振速通道和声压通道的噪声功率.依据式(2), AVS在
$ k $ 时刻的接收数据可表述为如下非线性量测方程:$$ \begin{equation} {\boldsymbol{Z}}_k={H}_k({\boldsymbol{\Theta}}_k)+{\boldsymbol{V}}_k \end{equation} $$ (3) 其中
$ {H}_k(\cdot) $ 为非线性量测函数. 量测模型是不同状态声目标的混合信号, 仅存在一组$ 4\times L $ 维的量测值, 因此量测集可表示为$ {\boldsymbol{Z}}_k=\{\boldsymbol{Z}_k\} $ . 对于给定声目标状态$ {\boldsymbol{\theta}}_k \in {{\boldsymbol{ \Theta}}}_k $ , 其被检测的概率为$ p_{D,k}({\boldsymbol{\theta}}_k) $ , 漏检概率为$ q_{D,k}({\boldsymbol{\theta}}_k)=1-p_{D,k}({\boldsymbol{\theta}}_k) $ . 量测噪声$ {\boldsymbol{V}}_k $ 为高斯噪声, 因此$ k $ 时刻量测$ {\boldsymbol{Z}}_k $ 服从复高斯分布$ {\boldsymbol{Z}}_k \sim {{CN}}({H}_k({{\boldsymbol{ \Theta}}}_k),{\boldsymbol{R}}_{Z,k}) $ , 即似然函数可表示为:$$ \begin{equation} \nonumber \begin{aligned} &g({{\boldsymbol{Z}}_k}|{{\boldsymbol{\theta }}_k})= \\ &\quad \frac{{\exp \left( { - {{\left( {{{\boldsymbol{Z}}_k} - {{{H}}_k}({{\boldsymbol{\Theta }}_k})} \right)}^{\mathrm H}}{\boldsymbol{R}}_{Z,k}^{ - 1}\left( {{{\boldsymbol{Z}}_k} - {{{H}}_k}({{\boldsymbol{\Theta }}_k})} \right)} \right)}}{{{\pi ^4}\det ({{\boldsymbol{R}}_{Z,k}})}} \end{aligned} \end{equation} $$ 其中
$ {\boldsymbol{R}}_{Z,k}={\mathrm E}\{\boldsymbol{Z}_k{\boldsymbol{Z}}_k^{\mathrm H}\}={{\boldsymbol{A}}_k}{{\boldsymbol{P}}_k}{\boldsymbol{A}}_k^{\mathrm H} + {{\boldsymbol{\varGamma }}_0} $ 为量测的协方差矩阵, 矩阵$ {\boldsymbol{P}}_k $ 用于存放声目标信号功率. 协方差矩阵可由${{\boldsymbol{\hat R}}_{Z,k}} \approx ({1}/{L}){{\boldsymbol{Z}}_k}{\boldsymbol{Z}}_k^{\mathrm H}$ 近似获取. 声目标信号通常未知, 因此无法获取$ {{{H}}_k}({{\boldsymbol{\Theta }}_k}) $ , 不便于求解上述似然函数. 然而, 空间谱易于获取, 当目标状态接近真实目标状态时才会达到最大值, 与粒子滤波似然函数具有相似效果[26, 30-34]. 因此后续算法研究中, 可依据空间谱这一特点设计新的似然函数.2. AVS虚拟扩展的多机动声目标跟踪算法
本部分提出AVS虚拟扩展的多机动声目标跟踪算法, 即Cum-AMMS-GLMB算法. 算法首先引入预处理方法, 设计新的高阶似然函数模型, 以抑制量测噪声、增加可跟踪目标数目. 然后在M
$ \delta $ -GLMB框架中引入IMM概念[37-40], 并将运动模型索引作为隐含状态扩展到声目标状态中, 进一步推导给出Cum-AMMS-GLMB算法的预测、更新步骤.2.1 AVS虚拟扩展预处理
声目标DOA估计通常由空间谱估计方法实现. 文献[41]指出阵列可估计声源数目为:
$$ \begin{equation} M^\prime\le\frac{N\cdot \mathrm{rank}({\boldsymbol{R}}_Z)}{p+q\cdot \mathrm{rank}({\boldsymbol{R}}_Z)} \end{equation} $$ (4) 其中
$\mathrm{rank}(\cdot)$ 表示矩阵求秩,$ {\boldsymbol{R}}_Z $ 为接收数据的协方差矩阵,$ p $ 为信号待估参量维度,$ q $ 为声信号维度,$ N $ 表示阵列通道数目. 对四通道AVS的二维DOA估计而言, 其可估计声目标数目满足$ M^\prime\le2.67 $ , 即最大可分辨声目标数目为$ 2 $ . 实际应用场景中, 通常期望使用更少的传感器跟踪更多的声目标. 本部分引入高阶累积量预处理方法, 并设计新的似然函数模型, 通过对AVS的虚拟扩展增加可跟踪声目标数目.2.1.1 高阶累积量预处理
在执行高阶累积量预处理前, 需对式(2)所示AVS接收数据模型做如下假设:
$ 1) $ 声目标信号$ {\boldsymbol{s}}_{k,m},m=1,\cdots,M(k) $ 均为零均值、平稳、非高斯随机过程, 且在统计上彼此独立;$ 2) $ 量测噪声$ {\boldsymbol{V}}_k $ 为零均值高斯噪声, 该高斯噪声与入射声目标信号统计上独立.根据以上假设, 预处理方法适用于高斯噪声背景下非高斯信号处理. 现实世界中, 声源通常为非高斯信号[11, 42-43]. 因此, 上述假设在实际应用中同样适用.
根据高阶累积量的独特性质[43-44], 对式(2)进行四阶累积量预处理:
$$ \begin{split} &{C_{4Z,k}}( {{i_1},{i_2},i_3 ,i_4 } ) =\frac{1}{L}\sum\limits_{l' = 1}^L cum\{ {{\boldsymbol{Z}}_k}( {{i_1},l'} ),\\ &{{\boldsymbol{Z}}_k^*}( {{i_2},l'} ),{{\boldsymbol{Z}}_k}( {{i_3},l'}),{{\boldsymbol{Z}}_k^*}( {{i_4},l'} )\} =\\ &\frac{1}{L}\sum\limits_{l' = 1}^L \Bigg\{ cum \Bigg\{ \sum\limits_{m = 1}^{M( k )} {{{\boldsymbol{A}}_k}( {{i_1},m} ){\boldsymbol{s}_{k,m}}( l' )} ,\cdots,\\ &\sum\limits_{n = 1}^{M( k )} {{{\boldsymbol{A}}_k^*}( {{i_4},n}){\boldsymbol{s}_{k,n}^*}( l' )} \Bigg\} + \\ &cum \left\{ {{\boldsymbol{V}_k(i_1,l')}, \cdots ,{\boldsymbol{V}_k^*(i_4,l')}} \right\} \Bigg\}= \\ &\frac{1}{L}\sum\limits_{l' = 1}^L cum \Bigg\{ \sum\limits_{m = 1}^{M( k )} {{{\boldsymbol{A}}_k}( {{i_1},m} ){\boldsymbol{s}_{k,m}}( l' )} ,\cdots,\\ &\sum\limits_{n = 1}^{M( k )} {{{\boldsymbol{A}}_k^*}( {{i_4},n}){\boldsymbol{s}_{k,n}^*}( l' )} \Bigg\}= \\ &\sum\limits_{m = 1}^{M( k)} {{{\boldsymbol{A}}_k}( {{i_1},m} ){{\boldsymbol{A}}_k^*}( {{i_2},m} ){{\boldsymbol{A}}_k}( {{i_3},m}){{\boldsymbol{A}}_k^*}( {{i_4},m}){\gamma _{4,{s_{k,m}}}}} \end{split} $$ (5) 其中
$ 1\le i_1,i_2,i_3,i_4\le4 $ 为AVS通道取值.$ {cum}\{\cdot\} $ 表示四阶累积量计算.$ (\cdot)^* $ 表示复数共轭.$ {\boldsymbol{s}}_{k,m}(l') $ 表示$ 1\times L $ 维向量$ {\boldsymbol{s}}_{k,m} $ 的第$ l' $ 个元素.$ {\boldsymbol{V}}_k(i,l') $ 表示$ 4\times L $ 维量测噪声$ {\boldsymbol{V}}_k $ 的第$ (i,l') $ 个元素.$ {\boldsymbol{A}}_k(i,m) $ 表示$ k $ 时刻阵列流型$ {\boldsymbol{A}}_k $ 的第$ (i,m) $ 个元素.$ {\gamma_{4,{s_{k,m}}}} $ 为$ k $ 时刻第$ m $ 个声目标的四阶累积量. 已知量测噪声为高斯白噪声, 预处理过程首先依据“高斯随机变量的高阶累积量恒等于0”的性质可得$({1}/ {L})\sum_{l' = 1}^L {{{cum}}\left\{ {{\boldsymbol{V}_k(i_1,l')}, \cdots ,{\boldsymbol{V}_k^*(i_4,l')}} \right\}} = 0$ , 从而达到抑制量测噪声的目的. 然后, 依据性质“多个随机变量的一个子集同其他部分独立, 则高阶累积量为0”可得, 当且仅当$ m=i=j=n $ 时,$({1}/{L})\sum_{l' = 1}^L {{cum}} \{ {{\boldsymbol{s}_{k,m}}(l'),{\boldsymbol{s}}_{k,i}^*(l'),{\boldsymbol{s}_{k,j}}(l'),{\boldsymbol{s}}_{k,n}^*(l')} \}$ 等于$ {\gamma _{4,{s_{k,m}}}} $ , 否则为0. 最终, 式(5)预处理结果中抑制了高斯量测噪声, 仅剩余目标信号的四阶累积量以及包含DOA信息的阵列流型的叠加信号.将式(5)预处理结果放入
$ 4^2\times4^2 $ 维累积量矩阵:$$ \begin{equation} { \begin{array}{l} {\boldsymbol{R}}_{Z,k}^{(\rm{cum})} = {{\boldsymbol{B}}_k}{{\boldsymbol{D}}_k}{\boldsymbol{B}}_k^{\mathrm H}=\\ \quad\quad\ \sum\limits_{m = 1}^{M(k)} {{\boldsymbol{b}}({\alpha _{k,m}},{\beta _{k,m}})} {\gamma _{4,{s_{k,m}}}}{\boldsymbol{b}}{({\alpha _{k,m}},{\beta _{k,m}})^{\mathrm H}} \end{array} } \end{equation} $$ (6) 其中
$ {\boldsymbol{B}}_k=[{\boldsymbol{b}}(\alpha_{k,1},\beta_{k,1}),\cdots,{\boldsymbol{b}}(\alpha_{k,M(k)},\beta_{k,M(k)})] $ ,${\boldsymbol{b}}(\alpha_{k,m},\beta_{k,m})\;=\;{\boldsymbol{a}}(\alpha_{k,m},\beta_{k,m})\otimes{\boldsymbol{a}}^\ast(\alpha_{k,m},\beta_{k,m})$ ,$ {\boldsymbol{D}}_k={\mathrm {diag}}\{{\gamma _{4,{s_{k,1}}}},\cdots,{\gamma _{4,{s_{k,M(k)}}}}\} $ . 四阶累积量矩阵$ {\boldsymbol{R}}_{Z,k}^{(\rm{cum})} $ 与协方差矩阵$ {{\boldsymbol{R}}_{Z,k}} $ 的区别在于: 1) 抑制了高斯量测噪声; 2) 预处理后的四通道AVS虚拟扩展为一个导向矢量形如$ {\boldsymbol{b}}(\alpha_{k,m},\beta_{k,m}) $ 的16通道冗余虚拟阵列. 去除冗余通道后, 虚拟阵列为一个10通道传声器. 根据式(4)可得, AVS虚拟扩展后的可估计声目标数目满足$ M^\prime\le8.33 $ , 即预处理后最大可估计8个目标.总结来说, 高阶累积量预处理的引入不仅能够抑制高斯量测噪声, 还可通过AVS的虚拟扩展以增加可跟踪声目标数目.
2.1.2 量测似然函数
基于高阶累积量预处理的噪声抑制和阵列虚拟扩展优势, 本部分建立高阶似然函数以增加AVS可跟踪声目标数目. 与MUSIC空间谱[7-8] 获取手段相同, 对式(6)累积量矩阵
$ {\boldsymbol{R}}_{Z,k}^{(\rm{cum})} $ 进行特征分解, 获取噪声子空间$ {\boldsymbol{U}}_{k,v} $ 和信号子空间$ {\boldsymbol{U}}_{k,s} $ . 基于信号、噪声子空间的正交性, 声目标DOA可由${\hat {\boldsymbol{\theta }}_k} = \arg \mathop {\min }\nolimits_{{{\boldsymbol{\theta }}_k}} {({\boldsymbol{b}}({\boldsymbol{\eta }}{{\boldsymbol{\theta }}_k}))^{\mathrm H}}{{\boldsymbol{U}}_{k,v}}{\boldsymbol{U}}_{k,v}^{\mathrm H}({\boldsymbol{b}}({\boldsymbol{\eta }}{{\boldsymbol{\theta }}_k}))$ 获取. 因此, 定义$ k $ 时刻AVS虚拟扩展的高阶似然函数为:$$ \begin{equation} g({\boldsymbol{Z}}_k|{\boldsymbol{\theta}}_k)=\frac{1} {({\boldsymbol{b}}({\boldsymbol{\eta}}{\boldsymbol{\theta}}_k))^{\mathrm H}{\boldsymbol{U}}_{k,v}{\boldsymbol{U}}_{k,v}^{\mathrm H} ({\boldsymbol{b}}({\boldsymbol{\eta}}{\boldsymbol{\theta}}_k))} \end{equation} $$ (7) 其中
$ {\boldsymbol{ \eta}}={\mathrm {diag}}\{1,0,0,1,0,0\} $ 用于提取声目标状态中的DOA信息,$ {\boldsymbol{b}}({\boldsymbol{\eta}}{\boldsymbol{\theta}}_k)={\boldsymbol{a}}({\boldsymbol{\eta}}{\boldsymbol{\theta}}_k)\otimes{\boldsymbol{a}}^\ast({\boldsymbol{\eta}}{\boldsymbol{\theta}}_k) $ . 高阶似然函数$ g({\boldsymbol{Z}}_k|{\boldsymbol{\theta}}_k) $ 是AVS虚拟扩展后的空间谱, 表示在量测数据为$ {\boldsymbol{Z}}_k $ 的条件下目标状态$ {\boldsymbol{\theta}}_k $ 的概率, 形式上可近似理解为概率密度函数.式(7)所示高阶似然函数与粒子滤波似然函数具有类似效果, 当目标状态接近真实声目标状态时, 似然函数很大, 否则很小. 尽管似然函数中并未直接体现式(3)量测模型中的量测
$ {\boldsymbol{Z}}_k $ 和噪声$ {\boldsymbol{V}}_k $ , 但用于构建高阶似然的$ {\boldsymbol{U}}_{k,v} $ 是通过$ {\boldsymbol{Z}}_k $ 的累积量矩阵$ {\boldsymbol{R}}_{Z,k}^{(\rm{cum})} $ 获取, 且其估计结果受$ {\boldsymbol{V}}_k $ 的影响. 图1为不同信噪比(Signal-to-noise ratio, SNR)情况下的高阶似然函数示例. 当SNR较高时谱峰尖锐, 容易获得更准确的DOA, 而SNR较低时谱峰主瓣过宽, 将导致DOA估计值偏离真实值.
图 1 不同信噪比下的归一化高阶似然函数示例Fig. 1 Example of normalized higher-order spatial spectrum under different SNR2.2 Cum-AMMS-GLMB滤波算法
GLMB滤波器可解决MeMBer、CPHD、PHD等滤波方法还需额外的航迹关联步骤才可解决的航迹区分问题[21-23]. 作为GLMB滤波器的特殊形式,
$ \delta $ -GLMB滤波器由于其能够降低GLMB滤波器的内存需求和计算成本, 对于目标跟踪更具优势. 一个具有状态空间$ {\boldsymbol{\Theta}} $ 和离散标签空间$ {L} $ 的$ \delta $ -GLMB RFS的概率密度可表示为[22]:$$ \begin{equation} \pi ({\boldsymbol{\Theta }}) = \Delta ({\boldsymbol{\Theta }})\sum\limits_{(I,\xi ) \in { {F}} ( L) \times \Xi } {{w^{(I,\xi )}}{\delta _I}({{L}} ({\boldsymbol{\Theta }})){{[{p^{(\xi )}}]}^{\boldsymbol{\Theta }}}} \end{equation} $$ (8) 其中
$ { {L}} ({\boldsymbol{\Theta}}) = \{{ {L}} ({\boldsymbol{\theta}}):{\boldsymbol{\theta}}\in {\boldsymbol{\Theta}}\} $ 表示$ {\boldsymbol{ \Theta}} $ 的标签集.$\Delta({\boldsymbol{\Theta}}) = \delta_{|{{\boldsymbol{ \Theta}}}|}(|{ {L}} ({{\boldsymbol{ \Theta}}})|)$ 表示不同的标签指示器,$ |\cdot| $ 表示集合的势. 每个$ I\in { {F}} ({L}) $ 表示$ k $ 时刻的一组轨迹标签.$ \Xi $ 为一离散空间, 每个$ \xi \in \Xi $ 表示$ k $ −1时刻航迹−量测的关联映射历史.$ (I,\xi) $ 表示一个假设, 即轨迹集合$ I $ 具有关联映射历史$ \xi $ .$ w^{({I,\xi})} $ 表示假设$ (I,\xi) $ 的权重,$ p^{(\xi)}(\cdot,l) $ 为关联映射历史$ \xi $ 中轨迹$ l $ 的运动状态的概率密度. M$ \delta $ -GLMB滤波器为一种$ \delta $ -GLMB滤波器近似方法, 相当于对$ \delta $ -GLMB中关联映射历史执行边缘化处理, 其概率密度可表示为[36]:$$ \begin{equation} \hat\pi ({\boldsymbol{\Theta }}) = \Delta ({\boldsymbol{\Theta }})\sum\limits_{I \in { {F}} ( L) } {{w^{(I)}}{\delta _I}({ {L}} ({\boldsymbol{\Theta }})){{[{p^{(I )}}]}^{\boldsymbol{\Theta }}}} \end{equation} $$ (9) 其中
$$ \begin{equation} \begin{aligned} &w^{(I)}=\sum_{\xi\in \Xi}w^{(I,\xi)}\\ &p^{(I)}({\boldsymbol{\theta}},l)=1_I(l)\frac{1}{w^{(I)}}\sum_{\xi\in \Xi} {w^{(I,\xi)}p^{(\xi)}({\boldsymbol{\theta}},l)}\\ \end{aligned} \end{equation} $$ $ \delta $ -GLMB RFS中权重$ w^{(I,\xi)} $ 和概率$ p^{(\xi)} $ 的存储容量分别为$ |{ {F}} ({L})\times \Xi| $ 和$ |\Xi | $ , 而M$ \delta $ -GLMB RFS中权重$ w^{(I)} $ 和概率$ p^{(I)} $ 的存储容量均为$ |{ {F}} ({L})| $ . 可见M$ \delta $ -GLMB滤波器在保留多目标后验关键统计量的同时, 能够降低存储容量和计算成本. 考虑到式(3)量测模型为所有声目标的混合信号, 任何一个轨迹假设集合均与该量测映射. 因此, 目标跟踪过程中无需引入航迹−量测的关联映射历史, M$ \delta $ -GLMB滤波框架更适用于基于混合量测的声目标跟踪算法设计.多目标跟踪问题中, 通常选取理想的匀速(Constant velocity, CV)、匀加速(Constant acceleration, CA)、协同转弯(Coordinate turn, CT) 等运动模型对目标运动方式进行建模, 然而单一运动模型无法准确描述真实目标运动状态. 为此, 本部分在M
$ \delta $ -GLMB滤波框架下引入IMM思想, 将包含轨迹标签的目标状态$ ({\boldsymbol{\theta}},l) $ 扩展为一额外增加运动模型索引的增广状态$ {\boldsymbol{\theta}}=({\boldsymbol{\theta}},l,o) $ , 设计并推导给出Cum-AMMS-GLMB算法. 算法通过不断更新多种运动模型下的加权混合目标状态, 实现对多机动声目标运动状态的拟合.2.2.1 Cum-AMMS-GLMB预测器
给定当前时刻的多目标状态
$ {\boldsymbol{\Theta}} $ , 及其标签空间$ {L} $ 和运动模型索引空间$ O $ . 每个声目标状态$ ({\boldsymbol{\theta}},l,o)\in {{\boldsymbol{ \Theta}}} $ 要么以$ p_{s}({\boldsymbol{\theta}},l,o) $ 的概率存活到下一时刻, 并以$ f^{(o)}({\boldsymbol{\theta}}_+|{\boldsymbol{\theta}},l)p(o_+|o)\delta_l(l_+) $ 的转移函数转移到新的状态$ ({\boldsymbol{\theta}}_+,l_+,o_+) $ , 要么以$ q_ {s}({\boldsymbol{\theta}},l,o)=1-p_ {s}({\boldsymbol{\theta}},l,o) $ 的概率消失. 幸存目标的状态集合可按照如下分布而得:$$ \begin{split} &{f_S}({{\boldsymbol{\Theta }}^{(s)}}|{\boldsymbol{\Theta }}) =\\ &\qquad \Delta ({{\boldsymbol{\Theta }}^{(s)}})\Delta ({\boldsymbol{\Theta }}){1_{{ {L}} ({\boldsymbol{\Theta }})}}({ {L}} ({{\boldsymbol{\Theta }}^{(s)}})){[\Phi ({{\boldsymbol{\Theta }}^{(s)}}; \cdot )]^{{{\boldsymbol{\Theta }}}}} \end{split} $$ (10) 其中
$$ \begin{equation} \nonumber \begin{split} &{{\Phi }}({{\boldsymbol{\Theta }}^{(s)}};{\boldsymbol{ \theta}} ,o,l) = \\ &\quad\sum\limits_{({\boldsymbol{\theta} _ + },{o_ + },{l_ + }) \in {{\boldsymbol{\Theta }}^{(s)}}} {\Big( {[1 - {1_{{ {L}} ({{\boldsymbol{\Theta }}^{(s)}})}}(l)]}\times }{q_s}({\boldsymbol{\theta}} ,o,l) + \\ &\quad{ {\delta _l}({l_ + }){p_s}({\boldsymbol{\theta}} ,o,l)f^{(o)}({\boldsymbol{\theta} _ + }\mid {\boldsymbol{\theta}} ,l)p({o_ + }\mid o)} \Big) \end{split} \end{equation} $$ 值得注意的是, 上述状态转移过程中仅将状态
$ {\boldsymbol{\theta}} $ 转移到状态$ {\boldsymbol{\theta}}_+ $ , 而标签被保留$ l_+=l $ , 运动模型索引$ \{o,o_+\} $ 作为隐含状态仅用于对多种运动模型下状态的加权混合求解, 并不显式体现. 对于标签空间为$ {B} $ 的新生目标状态集合$ {{\boldsymbol{\Theta }}^{(b)}} $ , 其分布如下:$$ \begin{equation} { \begin{aligned} {f_B}({{\boldsymbol{\Theta }}^{(b)}})= \Delta ({{\boldsymbol{\Theta }}^{(b)}}){\sum\limits_{J \in { {F}} ( B)} {{w_B}(J){\delta _J}({ {L}} ({{\boldsymbol{\Theta }}^{(b)}})){{[{p_B}]} ^{{{\boldsymbol{\Theta }}^{(b)}}}}}} \end{aligned} } \end{equation} $$ (11) 仅考虑声目标的幸存、新生、死亡. 下一时刻的多目标状态集
$ {{\boldsymbol{ \Theta}}}_+ $ 为幸存目标$ {{\boldsymbol{\Theta }}^{(s)}} $ 和新生目标$ {{\boldsymbol{\Theta }}^{(b)}} $ 的并集, 其中标签空间为$ {L}_+={L}\cup{B} $ , 运动模型索引空间始终为$ O $ .假设多目标先验分布为式(9)所示的M
$ \delta $ -GLMB形式, 则多目标预测同样为M$ \delta $ -GLMB形式:$$ \begin{equation} \pi ({{\boldsymbol{\Theta }}_ + }) = \Delta ({{\boldsymbol{\Theta }}_ + })\sum\limits_{{I_ + } \in { {F}} ({ L}_+) } {w_ + ^{({I } )}{\delta _{{I_ + }}}({ {L}} ({{\boldsymbol{\Theta }}_ + })){{[p_ + ^{(I )}]}^{{{\boldsymbol{\Theta }}_ + }}}} \end{equation} $$ (12) 其中,
$$ \begin{equation} \begin{split} &w_+^{(I)}=w_B(I_+\cap B)w_S^{(I)}(I_+\cap L)\\ &w_S^{(I )}(L) = {[\eta _S^{(I )}]^L}\sum\limits_{Y \subseteq L} {{1_Y}(L){{[1 - \eta _S^{(I )}]}^{Y - L}}{w^{(I )}}}\\ &p_+^{(I)}({\boldsymbol{\theta}}_+,o_+,l)=\\ &\qquad 1_{L}(l)p_S^{(I)}({\boldsymbol{\theta}}_+,o_+,l)+1_{B}(l)p_B(({\boldsymbol{\theta}}_+,l))\\ &p_S^{(I)}({\boldsymbol{\theta}}_+,o_+,l)=\\ &\ \frac{\sum\limits_{o \in O} {\langle {{p_s}( \cdot ,o,l)f^{(o)}({\boldsymbol{\theta}_ + }\mid \cdot,l)p({o_ + }\mid o),{p^{(I )}}( \cdot ,o,l)} \rangle } }{\eta_S^{(I)}(l)}\\ &\eta_S^{(I)}(l)=\sum\limits_{o \in O} {\langle {{p_s}( \cdot ,o,l),{p^{(I )}}( \cdot ,o,l)} \rangle } \\[-15pt]\end{split} \end{equation} $$ (13) 不同于M
$ \delta $ -GLMB算法, Cum-AMMS-GLMB算法预测器中的$ p_S^{(I)}({\boldsymbol{\theta}}_+,o_+,l) $ 和$ \eta_S^{(I)}(l) $ 则是通过多种运动模型概率和各模型下目标状态的加权混合获取.证明. 假设声目标之间独立进化, 新生目标与幸存目标无关, 则多目标预测为:
$$ \begin{equation} \nonumber \pi_+({{\boldsymbol{ \Theta}}}_+)=f_B({{\boldsymbol{ \Theta}}}^{(b)})\int {{f_S}({{\boldsymbol{ \Theta}}}^{(s)}\mid {{\boldsymbol{ \Theta}}})\pi ({{\boldsymbol{ \Theta}}})\mathrm{d} {{\boldsymbol{ \Theta}}}} \end{equation} $$ 根据式(10)可得,下一时刻幸存多目标的状态密度
$ {\pi _s}({{\boldsymbol{\Theta }}^{(s)}}) = \int {{f_S}({{\boldsymbol{\Theta }}^{(s)}}\mid {\boldsymbol{\Theta }})\pi ({\boldsymbol{\Theta }})\mathrm{d} {\boldsymbol{\Theta }}} $ 为:$$ \begin{equation} \nonumber { \begin{aligned} {\pi _s}({{\boldsymbol{\Theta }}^{(s)}}) =\;& \Delta ({{\boldsymbol{\Theta }}^{(s)}})\times \sum\limits_{I \in { {F}} (L) } \sum\limits_{Y \subseteq {L}} {{1_Y}({ {L}} ({{\boldsymbol{\Theta }}^{(s)}}))} {w^{(I )}}\times\\ &\prod_{l\in Y}\left \langle \Phi({\boldsymbol{\Theta} }^{(s)};\cdot,l)p^{(I)}(\cdot,l) \right \rangle \end{aligned} } \end{equation} $$ 由于增广状态值为离散索引值, 因此可得:
$$ \begin{equation} \nonumber { \begin{split} &{\pi _s}({{\boldsymbol{\Theta }}^{(s)}}) = \Delta ({{\boldsymbol{\Theta }}^{(s)}})\sum\limits_{I \in { {F}} (L) }\Big( { {\sum\limits_{Y \subseteq L} {{1_{Y}}({ {L}} ({{\boldsymbol{\Theta }}^{(s)}})){w^{(I )}}} } \times} \\ &\qquad{\prod\limits_{l \in Y} {\sum\limits_{o \in O} {\langle {\Phi ({{\boldsymbol{\Theta }}^{(s)}}; \cdot ,o,l),{p^{(I )}}( \cdot ,o,l)} \rangle } } } \Big)=\\ &\qquad\Delta ({{\boldsymbol{\Theta }}^{(s)}})\sum\limits_{I \in { {F}} ( L) } \Bigg({\sum\limits_{Y \subseteq L} {{1_{Y}}({ {L}} ({{\boldsymbol{\Theta }}^{(s)}})){w^{(I )}}} \times} \end{split} } \end{equation} $$ $$ \begin{equation} \nonumber { \begin{split} &\prod\limits_{l \in { {L}} ({\boldsymbol{\Theta }}^{(s)})} {\sum\limits_{({\boldsymbol{\theta}_ + },{l_ + }) \in {{\boldsymbol{\Theta }}^{(s)}}} {\sum\limits_{o \in O} \Big({{\delta _l}({l_ + }) } } \times} \\ &\langle {{p_s}( \cdot ,o,l)f^{(o)}({\boldsymbol{\theta} _ + }| \cdot ,l)p(o_+|o),{p^{(I )}}( \cdot ,o,l)} \rangle \Big)\times\\ &\prod\limits_{l \in Y - { {L}} ({{\boldsymbol{\Theta }}^{(s)}})} {\sum\limits_{o \in O} {\langle {1 - {p_s}( \cdot ,o,l),{p^{(I )}}( \cdot ,o,l)} \rangle } }\Bigg)= \\ &\Delta ({{\boldsymbol{\Theta }}^{(s)}})\sum\limits_{I \in { {F}} ( L) } \Big({\sum\limits_{Y \subseteq L} {{1_{Y}}({ {L}} ({{\boldsymbol{\Theta }}^{(s)}})){w^{(I)}}} }\times \\ &\prod\limits_{({\boldsymbol{\theta} _ + },l) \in {{\boldsymbol{\Theta }}^{(s)}}} {p_S^{(I )}({\boldsymbol{\theta} _ + },{o_ + },l)\eta _S^{(I )}(l)}\times\\ &\prod\limits_{l \in Y - { {L}} ({{\boldsymbol{\Theta }}^{(s)}})} {(1 - \eta _S^{(I )}(l))}\Big)= \\ &\Delta ({{\boldsymbol{\Theta }}^{(s)}})\sum\limits_{I \in { {F}} ( L)} {w_S^{(I )}({ {L}} ({{\boldsymbol{\Theta }}^{(s)}}){{[p_S^{(I )}]}^{{{\boldsymbol{\Theta }}^{(s)}}}}} \end{split} } \end{equation} $$ 其中
$ w_S^{(I)}(L) $ ,$ \eta_S^{(I)}(l) $ ,$ p_S^{(I)}({\boldsymbol{\theta}}_+,o_+,l) $ 如式(13)所示.记
$ {{\boldsymbol{ \Theta}}}^{(b)}={{\boldsymbol{ \Theta}}}_+-{{\boldsymbol{ \Theta}}}\times L $ ,$ {{\boldsymbol{ \Theta}}}^{(s)}={{\boldsymbol{ \Theta}}}_+ \cap({{\boldsymbol{ \Theta}}}\times L) $ , 则多目标预测$ {\pi _ + }({{\boldsymbol{\Theta }}_ + }) = {f_B}({{\boldsymbol{\Theta }}^{(b)}}){\pi _s}({{\boldsymbol{\Theta }}^{(s)}}) $ 可表示为:$$ \begin{equation} \nonumber { \begin{aligned} {\pi _ + }({{\boldsymbol{\Theta }}_ + })=\;& \Delta ({{\boldsymbol{\Theta }}_ + })\times\sum\limits_{J\in { {F}} ( B)}{w_B(J)\delta_J({ {L}} ({{\boldsymbol{ \Theta}}}^{(b)}))}\times\\ &\sum\limits_{I \in { {F}} ( L) } {w_S^{(I )}({ {L}} ({{\boldsymbol{\Theta }}^{(s)}})){{[p_ + ^{(I )}]}^{{{\boldsymbol{\Theta }}_ + }}}} \end{aligned} } \end{equation} $$ 其中
$[p_ + ^{(I)}] $ 如式(13)所示. 记$ { {L}} ({{\boldsymbol{\Theta }}^{(b)}})={ {L}} ({{\boldsymbol{ \Theta}}}_+)- L $ ,$ { {L}} ({{\boldsymbol{\Theta }}^{(s)}})={ {L}} ({{\boldsymbol{ \Theta}}}_+)\cap L $ . 式(13)中$ w_S^{(I)}({ {L}} ({{\boldsymbol{ \Theta}}}^{(s)})) $ 还可描述为$ \sum_{Y' \subseteq I} {{\delta _{Y'}}({ {L}} ({{\boldsymbol{\Theta }}^{(s)}})){{[\eta _S^{(I )}]}^{Y'}}{{[1 - \eta _S^{(I )}]}^{I - Y'}}{w^{(I)}}} $ , 因此可得:$$ \begin{equation} \nonumber \begin{aligned} &{\pi _ + }({{\boldsymbol{\Theta }}_ + }) =\Delta ({{\boldsymbol{\Theta }}_ + })\sum\limits_{I \in { {F}} ( L)} {w^{(I )}}\sum\limits_{J \in { {F}} ( B)} \sum\limits_{Y' \in { {F}} ( L)}\\ &\quad{\Big( {{1_I}(Y')} }{\delta _{Y' \cup J}}({ {L}} ({{\boldsymbol{\Theta }}_ + })){w_B}(J)\times\\ &\quad{{[\eta _S^{(I )}]}^{Y'}}{{[1 - \eta _S^{(I )}]}^{I - Y'}}{{[p_ + ^{(I )}]}^{{{\boldsymbol{\Theta }}_ + }}} \Big)=\\ &\quad\Delta ({{\boldsymbol{\Theta }}_ + })\sum\limits_{J \in { {F}} ( B)} \sum\limits_{Y' \in { {F}} ( L)}\Big({ {\delta _{Y' \cup J}}({ {L}} ({{\boldsymbol{\Theta }}_ + })) }\times\\ &\quad{w_B}(J){w_S^{(I )}(Y'){{[p_ + ^{(I )}]}^{{{\boldsymbol{\Theta }}_ + }}}}\Big)=\\ &\quad\Delta ({{\boldsymbol{\Theta }}_ + })\sum\limits_{{I_ + } \in { {F}} ({ L_ + }) } {{\delta _{{I_ + }}}({ {L}} ({{\boldsymbol{\Theta }}_ + }))w_ + ^{({I} )}{{[p_ + ^{(I )}]}^{{{\boldsymbol{\Theta }}_ + }}}} \end{aligned} \end{equation} $$ 其中
$ J=I_+\cap B $ ,$ Y'=I_+\cap L $ , 以及$ J\cup Y'=I_+ $ . 最终可得式(12)所示预测器. □2.2.2 Cum-AMMS-GLMB更新器
假设多目标预测为式(9)所示的M
$ \delta $ -GLMB形式, 则增广运动模型参数下的多目标后验具有如下M$ \delta $ -GLMB形式:$$ \begin{equation} \pi ({\boldsymbol{\Theta }}|{\boldsymbol{Z}}) = \Delta ({\boldsymbol{\Theta }})\sum\limits_{I \in { {F}} ( L) } {w_Z^{(I )}{\delta _I}({ {L}} ({\boldsymbol{\Theta }})){{[p_Z^{(I )}( \cdot |{\boldsymbol{Z}})]}^{\boldsymbol{\Theta }}}} \end{equation} $$ (14) 其中,
$$ \begin{equation} \begin{split} &w_Z^{(I)} \propto{{{w^{(I )}}{{[\eta _Z^{(I )}]}^I}}}\\ &p_Z^{(I )}({\boldsymbol{\theta }},o,l|{\boldsymbol{Z}}) = \frac{{{p^{(I )}}({\boldsymbol{\theta }},o,l){p_D}({\boldsymbol{\theta }},l)g({\boldsymbol{Z}}|{\boldsymbol{\theta }},l)}}{{\eta _Z^{(I)}(l)}}\\ &\eta _Z^{(I )}(l) = \sum\limits_{o \in O} {\langle {{p^{(I )}}({\boldsymbol{\theta }},o,l),{p_D}({\boldsymbol{\theta }},l)g({\boldsymbol{Z}}|{\boldsymbol{\theta }},l)} \rangle } \end{split} \end{equation} $$ (15) 上述更新器特点如下: 1) 算法量测模型为不同状态声目标的混合信号, 无需引入航迹−量测的关联映射; 2) 式(15)中的检测概率
$ {p_D}({\boldsymbol{\theta }},l) $ 和高阶似然函数$ g({\boldsymbol{Z}}|{\boldsymbol{\theta }},l) $ 已包含目标是否被检测, 无需再引入漏检概率$ {q_D}({\boldsymbol{\theta }},l) $ ; 3) 不同于M$ \delta $ -GLMB更新器,$ \eta _Z^{(I )}(l) $ 由多种运动模型的加权混合获取.证明. 量测
$ {\boldsymbol{Z}} $ 的叠加性导致其与目标状态具有较强相关性, 理论上无法分解多目标似然函数. 然而高阶似然基于噪声、信号子空间的正交性构建, 量测$ {\boldsymbol{Z}} $ 可等效为噪声子空间量测$ {\boldsymbol{U}}_{k,v} $ . 等效的量测可认为不具有“叠加性”, 因此多目标似然函数依然能够分解为$ g({\boldsymbol{Z}}|{{\boldsymbol{ \Theta}}}) =[p_D(\cdot)g({\boldsymbol{Z}}|\cdot)]^{\boldsymbol{\Theta}} $ . 虽然$ g({\boldsymbol{Z}}|\cdot) $ 中已包含杂波和目标状态的检测率, 但为增强高似然区域粒子, 抑制杂波粒子, 更新器中依然引入检测概率$ p_D(\cdot) $ . 对于$ g({\boldsymbol{Z}}|{{\boldsymbol{ \Theta}}})\pi({\boldsymbol{\Theta}}) $ , 有:$$ \begin{equation} \nonumber \begin{split} &g({\boldsymbol{Z}}|{\boldsymbol{\Theta }})\pi ({\boldsymbol{\Theta }}) = \\ &\qquad\Delta ({\boldsymbol{\Theta }})\sum\limits_{I\in { {F}} ( L) } {{w^{(I )}}{\delta _I}({ {L}} ({\boldsymbol{\Theta }})){{[\eta _Z^{(I )}]}^{{ {L}} ({\boldsymbol{\Theta }})}}{{[p_Z^{(I )}( \cdot |{\boldsymbol{Z}})]}^{\boldsymbol{\Theta }}}} \end{split} \end{equation} $$ 其中
$ \eta_Z^{(I)}(l) $ 和$ p_Z^{(I)}({\boldsymbol{\theta}},r,l|{\boldsymbol{Z}}) $ 如式(15)所示. 对上式求积分:$$ \begin{equation} \nonumber \begin{split} &\int {g({\boldsymbol{Z}}|{\boldsymbol{\Theta }}) \pi ({\boldsymbol{\Theta }})\mathrm{d} {\boldsymbol{\Theta }}}= \\ &\qquad\sum\limits_{I \in { {F}} ( L) } {\sum\limits_{Y \subseteq L} {{w^{(I )}}{\delta _I}(Y){{\left[ {\int {\eta _Z^{(I )}p_Z^{(I )}( \cdot |{\boldsymbol{Z}})\mathrm{d} {\boldsymbol{\Theta }}} } \right]}^Y}} }= \\ &\qquad\sum\limits_{I \in { {F}} ( L) } {\sum\limits_{Y \subseteq L} {{w^{(I )}}{\delta _I}(Y)} } \times\\ &\qquad{\left[ {\sum\limits_{o \in O} {\langle {{p^{(I )}}( \cdot ,o,l),{p_D}( \cdot ,l)g({\boldsymbol{Z}}| \cdot ,l)} \rangle } } \right]^Y}=\\ &\qquad\sum\limits_{I \in { {F}} ( L) } {\sum\limits_{Y \subseteq L} {{w^{(I)}}{\delta _I}(Y){{\left[ {\eta _Z^{(I)}} \right]}^Y}} } \end{split} \end{equation} $$ 联合上述两式, 可得多目标后验为:
$$ \begin{equation} \nonumber \begin{aligned} \pi ({\boldsymbol{\Theta }}|{\boldsymbol{Z}})=\; & \frac{{g({\boldsymbol{Z}}|{\boldsymbol{\Theta }})\pi ({\boldsymbol{\Theta }})}}{{\int {g({\boldsymbol{Z}}|{\boldsymbol{\Theta }})\pi ({\boldsymbol{\Theta }})\mathrm{d} {\boldsymbol{\Theta }}} }}=\\ &\Delta ({\boldsymbol{\Theta }})\sum\limits_{I \in { {F}} ( L) } {\frac{{{w^{(I )}}{{[\eta _Z^{(I )}]}^I}}}{{\sum\limits_{I \in { {F}} ( L) } {\sum\limits_{Y \subseteq L} {{w^{(I )}}{\delta _I}(Y){{\left[ {\eta _Z^{(I )}} \right]}^Y}} } }}}\times\\ &{\delta _I}({ {L}} ({\boldsymbol{\Theta }})){[p_Z^{(I )}( \cdot |{\boldsymbol{Z}})]^{\boldsymbol{\Theta }}}=\\ &\Delta ({\boldsymbol{\Theta }})\sum\limits_{I \in { {F}} ( L) } {w_Z^{(I )}{\delta _I}({{L}} ({\boldsymbol{\Theta }})){{[p_Z^{(I)}( \cdot |{\boldsymbol{Z}})]}^{\boldsymbol{\Theta }}}} \end{aligned} \end{equation} $$ 最终可得式(14)所示更新器, 其中
$ w_Z^{(I )} $ 如式(15)所示. 注意: 不同时刻检测概率不同. 后续算法实现过程中, 将$ {p_{D,k}}(\cdot) $ 表示为$ k $ 时刻的检测概率. □3. Cum-AMMS-GLMB算法实现
式(12) ~ 式(15)即为Cum-AMMS-GLMB算法的预测器和更新器. 本部分在基于AVS的声目标DOA跟踪场景下, 设计新生目标提取和检测概率拟合方法. 由于量测模型为非线性模型, 最终通过序贯蒙特卡洛(Sequential monte carlo, SMC)实现算法框架下的多机动声目标跟踪.
3.1 新生目标提取
$ k $ 时刻是否存在目标新生可依据$ k-1 $ 时刻和$ k $ 时刻量测值对应的空间谱判定. 设定判断阈值$ thr $ , 分别寻找$ k -1$ 时刻和$ k $ 时刻的归一化高阶空间谱中峰值集合$ {\boldsymbol{p}}_{\rm{peak}}^{(k-1)} $ 和$ {\boldsymbol{p}}_{\rm{peak}}^{(k)} $ 大于$ thr $ 对应的目标状态集合:$$ \begin{equation} \begin{split} &{{\boldsymbol{ \Theta}}}_{k-1}^{(P)}=\{\boldsymbol{\theta}^{(k-1)}\mid {\boldsymbol{p}}_{\rm{peak}}^{(k-1)}({\boldsymbol{\theta}}^{(k-1)})>thr\}\\ &{{\boldsymbol{ \Theta}}}_{k}^{(P)}=\{\boldsymbol{\theta}^{(k)}\mid {\boldsymbol{p}}_{\rm{peak}}^{(k)}({\boldsymbol{\theta}}^{(k)})>thr\} \end{split} \end{equation} $$ (16) 通过比较状态集合
$ {{\boldsymbol{ \Theta}}}_{k}^{(P)} $ 是否在$ {{\boldsymbol{ \Theta}}}_{k-1}^{(P)} $ 的基础上新增状态, 即可确认$ k $ 时刻是否产生新生目标$ {\boldsymbol{\theta}}_{B,k} $ .新生目标的新生概率可由该状态对应的归一化空间谱幅值表示. 为防止SNR较低情况下杂波产生的虚假源对声目标跟踪的影响, 引入接受系数
$ \epsilon_B $ , 并且令$ r_B=\epsilon_B\cdot {\boldsymbol{p}}_{\rm{peak}}^{(k)}({\boldsymbol{\theta}}_{B,k}) $ 作为新生概率.$ \epsilon_B $ 取一较小值即可, 本文设置为$ \epsilon_B=0.3 $ .3.2 检测概率模型拟合
不同声目标状态对应的归一化高阶空间谱幅值可体现检测概率. 空间谱类似于混合二维高斯分布, 因此
$ k $ 时刻的检测概率模型可由混合二维高斯分布拟合:$$ \begin{equation} p_{D,k}({\boldsymbol{\theta}})=\sum\limits_{i\in M^{(p)}(k)}{\boldsymbol{p}_{\rm{peak}}^{(k)}(i){\rm{N}}({\boldsymbol{\mu}}_\theta^i,{\boldsymbol{\varSigma}}_\theta^i)} \end{equation} $$ (17) 其中
$ M^{(p)}(k) $ 表示当前空间谱中大于阈值$ thr $ 的峰值数目,$ {\boldsymbol{p}}_{\rm{peak}}^{(k)}(i) $ 为$ k $ 时刻空间谱峰值中第$ i $ 个峰值对应的幅值,$ {\boldsymbol{\mu}}_\theta^{i} $ 和$ {\boldsymbol{\varSigma}}_\theta^i $ 分别表示第$ i $ 个峰值对应的声目标状态和协方差矩阵. 图2即为对图1(a)所示归一化空间谱通过式(17)拟合获得的检测概率函数.3.3 Cum-AMMS-GLMB算法的SMC实现
$ k-1 $ 时刻多目标概率密度参数集为$\{(L_{k-1}^{(h)}, w_{k-1}^{(h)}), p_{S,k-1}^{(h)}(\cdot,o,l)\}_{h=1}^{H_{k-1}}$ , 其中$ L_{k-1}^{(h)} $ 为第$ h $ 个假设的标签集合,$ H_{k-1} $ 为总假设数目. 标签$ l\in L_{k-1}^{(h)} $ 的单目标状态$ p_{S,k-1}^{(h)}(\cdot,o,l) $ 由一组加权粒子集$\{w_{k-1,i}(l), {\boldsymbol{\theta}}_{k-1,i}(o,l)\}_{i=1}^{J_{k-1}^{(l)}}$ 给出, 其中$ J_{k-1}^{(l)} $ 为总粒子数目. 以下列出Cum-AMMS-GLMB算法的SMC实现过程.3.3.1 初始化输入
$ k -1$ 时刻样本$\big\{w_{k-1,i}(l),{\boldsymbol{\theta}}_{k-1,i}(o,l)\big\}_{i=1}^{J_{k-1}^{(l)}}, l\in {L}_{k-1}$ ,$ k $ 时刻量测数据$ {\boldsymbol{Z}}_k $ , 判定阈值$ thr $ , 各运动模型$ {\boldsymbol{F}}^{(o)} $ 和$ {\boldsymbol{G}}^{(o)} $ , 运动模型概率$ \left\{{{P_{k-1,i}(o)}}\right\}_{i=1}^{J_{k-1}^{(l)}} $ , 及运动模型转移概率$p(o=i|o=j),i,j \in O$ .3.3.2 预处理
对量测数据进行高阶累积量预处理, 抑制量测高斯噪声和虚拟扩展AVS. 具体步骤如下:
步骤 1. 通过式(5)、(6)对
$ {\boldsymbol{Z}}_k $ 进行高阶累积量预处理获取累积量矩阵$ {\boldsymbol{R}}_{Z,k}^{(\rm{cum})} $ , 并对$ {\boldsymbol{R}}_{Z,k}^{(\rm{cum})} $ 特征分解获取噪声子空间$ {\boldsymbol{U}}_{k,v} $ .步骤 2. 通过式(7)求解归一化高阶空间谱
$ {\boldsymbol{P}}_{\rm{cum}}^{(k)} $ , 寻找$ {\boldsymbol{P}}_{\rm{cum}}^{(k)} $ 中大于阈值$ thr $ 的峰值$ {\boldsymbol{p}}_{\rm{peak}}^{(k)} $ 用于式(17) 所示检测概率函数$ p_{D,k}({\boldsymbol{\theta}}) $ 的拟合.3.3.3 预测器
预测器输出
$ \left\{ {w_{k|k - 1,i}(l),{\boldsymbol{\theta }}_{k|k - 1,i}(o,l)} \right\}_{i = 1}^{J_{k|k - 1}^{(l )}} $ 为新生、幸存粒子的并集, 其中$ J_{k|k - 1}^{(l)} = J_{B,k}^{(l)} + J_{k - 1}^{(l)} $ . 具体步骤如下:步骤 3. 新生粒子
$ \{{w_{B,k,i}(l),{\boldsymbol{\theta}}_{B,k,i}(l)}\}_{i=1}^{J_{B,k}^{(l)}} $ :根据(16)判定是否存在新增目标状态
$ {\boldsymbol{\theta}}_{B,k} $ , 并生成粒子$ {\boldsymbol{\theta}}_{B,k,i} \sim {\boldsymbol{\theta}}_{B,k} $ 和权重$ w_{B,k,i} \sim {1}/{J_{B,k}^{(l)}} $ , 其中$ {J_{B,k}^{(l)}} $ 为采样粒子个数.步骤 4. 各运动模型下目标状态预测值
$ {\boldsymbol{\theta}}_{k|k-1,i}(o,l) $ :$$ {\boldsymbol{\theta}}_{k|k-1,i}(o,l)={\boldsymbol{F}}^{(o)}{\boldsymbol{\theta}}_{k-1,i}(o,l)+{\boldsymbol{G}}^{(o)}{\boldsymbol{w}}_k,o\in O $$ 步骤 5. 运动模型概率预测
$ P_{k|k-1,i} $ 和更新$ P_{k,i} $ :$$ \begin{aligned} &{P_{k|k - 1,i}}(o = j) = \sum\limits_{h \in O} {p (o = j|o = h){P_{k - 1,i}}(o= h)} \\ &{P_{k,i}}(o = j) = \frac{{g({{\boldsymbol{Z}}_k}|{\boldsymbol{\theta}}_{k|k - 1,i}(o,l)){P_{k|k - 1,i}}(o = j)}}{{\sum\limits_{o \in O} {g({{\boldsymbol{Z}}_k}|{\boldsymbol{\theta}}_{k|k - 1,i}(o,l)){P_{k|k - 1,i}}(o = h)} }} \end{aligned} $$ 步骤 6. 幸存粒子
$\{w_{S,k|k-1,i}(l),{\boldsymbol{\theta}}_{S,k|k-1,i}(o, l)\}_{i=1}^{J_{k-1}^{(l)}}$ :$$ \begin{aligned} &{\boldsymbol{\theta }}_{S,k|k - 1,i}(o,l) = \sum\limits_{j \in O} {{P_{k,i}}(o = j){\boldsymbol{\theta }}_{k|k - 1,i}(o,l)} \\ &w_{S,k|k - 1,i}(o,l) = \frac{{{p_{s,k}}w_{k - 1,i}(l)}}{{\sum\limits_{i = 1}^{J_{k - 1}^{(l)}} {{p_{s,k}}w_{k - 1,i}(l)} }} \end{aligned} $$ 根据以上预测器, 可获得
$ k $ 时刻的预测概率密度参数集为$ \{ {( {L_{k|k - 1}^{(h)},w_{k|k - 1}^{(h)},p_{k|k - 1}^{(h)}( \cdot ,o,l)} )} \}_{h = 1}^{{H_{k|k - 1}}} $ , 其中$ L_{k|k-1}^{(h)} $ 和$ w_{k|k-1}^{(h)} $ 分别表示预测器标签集和权重,$ H_{k|k-1} $ 表示预测器假设数目.3.3.4 更新器
更新器输出结果为
$ \left\{ {{w_{k,i}}(l),{{\boldsymbol{\theta }}_{k,i}}(o,l)} \right\}_{i = 1}^{J_k^{(l)}},l \in { {L}_k} $ . 具体计算步骤如下:步骤 7. 高阶似然函数
$ g({\boldsymbol{Z}}_k|{\boldsymbol{\theta}}_{k|k-1,i}(o,l)) $ 和检测概率$ p_{D,k}({\boldsymbol{\theta}}_{k| k-1,i}(o,l)) $ 分别由式(7)、(17)获取.步骤 8. 计算
$ \eta _Z(l) $ :$$ \begin{split} \eta _Z(l) =\;& \sum\limits_{o \in O} {{P_{k,i}}(o)\sum\limits_{i = 1}^{J_{k|k - 1}^{(l)}} {\left( {w_{k|k - 1,i}(l)} \right.} }\times \\ &\left. {{p_{D,k}}({\boldsymbol{\theta }}_{k|k - 1,i}(l))g({{\boldsymbol{Z}}_k}|{\boldsymbol{\theta }}_{k|k - 1,i}(o,l))} \right) \end{split} $$ 步骤 9. 更新权重
$ w_{k,i}(l) $ :$$\begin{aligned}&w_{k ,i}(l)= \\ &\quad\quad\frac{{w_{k|k - 1,i}(l){p_{D,k}}({\boldsymbol{\theta }}_{k|k - 1,i}(l))g({{\boldsymbol{Z}}_k}|{\boldsymbol{\theta }}_{k|k - 1,i}(o,l))}}{{\eta _Z(l)}}\end{aligned} $$ 最终可得
$ k $ 时刻多目标后验概率密度参数集为$ \{ {( {L_k^{(h)},w_k^{(h)},p_{Z,k}^{(h)}( \cdot ,o,l|{{\boldsymbol{Z}}_k})} )} \}_{h = 1}^{{H_k}} $ , 其中$ H_k $ 为更新器假设数目. 注意: 更新器中后验假设的标签集合依据预测假设的标签集合给定.3.3.5 剪枝、合并及状态估计
算法的剪枝、合并以及声目标状态估计与常规
$ \delta $ -GLMB算法[21-23]一致. 具体处理过程如下:步骤 10. 剪枝与合并: 修剪粒子权重低于指定剪枝阈值的分量, 并将粒子总数控制在最大限度内.
步骤 11. 目标数目、状态估计: 确定势分布中概率最高的势(对应于目标数目), 然后在该势下寻找权值比设定阈值高的分量进行目标状态提取.
依据以上Cum-AMMS-GLMB算法步骤即可实现基于AVS的多机动声目标跟踪. 由于本文重点在于高阶似然函数的设计和添加表征不同运动模型索引的增广参数的M
$ \delta $ -GLMB算法推导, 因此预测、更新器中假设截断方法参考文献[21]实现, 本部分不多加赘述.3.4 算法计算复杂度
根据上述实现步骤, Cum-AMMS-GLMB算法的计算复杂度主要增加在以下部分:
$ 1) $ 四通道AVS的高阶累积量预处理, 计算复杂度为$ {\rm{O}}(N^4) $ , 其中$ N=4 $ 为AVS通道数目.$ 2) $ 预测部分增加了对运动模型概率的预测、更新. 给定预测器中$ H $ 个假设下的标签集合$ \{L^{(h)}\}_{h=1}^{H} $ 、粒子数目$ J $ 以及运动模型数目$ O $ , 则运动模型概率预测、更新的总计算复杂度为$ {\rm{O}}(O\cdot \sum_{h=1}^{H}{ ( | L^{(h)} | J )}) $ , 其中$ | L^{(h)} | $ 表示第$ h $ 个假设下的目标数目.$ 3) $ 更新部分增加了高阶似然的求解. 高阶似然通过$ N^2 \times N^2 $ 维累积量矩阵$ {\boldsymbol{R}}_{Z,k}^{(\rm{cum})} $ 的特征分解获取, 计算复杂度为$ {\rm{O}}((N^2)^3) $ .因此, Cum-AMMS-GLMB算法的主要计算复杂度为
$ {\rm{O}}(N^4+O\cdot \sum_{h=1}^{H}{ ( | L^{(h)} | J )}+N^6) $ .Cum-AMMS-GLMB算法滤波主体部分无需引入IMM-GLMB算法中的航迹−量测关联映射, 但由于执行高阶累积量预处理, 算法计算复杂度依然高于IMM-GLMB算法. 然而, 后续仿真验证了该算法可增加AVS的可跟踪声目标数目、提高跟踪性能. 算法的跟踪优势可一定程度上补偿略高的计算复杂度. 为降低计算复杂度, 可在预处理后去除虚拟AVS的冗余通道, 从而降低矩阵特征分解的计算复杂度. 除此之外, 还可借鉴GLMB快速算法[45-46], 通过联合预测和更新步骤以降低计算复杂度.
4. 声目标跟踪性能评价
本部分引入后验克拉美罗下界PCRLB[47]和最优子模型分配距离(Optimal sub-pattern assignment, OSPA)[48]两种指标, 从不同角度评估算法性能.
4.1 后验克拉美罗界
PCRLB能够提供随机变量估计所能达到的最低误差界, 可通过比较不同条件下算法(Root mean square error, RMSE)与PCRLB的逼近程度以评估算法跟踪性能. 由于本文引入了IMM思想, 已有PCRLB通解方法已不再适用, 因此本部分推导求解多运动模型交互的AVS目标跟踪的PCRLB.
4.1.1 基于多运动模型交互的PCRLB推导
式(1)的过程噪声协方差矩阵为奇异矩阵, 因此需将状态向量
$ {\boldsymbol{\theta}}_k $ 分割为目标DOA向量$ {\boldsymbol{\theta}}_k^{(1)}= \{\alpha_k,\beta_k\} $ , 速度向量$ {\boldsymbol{\theta}}_k^{(2)}=\{{\dot{\alpha}}_k,{\dot{\beta}}_k\} $ 和加速度向量$ {\boldsymbol{\theta}}_k^{(3)}=\{{\ddot{\alpha}}_k,{\ddot{\beta}}_k\} $ 三部分. 以CV、CA和CT三种运动模型为例, 实际上算法将三种运动模型通过一定的概率权重进行混合以实现声目标状态拟合. 记$ \{p^{(o)}\}_{o=1}^3 $ 分别为算法对以上三种运动模型更新的模型概率, 则$ k $ 时刻目标状态可表示为:$$ \begin{equation} { \begin{aligned} &{\boldsymbol{\theta}}_k^{(1)}={\boldsymbol{D}}^{(1)}{\boldsymbol{\theta}}_{k-1}^{(1)} +{\boldsymbol{D}}^{(2)}{\boldsymbol{\theta}}_{k-1}^{(2)}+{\boldsymbol{D}}^{(3)}{\boldsymbol{\theta}}_{k-1}^{(3)} + {\boldsymbol{D}}^{(4)}{\boldsymbol{w}}_{k-1}\\ &{\boldsymbol{\theta}}_k^{(2)}={\boldsymbol{R}}^{(1)}{\boldsymbol{\theta}}_{k-1}^{(1)} +{\boldsymbol{R}}^{(2)}{\boldsymbol{\theta}}_{k-1}^{(2)}+{\boldsymbol{R}}^{(3)}{\boldsymbol{\theta}}_{k-1}^{(3)} +{\boldsymbol{R}}^{(4)}{\boldsymbol{\theta}}_k^{(1)}\\ &{\boldsymbol{\theta}}_k^{(3)}={\boldsymbol{U}}^{(1)}{\boldsymbol{\theta}}_{k-1}^{(1)} +{\boldsymbol{U}}^{(2)}{\boldsymbol{\theta}}_{k-1}^{(2)}+{\boldsymbol{U}}^{(3)}{\boldsymbol{\theta}}_{k-1}^{(3)} +{\boldsymbol{U}}^{(4)}{\boldsymbol{\theta}}_k^{(1)} \end{aligned} } \end{equation} $$ (18) 其中
$$ \begin{equation} \nonumber \begin{split} &{\boldsymbol{D}}^{(1)}={\boldsymbol{I}}_2, {\boldsymbol{D}}^{(2)}=\left((p^{(1)}+p^{(3)})T{\boldsymbol{I}}_2+p^{(2)}{\boldsymbol{G}}^{(1)}\right)\\ &{\boldsymbol{D}}^{(3)}=\frac{p^{(3)}T^2}{2}{\boldsymbol{I}}_2, {\boldsymbol{D}}^{(4)}=\Big(\frac{(p^{(1)}+p^{(2)})T^2}{2}+\frac{p^{(3)}T^3}{6}\Big){\boldsymbol{I}}_2\\ &{\boldsymbol{R}}^{(1)}=-\left(\frac{2(p^{(1)}+p^{(2)})+3p^{(3)}}{T}\right){\boldsymbol{I}}_2 \end{split} \end{equation} $$ $$ \begin{equation} \nonumber \begin{split} &{\boldsymbol{R}}^{(2)}=-\left((p^{(1)}+2p^{(3)}){\boldsymbol{I}}_2+p^{(2)}(\frac{2}{T}{\boldsymbol{G}}^{(1)}-{\boldsymbol{G}}^{(2)})\right)\\ &{\boldsymbol{R}}^{(3)}=-\frac{p^{(3)}T}{2}{\boldsymbol{I}}_2, {\boldsymbol{R}}^{(4)}=\left(\frac{2(p^{(1)}+p^{(2)})+3p^{(3)}}{T}\right){\boldsymbol{I}}_2\\ &{\boldsymbol{U}}^{(1)}=-\frac{6p^{(3)}}{T^2}{\boldsymbol{I}}_2,{\boldsymbol{U}}^{(2)}=-\frac{6p^{(3)}}{T}{\boldsymbol{I}}_2\\ &{\boldsymbol{U}}^{(3)}=-2p^{(3)}{\boldsymbol{I}}_2,{\boldsymbol{U}}^{(4)}=\frac{6p^{(3)}}{T^2}{\boldsymbol{I}}_2\\ &{\boldsymbol{G}}^{(1)}=\left[\begin{matrix}\dfrac{\sin{\omega T}}{\omega}&\dfrac{\cos{\omega T}-1}{\omega}\\\dfrac{1-\cos{\omega T}}{\omega}&\dfrac{\sin{\omega T}}{\omega}\\\end{matrix}\right]\\ &{\boldsymbol{G}}^{(2)}=\left[\begin{matrix}\cos{\omega T}&-\sin{\omega T}\\\sin{\omega T}&\cos{\omega T}\\\end{matrix}\right] \end{split} \end{equation} $$ 记
$ \{\boldsymbol{\theta}_{1:k-2}^{(1)},{\boldsymbol{\theta}}_{k-1}^{(1)}, {\boldsymbol{\theta}}_{k-1}^{(2)},{\boldsymbol{\theta}}_{k-1}^{(3)},{\boldsymbol{Z}}_{k-1}\} $ 的联合概率密度为:$ P_{k-1}=P({\boldsymbol{\theta}}_{1:k-2}^{(1)},{\boldsymbol{\theta}}_{k-1}^{(1)},{\boldsymbol{\theta}}_{k-1}^{(2)},{\boldsymbol{\theta}}_{k-1}^{(3)},{\boldsymbol{Z}}_{k-1}) $ , 其信息矩阵为:$$ \begin{equation} \begin{aligned}[b] &{\tilde{{\boldsymbol{ J}}}}({\boldsymbol{\theta }}_{1:k - 2}^{(1)},{\boldsymbol{\theta }}_{k-1}^{(1)},{\boldsymbol{\theta }}_{k-1}^{(2)},{\boldsymbol{\theta }}_{k-1}^{(3)}) = \left[ {\begin{array}{*{20}{c}} {{{{\tilde{{\boldsymbol{ J}}}}}^{(11)}}}&{{{{\tilde{{\boldsymbol{ J}}}}}^{(12)}}}\\ {{{{\tilde{{\boldsymbol{ J}}}}}^{(21)}}}&{{{{\tilde{{\boldsymbol{ J}}}}}^{(22)}}} \end{array}} \right]=\\ &\qquad \left[ {\begin{array}{*{20}{c|ccc}} {{\tilde{{\boldsymbol{ J}}}}^{11}}&{{\tilde{{\boldsymbol{ J}}}}^{12}}&{{\tilde{{\boldsymbol{ J}}}}^{13}}&{{\tilde{{\boldsymbol{ J}}}}^{14}}\\ \hline\\[-10pt] {{{({\tilde{{\boldsymbol{ J}}}}^{12})}^{\mathrm T}}}^{\;}&{{\tilde{{\boldsymbol{ J}}}}^{22}}&{{\tilde{{\boldsymbol{ J}}}}^{23}}&{{\tilde{{\boldsymbol{ J}}}}^{24}}\\ {{{({\tilde{{\boldsymbol{ J}}}}^{13})}^{\mathrm T}}}&{{{({\tilde{{\boldsymbol{ J}}}}^{23})}^{\mathrm T}}}&{{\tilde{{\boldsymbol{ J}}}}^{33}}&{{\tilde{{\boldsymbol{ J}}}}^{34}}\\ {{{({\tilde{{\boldsymbol{ J}}}}^{14})}^{\mathrm T}}}&{{{({\tilde{{\boldsymbol{ J}}}}^{24})}^{\mathrm T}}}&{{{({\tilde{{\boldsymbol{ J}}}}^{34})}^{\mathrm T}}}&{{\tilde{{\boldsymbol{ J}}}}^{44}} \end{array}} \right]\\[-40pt] \end{aligned} \end{equation} $$ (19) 则
$ k -1$ 时刻目标状态$ {\boldsymbol{\theta}}_{k-1} $ 的信息矩阵$ {\boldsymbol{ J}}_{k-1} $ 即为$ \tilde{{\boldsymbol{ J}}}({\boldsymbol{\theta}}_{1:k-2}^{(1)},{\boldsymbol{\theta}}_{k-1}^{(1)}, {\boldsymbol{\theta}}_{k-1}^{(2)},{\boldsymbol{\theta}}_{k-1}^{(3)})^{-1} $ 右下分块矩阵的逆:$$ \begin{split} {\boldsymbol{ J}}_{k-1}=\;&{\tilde{\boldsymbol{J}}}^{({21})}-{\tilde{\boldsymbol{J}}}^{({22})} ({\tilde{\boldsymbol{J}}}^{({11})})^{-1}{\tilde{\boldsymbol{J}}}^{({12})}=\\ &\left[\begin{matrix}{\boldsymbol{ J}}^{11}&{\boldsymbol{ J}}^{12}&{\boldsymbol{ J}}^{13}\\ {\boldsymbol{ J}}^{21}&{\boldsymbol{ J}}^{22}&{\boldsymbol{ J}}^{23}\\ {\boldsymbol{ J}}^{31}&{\boldsymbol{ J}}^{32}&{\boldsymbol{ J}}^{33}\\ \end{matrix}\right] \end{split} $$ (20) 根据式(18)可得:
$$ \begin{equation} \begin{split} \left[ {\begin{array}{*{20}{c}} {{\boldsymbol{\theta }}_{1:k - 2}^{(1)}}\\ {{\boldsymbol{\theta }}_{k-1}^{(1)}}\\ {{\boldsymbol{\theta }}_{k }^{(1)}}\\ {{\boldsymbol{\theta }}_{k }^{(2)}}\\ {{\boldsymbol{\theta }}_{k }^{(3)}} \end{array}} \right] = \;&\left[ {\begin{array}{*{20}{c|cccc}} {{{\boldsymbol{I}}_2}}&{\bf{0}}&{\bf{0}}&{\bf{0}}&{\bf{0}}\\ \hline {\bf{0}}&{{{\boldsymbol{I}}_2}}&{\bf{0}}&{\bf{0}}&{\bf{0}}\\ {\bf{0}}&{\bf{0}}&{\bf{0}}&{\bf{0}}&{{{\boldsymbol{I}}_2}}\\ {\bf{0}}&{{\boldsymbol{R}}^{(1)}}&{{\boldsymbol{R}}^{(2)}}&{{\boldsymbol{R}}^{(3)}}&{{\boldsymbol{R}}^{(4)}}\\ {\bf{0}}&{{\boldsymbol{U}}^{(1)}}&{{\boldsymbol{U}}^{(2)}}&{{\boldsymbol{U}}^{(3)}}&{{\boldsymbol{U}}^{(4)}} \end{array}} \right]\times\\ &\left[ {\begin{array}{*{20}{c}} {{\boldsymbol{\theta }}_{1:k - 2}^{(1)}}\\ {{\boldsymbol{\theta }}_{k-1}^{(1)}}\\ {{\boldsymbol{\theta }}_{k-1}^{(2)}}\\ {{\boldsymbol{\theta }}_{k-1}^{(3)}}\\ {{\boldsymbol{\theta }}_{k }^{(1)}} \end{array}} \right] \\[-50pt]\end{split} \end{equation} $$ (21) 对于函数
$ y=Mx $ ,$ x $ 与$ y $ 的信息矩阵关系为$ J(y)= M^{-\mathrm T}J(x)M $ [47]. 记式(21)等式右侧右下分块矩阵为$ \tilde{{\boldsymbol{ M}}} $ , 则有:$$ \begin{equation} \begin{split} &{\boldsymbol{ J}}({\boldsymbol{\theta}}_{k-1}^{(1)},{\boldsymbol{\theta}}_{k}^{(1)},{\boldsymbol{\theta}}_{k}^{(2)}, {\boldsymbol{\theta}}_{k}^{(3)})=\\ &\qquad{\tilde{\boldsymbol{M}}}^{-\mathrm T}{\boldsymbol{ J}}({\boldsymbol{\theta}}_{k-1}^{(1)},{\boldsymbol{\theta}}_{k-1}^{(2)},{\boldsymbol{\theta}}_{k-1}^{(3)}, {\boldsymbol{\theta}}_{k}^{(1)}){\tilde{\boldsymbol{M}}} \end{split} \end{equation} $$ (22) 为求解
$ {\boldsymbol{\theta}}_{k} $ 的信息矩阵$ {{\boldsymbol{ J}}}_k $ , 首先给出$ \{\boldsymbol{\theta}_{1:{k-2}}^{(1)}, $ $ {\boldsymbol{\theta}}_{k-1}^{(1)},{\boldsymbol{\theta}}_{k-1}^{(2)},{\boldsymbol{\theta}}_{k-1}^{(3)}, {\boldsymbol{\theta}}_{k}^{(1)},{\boldsymbol{Z}}_{k}\} $ 的联合概率密度函数为:$$ \begin{equation} \nonumber P_{k}=P_{k-1} \cdot p({\boldsymbol{\theta}}_{k}^{(1)}|{\boldsymbol{\theta}}_{k-1})p({\boldsymbol{Z}}_{k}|{\boldsymbol{\theta}}_{k}^{(1)})=P_{k-1} \cdot\tilde P_{k} \end{equation} $$ 其中
$ \tilde P_{k} = p({\boldsymbol{\theta}}_{k}^{(1)}|{\boldsymbol{\theta}}_{k-1}) p({\boldsymbol{Z}}_{k}|{\boldsymbol{\theta}}_{k}^{(1)}) $ . 结合式(19), 可得$ \{\boldsymbol{\theta}_{1:k-2}^{(1)}, {\boldsymbol{\theta}}_{k-1}^{(1)}, {\boldsymbol{\theta}}_{k-1}^{(2)},{\boldsymbol{\theta}}_{k-1}^{(3)},{\boldsymbol{\theta}}_{k}^{(1)}\} $ 的信息矩阵为:$$ \begin{equation} \begin{split} &{\boldsymbol{J}}({\boldsymbol{\theta }}_{1:k-2}^{(1)},{\boldsymbol{\theta }}_{k-1}^{(1)},{\boldsymbol{\theta}}_{k-1}^{(2)},{\boldsymbol{\theta }}_{k-1}^{(3)},{\boldsymbol{\theta }}_{k}^{(1)})=\\ &{ \left[ {\begin{array}{*{20}{c|cccc}} {{\tilde{{\boldsymbol{ J}}}}^{11}} & {{\tilde{{\boldsymbol{ J}}}}^{12}} & {{\tilde{{\boldsymbol{ J}}}}^{13}} & {{\tilde{{\boldsymbol{ J}}}}^{14}} & {\bf{0}}\\ \hline\\[-10pt] {{{({\tilde{{\boldsymbol{ J}}}}^{12})}^{\mathrm T}}} & {{\tilde{{\boldsymbol{ J}}}}^{22} + {\boldsymbol{L}}_k^{11}} & {{\tilde{{\boldsymbol{ J}}}}^{23} + {\boldsymbol{L}}_k^{12}} & {{\tilde{{\boldsymbol{ J}}}}^{24} + {\boldsymbol{L}}_k^{13}} & {{\boldsymbol{L}}_k^{14}}\\ {{{({\tilde{{\boldsymbol{ J}}}}^{13})}^{\mathrm T}}} & {{{({\tilde{{\boldsymbol{ J}}}}^{23} + {\boldsymbol{L}}_k^{12})}^{\mathrm T}}} & {{\tilde{{\boldsymbol{ J}}}}^{33} + {\boldsymbol{L}}_k^{22}} & {{\tilde{{\boldsymbol{ J}}}}^{34} + {\boldsymbol{L}}_k^{23}} & {{\boldsymbol{L}}_k^{24}}\\ {{{({\tilde{{\boldsymbol{ J}}}}^{14})}^{\mathrm T}}} & {{{({\tilde{{\boldsymbol{ J}}}}^{24} + {\boldsymbol{L}}_k^{13})}^{\mathrm T}}} & {{{({\tilde{{\boldsymbol{ J}}}}^{34} + {\boldsymbol{L}}_k^{23})}^{\mathrm T}}} & {{\tilde{{\boldsymbol{ J}}}}^{44} + {\boldsymbol{L}}_k^{33}} & {{\boldsymbol{L}}_k^{34}}\\ {\bf{0}} & {{{({\boldsymbol{L}}_k^{14})}^{\mathrm T}}} & {{{({\boldsymbol{L}}_k^{24})}^{\mathrm T}}} & {{{({\boldsymbol{L}}_k^{34})}^{\mathrm T}}} & {{\boldsymbol{L}}_k^{44}} \end{array}} \right] } \end{split} \end{equation} $$ (23) 其中
$$ \begin{equation} \nonumber { \begin{array}{l} {\boldsymbol{L}}_k^{11} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(1)}}^{{\boldsymbol{\theta }}_{k-1}^{(1)}}\ln {\tilde P_{k }}} \Big\}, {\boldsymbol{L}}_k^{12} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(1)}}^{{\boldsymbol{\theta }}_{k-1}^{(2)}}\ln {\tilde P_{k }}} \Big\}\\ {\boldsymbol{L}}_k^{13} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(1)}}^{{\boldsymbol{\theta }}_{k-1}^{(3)}}\ln {\tilde P_{k }}} \Big\}, {\boldsymbol{L}}_k^{14} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(1)}}^{{\boldsymbol{\theta }}_{k}^{(1)}}\ln {\tilde P_{k }}} \Big\}\\ {\boldsymbol{L}}_k^{22} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(2)}}^{{\boldsymbol{\theta }}_{k-1}^{(2)}}\ln {\tilde P_{k }}} \Big\}, {\boldsymbol{L}}_k^{23} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(2)}}^{{\boldsymbol{\theta }}_{k-1}^{(3)}}\ln {\tilde P_{k }}} \Big\}\\ {\boldsymbol{L}}_k^{24} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(2)}}^{{\boldsymbol{\theta }}_{k}^{(1)}}\ln {\tilde P_{k }}} \Big\}, {\boldsymbol{L}}_k^{33} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(3)}}^{{\boldsymbol{\theta }}_{k-1}^{(3)}}\ln {\tilde P_{k }}} \Big\}\\ {\boldsymbol{L}}_k^{34} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k-1}^{(3)}}^{{\boldsymbol{\theta }}_{k}^{(1)}}\ln {\tilde P_{k }}} \Big\}, {\boldsymbol{L}}_k^{44} = \mathrm E\Big\{ { - \Delta _{{\boldsymbol{\theta }}_{k}^{(1)}}^{{\boldsymbol{\theta }}_{k}^{(1)}}\ln {\tilde P_{k }}} \Big\} \end{array} } \end{equation} $$ 因此信息矩阵
$ {\tilde{\boldsymbol{S}}_{k}}={\boldsymbol{ J}}\{\boldsymbol{\theta}_{k-1}^{(1)},{\boldsymbol{\theta}}_{k-1}^{(2)}, {\boldsymbol{\theta}}_{k-1}^{(3)},{\boldsymbol{\theta}}_{k}^{(1)}\} $ 即为式(23)逆的右下分块矩阵的逆:$$ \begin{equation} \nonumber { {\tilde{\boldsymbol{S}}_{k}} = \left[ {\begin{array}{*{20}{c}} {{\boldsymbol{J}}^{11} + {\boldsymbol{L}}_k^{11}} & {{\boldsymbol{J}}^{12} + {\boldsymbol{L}}_k^{12}} & {{\boldsymbol{J}}^{13} + {\boldsymbol{L}}_k^{13}} & {{\boldsymbol{L}}_k^{14}}\\ {{{({\boldsymbol{J}}^{12} + {\boldsymbol{L}}_k^{12})}^{\mathrm T}}} & {{\boldsymbol{J}}^{22} + {\boldsymbol{L}}_k^{22}} & {{\boldsymbol{J}}^{23} + {\boldsymbol{L}}_k^{23}} & {{\boldsymbol{L}}_k^{24}}\\ {{{({\boldsymbol{J}}^{13} + {\boldsymbol{L}}_k^{13})}^{\mathrm T}}} & {{{({\boldsymbol{J}}^{23} + {\boldsymbol{L}}_k^{23})}^{\mathrm T}}} & {{\boldsymbol{J}}^{33} + {\boldsymbol{L}}_k^{33}} & {{\boldsymbol{L}}_k^{34}}\\ {{{({\boldsymbol{L}}_k^{14})}^{\mathrm T}}} & {{{({\boldsymbol{L}}_k^{24})}^{\mathrm T}}} & {{{({\boldsymbol{L}}_k^{34})}^{\mathrm T}}} & {{\boldsymbol{L}}_k^{44}} \end{array}} \right] } \end{equation} $$ 依据式(22)可得,
$ {{\boldsymbol{S}}_{k }}={\boldsymbol{ J}}\{\boldsymbol{\theta}_{k-1}^{(1)},{\boldsymbol{\theta}}_{k}^{(1)}, {\boldsymbol{\theta}}_{k}^{(2)},{\boldsymbol{\theta}}_{k}^{(3)}\} $ 为:$$ \begin{equation} \begin{split} {{\boldsymbol{S}}_{k }} =\;& {\tilde{{\boldsymbol{ M}}}}^{ - \mathrm T}{{\tilde{{\boldsymbol{ S}}}}_{k}}{{\tilde{{\boldsymbol{ M}}}}} =\left[ {\begin{array}{*{20}{c}} {{\boldsymbol{S}}_{k }^{(11)}}&{{\boldsymbol{S}}_{k }^{(12)}}\\ {{\boldsymbol{S}}_{k }^{(21)}}&{{\boldsymbol{S}}_{k }^{(22)}} \end{array}} \right] =\\ & \left[ {\begin{array}{*{20}{c|ccc}} {{\boldsymbol{s}}_{k}^{11}}&{{\boldsymbol{s}}_{k}^{12}}&{{\boldsymbol{s}}_{k}^{13}}&{{\boldsymbol{s}}_{k}^{14}}\\ \hline\\[-10pt] {{\boldsymbol{s}}_{k}^{21}}&{{\boldsymbol{s}}_{k}^{22}}&{{\boldsymbol{s}}_{k}^{23}}&{{\boldsymbol{s}}_{k}^{24}}\\ {{\boldsymbol{s}}_{k}^{31}}&{{\boldsymbol{s}}_{k}^{32}}&{{\boldsymbol{s}}_{k}^{33}}&{{\boldsymbol{s}}_{k}^{34}}\\ {{\boldsymbol{s}}_{k}^{41}}&{{\boldsymbol{s}}_{k}^{43}}&{{\boldsymbol{s}}_{k}^{43}}&{{\boldsymbol{s}}_{k}^{44}} \end{array}} \right] \end{split} \end{equation} $$ (24) 因此
$ {\boldsymbol{S}}_{k}^{-1} $ 的右下分块矩阵的逆即为$ {\boldsymbol{\theta}}_{k} $ 的信息矩阵$ {{\boldsymbol{ J}}}_k $ :$$ \begin{equation} {{\boldsymbol{J}}_{k}} = {\boldsymbol{S}}_{k}^{(22)} - {\boldsymbol{S}}_{k}^{(21)}{( {{\boldsymbol{S}}_{k}^{(11)}} )^{ - 1}}{\boldsymbol{S}}_{k}^{(12)} \end{equation} $$ (25) 根据PCRLB与信息矩阵的关系, 最终可得
$ k $ 时刻目标状态分量估计的PCRLB为:$ {{\rm{PCRLB}}}_{k} ={{\boldsymbol{J}}_{k}^{-1}} $ .4.1.2 AVS的目标跟踪PCRLB
根据式(18)可得声目标状态分量的转移概率为:
$$ \begin{equation} p({\boldsymbol{\theta }}_{k}^{(1)}|{{\boldsymbol{\theta }}_{k-1}}) = \frac{{\exp \left\{ { - \frac{1}{2}{\boldsymbol{\theta }}_{\rm{new}}^{\mathrm T}{\boldsymbol{Q}}_w^{ - 1}{{\boldsymbol{\theta }}_{\rm{new}}}} \right\}}}{{2\pi \sqrt {\det ({{\boldsymbol{Q}}_w})} }} \end{equation} $$ (26) 其中
${\boldsymbol{\theta }}_{\rm{new}}={\boldsymbol{\theta}}_{k}^{(1)} - {\boldsymbol{D}}^{(1)}{\boldsymbol{\theta}}_{k-1}^{(1)}-{\boldsymbol{D}}^{(2)}{\boldsymbol{\theta}}_{k-1}^{(2)} - {\boldsymbol{D}}^{(3)}{\boldsymbol{\theta}}_{k-1}^{(3)}$ ,$ {\boldsymbol{Q}}_w={\boldsymbol{D}}^{(4)}{\boldsymbol{Q}}_{k-1}({\boldsymbol{D}}^{(4)})^{\mathrm T} $ . 因此式(23)中$ \{\boldsymbol{L}_k^{ij}\}_{i,j=1}^4 $ 为:$$ \begin{equation} \begin{split} &{\boldsymbol{L}}_k^{11} = {\boldsymbol{Q}}_w^{ - 1},{\boldsymbol{L}}_k^{12} = {({{\boldsymbol{D}}^{(2)}})^{\mathrm T}}{\boldsymbol{Q}}_w^{ - 1}\\ &{\boldsymbol{L}}_k^{13} = {({{\boldsymbol{D}}^{(3)}})^{\mathrm T}}{\boldsymbol{Q}}_w^{ - 1},{\boldsymbol{L}}_k^{14} = - {\boldsymbol{Q}}_w^{ - 1}\\ &{\boldsymbol{L}}_k^{22} = {({{\boldsymbol{D}}^{(2)}})^{\mathrm T}}{\boldsymbol{Q}}_w^{ - 1}{{\boldsymbol{D}}^{(2)}}\\ &{\boldsymbol{L}}_k^{23} = {({{\boldsymbol{D}}^{(3)}})^{\mathrm T}}{\boldsymbol{Q}}_w^{ - 1}{{\boldsymbol{D}}^{(2)}},{\boldsymbol{L}}_k^{24} = - {\boldsymbol{Q}}_w^{ - 1}{{\boldsymbol{D}}^{(2)}}\\ &{\boldsymbol{L}}_k^{33} = {({{\boldsymbol{D}}^{(3)}})^{\mathrm T}}{\boldsymbol{Q}}_w^{ - 1}{{\boldsymbol{D}}^{(3)}},{\boldsymbol{L}}_k^{34} = - {\boldsymbol{Q}}_w^{ - 1}{{\boldsymbol{D}}^{(3)}}\\ &{\boldsymbol{L}}_k^{44} = {\boldsymbol{Q}}_w^{ - 1} + {\mathrm E}\left\{ { - \Delta _{{\boldsymbol{\theta }}_{k}^{(1)}}^{{\boldsymbol{\theta }}_{k}^{(1)}}\ln p({{\boldsymbol{Z}}_{k}}|{\boldsymbol{\theta }}_{k}^{(1)})} \right\} \end{split} \end{equation} $$ (27) 其中
$ {\boldsymbol{L}}_k^{44} $ 右侧第二项即为AVS的DOA信息矩阵.为便于分析, 将
$ 4\times L $ 维量测$ {{\boldsymbol{ Z}}}_k $ 记为$ 4L\times 1 $ 维列向量$ {\tilde{\boldsymbol{ z}}_k} = {\mathrm {vec}}({{\boldsymbol{Z}}_k}) $ , 则后验概率密度函数为:$$ \begin{equation} p({\tilde{\boldsymbol{ z}}_k}|{\boldsymbol{\theta }}_{k }^{(1)}) = \frac{{ {\exp \left( { - {{({\tilde{\boldsymbol{ z}}_k} - {\boldsymbol{\mu }})}^{\mathrm H}}{{\boldsymbol{\varGamma }}^{ - 1}}({\tilde{\boldsymbol{ z}}_k} - {\boldsymbol{\mu }})} \right)} }}{{{\pi ^{4L}}\det ({\boldsymbol{\varGamma }})}} \end{equation} $$ (28) 其中
$ {\boldsymbol{\varGamma}}={\boldsymbol{\varGamma}}_0\otimes {\boldsymbol{I}}_L $ ,$ {\boldsymbol{\mu }} = \sum_{m = 1}^{M(k)} {{\boldsymbol{a}}({\boldsymbol{\theta }}_{k,m}^{(1)}) \otimes {{\boldsymbol{s}}_{k,m}^{\mathrm T}}} $ . 文献[11]推导了AVS的DOA信息矩阵为:$$ \begin{equation} \begin{split} {\boldsymbol{ J}}_{\rm{DOA}}=\;&{\rm{E}}\left\{ { - \Delta _{{\boldsymbol{\theta }}_{k}^{(1)}}^{{\boldsymbol{\theta }}_{k }^{(1)}}\ln {p}({{\boldsymbol{Z}}_{k}}|{\boldsymbol{\theta }}_{k}^{(1)})} \right\}=\\ &\left[ {\begin{array}{*{20}{c}} {2L \cdot {{\rm{SNR}}}\cdot\cos{\beta}_{k}}&{}\\ {}&{2L \cdot {{\rm{SNR}}}} \end{array}} \right] \end{split} \end{equation} $$ (29) 其中
$ {\beta}_{k} $ 为$ {\boldsymbol{\theta }}_{k}^{(1)} $ 的第二项, 表示$ k $ 时刻声目标的俯仰角.$ L $ 表示每一时刻的采样快拍数,${\rm{SNR}}$ 为信噪比.依据上述推导, 以下给出基于AVS的目标跟踪PCRLB的求解步骤:
$ 1) $ 初始化信息矩阵$ {\boldsymbol{ J}}_0 $ , 可依据先验信息设置, 无先验信息时可设为$ {\bf{0}} $ ;$ 2) $ 通过式(27)、(29)计算$ \{\boldsymbol{L}_k^{ij},i,j= 1,\cdots, 4\} $ , 并将其与$ {\boldsymbol{ J}}_{k-1} $ 结合获得$ {\tilde{\boldsymbol{S}}_{k}} $ , 最后通过式(24)计算$ {\boldsymbol{S}}_{k} $ ;$ 3) $ 根据式(25)求解$ {\boldsymbol{\theta}}_{k} $ 的信息矩阵$ {\boldsymbol{ J}}_{k} $ , 最终获得$ k $ 时刻的PCRLB为$ {\boldsymbol{ J}}_{k}^{-1} $ .后续仿真实验中, 将目标跟踪的时间平均PCRLB:
$\overline {{{PCRLB}}} =({\sum_{k=1}^{K}{{PCRLB}}_{k}})/{K}$ 作为均方根误差的参考下界, 通过比较不同SNR和过程噪声情况下的RMSE来评估算法跟踪性能.4.2 OSPA距离指标
OSPA距离能够直观反映算法对于每一时刻的目标个数和状态估计误差. 记
$ {\boldsymbol{\Theta}}=\{\boldsymbol{\theta}_1,\cdots,{\boldsymbol{\theta}}_m\} $ 和$ \widehat{\boldsymbol{\Theta}}=\{\hat{\boldsymbol{\theta}}_1,\cdots,\hat{\boldsymbol{\theta}}_n\} $ 分别表示真实和估计的多目标状态集, 则OSPA距离定义如下:$$ \begin{equation} \begin{split} &{D_{p,c}}({\boldsymbol{\Theta }},\widehat {\boldsymbol{\Theta }}) = \\ &\left\{ \begin{aligned} &{\left[ {\frac{1}{n}\left( {\mathop {\min }\limits_{\pi \in {\Pi _n}} \sum\limits_{i = 1}^m {{{\left( {{d_c}({{\boldsymbol{\theta }}_i},{{\hat {\boldsymbol{\theta }}}_{\pi (i)}})} \right)}^p}} + (n - m) \cdot {c^p}} \right)} \right]^{\frac{1}{p}}},\\ &\qquad\qquad\qquad\qquad m \le n\\ &{D_{p,c}}({\boldsymbol{\Theta }},\widehat {\boldsymbol{\Theta }}),\qquad\; m > n \end{aligned} \right. \end{split} \end{equation} $$ (30) 其中
$ {\Pi _n} $ 表示从集合$ \{1,2,\cdots\} $ 中取出$ m $ 个元素进行排列的集合.$ c>0 $ 表示最小截止距离, 表示$ {\boldsymbol{\theta}} $ 和$ \hat{\boldsymbol{\theta}} $ 之间的距离在$ c $ 处截断.$ p $ 为阶数, 用于惩罚多目标状态的估计偏差. 在后续仿真实验中取$ c=5 $ ,$ p=1 $ , 其中OPSA距离计算过程中通过匈牙利算法[49]实现最优分配问题的求解.5. 实验验证
5.1 多声源跟踪仿真
单个AVS作为接收传感器, 并放置5个运动声目标于整个实验场景中. 设置SNR为10 dB, 加速度扰动噪声标准差
$ \sigma_w $ 为$ 0.5^\circ/{\rm s}^2 $ , 采样率为10 kHz, 共采集50 s数据. 表1、表2分别给出了各时间段所对应的运动模型和声目标的运动状态.表 1 各时间段对应的运动模型Table 1 Movement model corresponding to each time period时间 1~15 s 16~30 s 31~40 s 41~50 s 运动模型 CV模型 CA模型
$0.1^\circ/{\rm s}^2$CT模型
$ \omega=-2\pi/180$CV模型 表 2 声目标的运动状态和幸存时间Table 2 Motion state and survival time of acoustic targets声目标 初始状态 存在时间 目标1 $\mathrm{DOA}:\ \{30^\circ,30^\circ\}$, 速度$:\ \{1^\circ/{\rm s},0.05^\circ/{\rm s}\}$ 1~50 s 目标2 $\mathrm{DOA}:\ \{300^\circ,80^\circ\}$, 速度$:\ \{-1^\circ/{\rm s},1^\circ/{\rm s}\}$ 1~50 s 目标3 $\mathrm{DOA}:\ \{200^\circ,100^\circ\}$, 速度$:\ \{0^\circ/{\rm s},0.1^\circ/{\rm s}\}$ 10~20 s 目标4 $\mathrm{DOA}:\ \{150^\circ,40^\circ\}$, 速度$:\ \{1^\circ/{\rm s},0^\circ/{\rm s}\}$ 1~50 s 目标5 $\mathrm{DOA}:\ \{90^\circ,80^\circ\}$, 速度$:\ \{1^\circ/{\rm s},0.5^\circ/{\rm s}\}$ 25~50 s 使用Cum-AMMS-GLMB算法进行声目标跟踪, 可得图3所示的声目标跟踪、目标数目估计和OSPA距离结果. 跟踪结果指出Cum-AMMS-GLMB算法在声目标数目大于AVS通道数目的情况下依然能够准确跟踪声目标. 这主要是由于该算法设计的高阶似然函数可对AVS进行虚拟扩展, 从而提高AVS的可跟踪声目标数目.
图 3 Cum-AMMS-GLMB算法的多声目标跟踪结果Fig. 3 Multiple acoustic target tracking results of Cum-AMMS-GLMB algorithm5.2 算法对比仿真
本部分对表3所示4种算法的跟踪性能进行比较. 算法均执行IMM处理, 不同之处在于算法滤波器和似然函数的使用. 值得注意的是, Cum-AMMS-GLMB算法基于M
$ \delta $ -GLMB滤波框架实现, 且其量测$ {\boldsymbol{Z}}_k $ 为所有声目标的叠加信号, 因此无需引入Cov-GLMB算法中的航迹−量测的关联映射历史. CBMeMBer算法和Cum-CBMeMBer算法并未设计对声目标状态进行航迹区分的数据关联方法, 但不影响目标跟踪的实现.表 3 各对比算法的似然函数、滤波器的区别Table 3 The difference between the likelihood function and the filter of each comparison algorithm5.2.1 场景一
仿真场景如下: AVS共采集50 s数据, 采样率为10 kHz, 信噪比为10 dB,
$ \sigma_w $ 设置为$ 0.5^\circ/{\rm s}^2 $ . 两个声目标初始DOA均为$ \{180^\circ,30^\circ\} $ , 目标1在第1 s入场, 初始速度为$ \{1^\circ/{\rm s},0.05^\circ/{\rm s}\} $ , 目标2在第20 s入场, 初始速度为$ \{-1^\circ/{\rm s},1^\circ/{\rm s}\} $ , 其中不同时刻目标对应的运动模型参考表1. 使用表3所示4种算法进行仿真, 可得图4(a) ~ 4(b)的声目标跟踪结果. 为避免实验结果的特殊性, 还执行50次蒙特卡洛实验获得图4(c) OSPA距离以评估算法性能.
图 4 双声目标情况下不同算法的估计结果Fig. 4 The estimation results of different algorithms in the case of two acoustic targets与对比算法相比, Cum-AMMS-GLMB算法能够更准确的跟踪声目标轨迹和估计目标数目, 并且具有更小的OSPA距离. CBMeMBer、以及Cum-CBMeMBer算法还需额外的数据关联才可区分声目标轨迹, 而基于GLMB框架的Cum-AMMS-GLMB算法和Cov-GLMB算法由于航迹标签参数的引入, 能够直接分辨目标轨迹. 通过图4可得, 第20 s目标新生时, 引入高阶似然的Cum-AMMS-GLMB算法和Cum-CBMeMBer算法能够准确估计目标数目, 区分两个声源, 而Cov-GLMB算法和CBMeMBer算法并未准确估计目标数目. 这主要是由于AVS本身孔径较小, 且新生、幸存目标之间的角度间隔较小, 因此使用传统方法对角度间隔较近目标进行角度分辨时无法准确区分目标, 从而导致算法性能下降.
5.2.2 场景二
仿真场景如下: 两个声目标的初始DOA分别为
$ \{180^\circ,30^\circ\} $ 和$ \{100^\circ,30^\circ\} $ , 初始速度分别为$\{1^\circ/{\rm s}, 0.05^\circ/{\rm s}\}$ 和$ \{-1^\circ/{\rm s},1^\circ/{\rm s}\} $ . 目标在不同时刻对应的运动模型参考表1.首先评估不同信噪比下表3所示各算法的跟踪性能.
$ \sigma_w $ 设为$ 0.5^\circ/{\rm s}^2 $ , SNR分别从0dB以5dB的间隔逐渐增加到20 dB, 进行蒙特卡洛实验获取图5 所示RMSE估计结果. 根据图5可得, 较低SNR情况下Cum-AMMS-GLMB算法的RMSE低于对比算法, 且随着SNR的增加, RMSE逐渐逼近PCRLB. 除此之外, 由于使用了式(7)所示的高阶似然, Cum-CBMeMBer算法的跟踪性能同样优于相同滤波器框架的CBMeMBer算法.
图 5 不同信噪比下PCRLB和各算法的RMSE估计结果Fig. 5 PCRLB and RMSE estimation results of each algorithm under different SNR对于不同过程噪声下各算法性能, SNR设置为5 dB, 加速度扰动噪声标准差
$ \sigma_w $ 分别从$ 0^\circ/{\rm s}^2 $ 以$ 0.5^\circ/{\rm s}^2 $ 的间隔逐渐增加到$ 4.5^\circ/{\rm s}^2 $ 进行蒙特卡洛实验. 根据图6所示RMSE估计结果可得, 随着$ \sigma_w $ 的增大, 各算法的RMSE和PCRLB均呈上升趋势. Cum-AMMS-GLMB算法的RMSE性能指标优于对比算法, 更逼近于PCRLB.
图 6 不同$\sigma_w$ 下PCRLB和各算法的RMSE估计结果Fig. 6 PCRLB and RMSE estimation results of each algorithm under different$\sigma_w$ 5.3 真实实验
为验证理论算法可行性, 使用实验室自制AVS于半消声室进行声目标跟踪实验. 实验场景如图7(a)所示, AVS与三个声目标约在同一水平面, 以保证各声目标的俯仰角近似为
$ 0^\circ $ . 三个声目标分别以图7(b)的轨迹运动, 静止声目标1的方位角约为$ 8.1^\circ $ , 运动声目标2、3中, 目标2的方位角从$ 60.95^\circ $ 逐渐增加到$ 113.96^\circ $ , 目标3的方位角从$ 240.95^\circ $ 逐渐增加到$ 299.05^\circ $ . 运动过程中目标1、2一直存在, 目标3分别于一直存在、突然出现和突然消失三种情况下进行实验. 实验室自制采集系统进行信号采集, 采样率为3 kHz, 使用Cum-AMMS-GLMB算法进行声目标跟踪, 每一秒进行一次跟踪, 可得图8~10所示的声目标跟踪结果.根据图8所示声目标一直存在情况下的跟踪结果可得, 静止目标1在
$ \{7.13^\circ, 0.82^\circ\} $ 附近波动; 运动目标2的方位角跟踪轨迹从$ 61.88^\circ $ 逐渐增加到$ 109.91^\circ $ , 俯仰角约在$ 1.02^\circ $ 附近波动; 目标3的方位角跟踪轨迹从$ 233.38^\circ $ 逐渐增加到$ 300.95^\circ $ , 俯仰角约在$ 0.66^\circ $ 附近波动. 图9、图10所示的声目标跟踪结果可得, 在声目标突然出现情况下, 目标3在第19 s突然出现, 其方位角跟踪轨迹从$ 239.01^\circ $ 逐渐增加到$ 294.06^\circ $ . 对于声目标突然消失的情况, 目标3的方位角跟踪轨迹从$ 227.52^\circ $ 逐渐增加到$ 267.21^\circ $ , 随后在第38 s消失.以上三种不同场景的跟踪结果中, 声目标2、3的跟踪结果存在波动. 这主要是由于真实DOA轨迹为测量声目标在起始、终止位置处的连线, 而在匀速移动声目标2、3的过程中会产生抖动, 从而导致运动目标的估计轨迹存在角度波动. 不考虑这一不可避免的因素和传声器自身的测量误差, 图8~10所示的声目标跟踪结果可得, Cum-AMMS-GLMB算法估计的目标跟踪轨迹与真实运动轨迹基本一致, 该算法是可行、有效的.
6. 结论
本文提出的基于AVS的Cum-AMMS-GLMB算法, 解决了AVS可跟踪目标数目少、跟踪性能较差的问题. 算法中高阶累积量预处理的引入和高阶似然函数的建立, 不仅有效抑制了高斯量测噪声, 还通过AVS的虚拟扩展, 提高了可跟踪声目标数目. Cum-AMMS-GLMB算法更适配于声目标运动模型多样的实际场景. 在M
$ \delta $ -GLMB滤波框架下, 算法通过多种运动模型的加权混合以逼近目标的真实运动状态, 获得了优于单一运动模型的跟踪效果. 该算法无需引入航迹−量测的关联映射, 能够降低滤波器的存储容量和计算成本, 从而抵消了预处理过程的部分计算负担. 除此之外, 算法通过归一化高阶空间谱实现对检测概率函数的拟合, 增强了高似然区域粒子, 抑制了杂波向可用粒子的扩散. 实验验证了算法的可行性和有效性, 并且较已有基于AVS的声目标跟踪算法的性能更优. 显然所提小尺寸、便携AVS的多目标DOA跟踪, 对未来机器人听觉跟踪系统的设计具有重要参考价值. -
图 1 不同信噪比下的归一化高阶似然函数示例
Fig. 1 Example of normalized higher-order spatial spectrum under different SNR
图 3 Cum-AMMS-GLMB算法的多声目标跟踪结果
Fig. 3 Multiple acoustic target tracking results of Cum-AMMS-GLMB algorithm
图 4 双声目标情况下不同算法的估计结果
Fig. 4 The estimation results of different algorithms in the case of two acoustic targets
图 5 不同信噪比下PCRLB和各算法的RMSE估计结果
Fig. 5 PCRLB and RMSE estimation results of each algorithm under different SNR
图 6 不同
$\sigma_w$ 下PCRLB和各算法的RMSE估计结果Fig. 6 PCRLB and RMSE estimation results of each algorithm under different
$\sigma_w$ 表 1 各时间段对应的运动模型
Table 1 Movement model corresponding to each time period
时间 1~15 s 16~30 s 31~40 s 41~50 s 运动模型 CV模型 CA模型
$0.1^\circ/{\rm s}^2$CT模型
$ \omega=-2\pi/180$CV模型 
下载: 导出CSV
表 2 声目标的运动状态和幸存时间
Table 2 Motion state and survival time of acoustic targets
声目标 初始状态 存在时间 目标1 $\mathrm{DOA}:\ \{30^\circ,30^\circ\}$, 速度$:\ \{1^\circ/{\rm s},0.05^\circ/{\rm s}\}$ 1~50 s 目标2 $\mathrm{DOA}:\ \{300^\circ,80^\circ\}$, 速度$:\ \{-1^\circ/{\rm s},1^\circ/{\rm s}\}$ 1~50 s 目标3 $\mathrm{DOA}:\ \{200^\circ,100^\circ\}$, 速度$:\ \{0^\circ/{\rm s},0.1^\circ/{\rm s}\}$ 10~20 s 目标4 $\mathrm{DOA}:\ \{150^\circ,40^\circ\}$, 速度$:\ \{1^\circ/{\rm s},0^\circ/{\rm s}\}$ 1~50 s 目标5 $\mathrm{DOA}:\ \{90^\circ,80^\circ\}$, 速度$:\ \{1^\circ/{\rm s},0.5^\circ/{\rm s}\}$ 25~50 s
下载: 导出CSV
-
[1] Leslie C B, Kendall J M, Jones J L. Hydrophone for measuring particle velocity. The Journal of the Acoustical Society of America, 1956, 28(4): 711-715 doi: 10.1121/1.1908455 [2] Nehorai A, Paldi E. Acoustic vector-sensor array processing. IEEE Transactions on Signal Processing, 1994, 42(9): 2481-2491 doi: 10.1109/78.317869 [3] Bereketli A, Guldogan M B, Kolcak T, Gudu T, Avsar A L. Experimental results for direction of arrival estimation with a single acoustic vector sensor in shallow water. Journal of Sensors, 2015, 2015: Article No. 401353 [4] Zhang W D, Guan L G. Zhang G J, Xue C Y. Zhang K R, Wang J P. Research of DOA estimation based on single MEMS vector hydrophone. Sensors, 2009, 9(9): 6823-6831 doi: 10.3390/s90906823 [5] Levin D, Habets E A P, Gannot S. Maximum likelihood estimation of direction of arrival using an acoustic vector-sensor. The Journal of the Acoustical Society of America, 2012, 131(2): 1240-1248 doi: 10.1121/1.3676699 [6] Hawkes M, Nehorai A. Acoustic vector-sensor beamforming and Capon direction estimation. IEEE Transactions on Signal Processing, 1998, 46(9): 2291-2304 doi: 10.1109/78.709509 [7] Wong K T, Zoltowski M D. Root-MUSIC-based azimuth-elevation angle-of-arrival estimation with uniformly spaced but arbitrarily oriented velocity hydrophones. IEEE Transactions on Signal Processing, 1999, 47(12): 3250-3260 doi: 10.1109/78.806070 [8] Miron S, Bihan N L, Mars J I. Quaternion-MUSIC for vector-sensor array processing. IEEE Transactions on Signal Processing, 2006, 54(4): 1218-1229 doi: 10.1109/TSP.2006.870630 [9] Tichavsky P, Wong K T, Zoltowski M D. Near-field/far-field azimuth and elevation angle estimation using a single vector hydrophone. IEEE Transactions on Signal Processing, 2001, 49(11): 2498-2510 doi: 10.1109/78.960397 [10] Wong K T, Zoltowski M D. Extended-aperture underwater acoustic multisource azimuth/elevation direction-finding using uniformly but sparsely spaced vector hydrophones. IEEE Journal of Oceanic Engineering, 1997, 22(4): 659-672 doi: 10.1109/48.650832 [11] Zhang J, Xu X Y, Chen Z F, Bao M, Zhang X P, Yang J H. High-resolution DOA estimation algorithm for a single acoustic vector sensor at low SNR. IEEE Transactions on Signal Processing, 2020, 68: 6142-6158 doi: 10.1109/TSP.2020.3021237 [12] Mahler R P S. Statistical Multisource-multitarget Information Fusion. Boston: Artech House, 2007. [13] Mahler R P S. Multitarget Bayes filtering via first-order multitarget moments. IEEE Transactions on Aerospace and Electronic Systems, 2003, 39(4): 1152-1178 doi: 10.1109/TAES.2003.1261119 [14] Vo B N, Ma W K. The Gaussian mixture probability hypothesis density filter. IEEE Transactions on Signal Processing, 2006, 54(11): 4091-4104 doi: 10.1109/TSP.2006.881190 [15] Vo B N, Singh S, Doucet A. Sequential Monte Carlo methods for multitarget filtering with random finite sets. IEEE Transactions on Aerospace and Electronic Systems, 2005, 41(4): 1224-1245 doi: 10.1109/TAES.2005.1561884 [16] Mahler R. PHD filters of higher order in target number. IEEE Transactions on Aerospace and Electronic Systems, 2007, 43(4): 1523-1543 doi: 10.1109/TAES.2007.4441756 [17] Vo B T, Vo B N, Cantoni A. The cardinality balanced multi-target multi-Bernoulli filter and its implementations. IEEE Transactions on Signal Processing, 2009, 57(2): 409-423 doi: 10.1109/TSP.2008.2007924 [18] Ristic B, Vo B T, Vo B N, Farina A. A tutorial on Bernoulli filters: Theory, implementation and applications. IEEE Transactions on Signal Processing, 2013, 61(13): 3406-3430 doi: 10.1109/TSP.2013.2257765 [19] Wong S, Vo B T, Papi F. Bernoulli forward-backward smoothing for track-before-detect. IEEE Signal Processing Letters, 2014, 21(6): 727-731 doi: 10.1109/LSP.2014.2310137 [20] Vo B T, See C M, Ma N, Ng W T. Multi-sensor joint detection and tracking with the Bernoulli filter. IEEE Transactions on Aerospace and Electronic Systems, 2012, 48(2): 1385-1402 doi: 10.1109/TAES.2012.6178069 [21] Vo B N, Vo B T, Phung D. Labeled random finite sets and the Bayes multi-target tracking filter. IEEE Transactions on Signal Processing, 2014, 62(24): 6554-6567 doi: 10.1109/TSP.2014.2364014 [22] Vo B N, Vo B T. Labeled random finite sets and multi-object conjugate priors. IEEE Transactions on Signal Processing, 2013, 61(13): 3460-3475 doi: 10.1109/TSP.2013.2259822 [23] Reuter S, Vo B T, Vo B N, Dietmayer K. The labeled multi-Bernoulli filter. IEEE Transactions on Signal Processing, 2014, 62(12): 3246-3260 doi: 10.1109/TSP.2014.2323064 [24] 陈辉, 邓东明, 韩崇昭. 基于区间箱粒子多伯努利滤波器的传感器控制策略. 自动化学报, 2021, 47(6): 1428-1443 doi: 10.16383/j.aas.c180541Chen Hui, Deng Dong-Ming, Han Chong-Zhao. Sensor control based on interval box-particle multi-Bernoulli filter. Acta Automatica Sinica, 2021, 47(6): 1428-1443 doi: 10.16383/j.aas.c180541 [25] 侯利明, 连峰, 谭顺成, 徐从安. 闪烁噪声统计特性未知情况下的鲁棒广义标签多伯努利滤波器. 电子学报, 2021, 49(7): 1346-1353 doi: 10.12263/DZXB.20200960Hou Li-Ming, Lian Feng, Tan Shun-Cheng, Xu Cong-An. Robust generalized labeled multi-Bernoulli filter for multi-target tracking with unknown statistical characteristics of glint noise. Acta Electronica Sinica, 2021, 49(7): 1346-1353 doi: 10.12263/DZXB.20200960 [26] Zhao J, Gui R Z, Dong X D, Wu S Y. Time-varying DOA tracking algorithm based on generalized labeled multi-Bernoulli. IEEE Access, 2021, 9: 5943-5950 doi: 10.1109/ACCESS.2020.3048952 [27] Sathyan T, Chin T J, Arulampalam S, Suter D. A multiple hypothesis tracker for multitarget tracking with multiple simultaneous measurements. IEEE Journal of Selected Topics in Signal Processing, 2013, 7(3): 448-460 doi: 10.1109/JSTSP.2013.2258322 [28] Habtemariam B, Tharmarasa R, Thayaparan T, Mallick M, Kirubarajan T. A multiple-detection joint probabilistic data association filter. IEEE Journal of Selected Topics in Signal Processing, 2013, 7(3): 461-471 doi: 10.1109/JSTSP.2013.2256772 [29] Zhong X H, Premkumar A B, Madhukumar A S, Tong L C. Multi-modality likelihood based particle filtering for 2-D direction of arrival tracking using a single acoustic vector sensor. In: Proceedings of IEEE International Conference on Multimedia and Expo. Barcelona, Spain: IEEE, 2011. 1−6 [30] Zhong X H, Premkumar A B. Particle filtering approaches for multiple acoustic source detection and 2-D direction of arrival estimation using a single acoustic vector sensor. IEEE Transactions on Signal Processing, 2012, 60(9): 4719-4733 doi: 10.1109/TSP.2012.2199987 [31] Zhong X H, Hari V N, Premkumar A B, Madhukumar A S. Particle filtering with enhanced likelihood model for underwater acoustic source DOA tracking. In: Proceedings of the OCEANS IEEE-Spain. Santander, Spain: IEEE, 2011. 1−6 [32] Zhong X H, Premkumar A B, Madhukumar A S. Particle filtering and posterior cramér-rao bound for 2-D direction of arrival tracking using an acoustic vector sensor. IEEE Sensors Journal, 2012, 12(2): 363-377 doi: 10.1109/JSEN.2011.2168204 [33] Gunes A, Guldogan M B. Multi-target bearing tracking with a single acoustic vector sensor based on multi-Bernoulli filter. In: Proceedings of the OCEANS Genova. Genova, Italy: IEEE, 2015. 1−5 [34] Dong X D, Zhang X F, Zhao J, Sun M, Wu Q H. Multi-maneuvering sources DOA tracking with improved interactive multi-model multi-Bernoulli filter for acoustic vector sensor (AVS) array. IEEE Transactions on Vehicular Technology, 2021, 70(8): 7825-7838 doi: 10.1109/TVT.2021.3093063 [35] Fantacci C, Papi F. Scalable multisensor multitarget tracking using the marginalized δ-GLMB density. IEEE Signal Processing Letters, 2016, 23(6): 863-867 doi: 10.1109/LSP.2016.2557078 [36] Fantacci C, Vo B T, Papi F, Vo B N. The marginalized δ-GLMB filter. arXiv preprint arXiv: 1501.00926, 2015. [37] Wu S Y, Dong X D, Zhao J, Sun X Y, Cai R H. A fast implementation of interactive-model generalized labeled multi-Bernoulli filter for interval measurements. Signal Processing, 2019, 164: 345-353 doi: 10.1016/j.sigpro.2019.05.028 [38] Yi W, Jiang M, Hoseinnezhad R. The multiple model Vo–Vo filter. IEEE Transactions on Aerospace and Electronic Systems, 2017, 53(2): 1045-1054 doi: 10.1109/TAES.2017.2667300 [39] Mazor E, Averbuch A, Bar-Shalom Y, Dayan J. Interacting multiple model methods in target tracking: A survey. IEEE Transactions on Aerospace and Electronic Systems, 1998, 34(1): 103-123 doi: 10.1109/7.640267 [40] 李昊润, 卜凡康, 周剑雄. 修正的马尔科夫转移矩阵自适应IMM算法. 火力与指挥控制, 2021, 46(9): 118-124, 132 doi: 10.3969/j.issn.1002-0640.2021.09.021Li Hao-Run, Bu Fan-Kang, Zhou Jian-Xiong. Improved adaptive Markov matrix IMM algorithm. Fire Control & Command Control, 2021, 46(9): 118-124, 132 doi: 10.3969/j.issn.1002-0640.2021.09.021 [41] Agarwal A, Kumar A, Agrawal M, Fauziya F. Higher order statistics based direction of arrival estimation with single acoustic vector sensor in the under-determined case. In: Proceedings of the OCEANS MTS/IEEE Monterey. Monterey, USA: IEEE, 2016. 1−9 [42] Bao M, Zheng C S, Li X D, Yang J, Tian J. Acoustical vehicle detection based on bispectral entropy. IEEE Signal Processing Letters, 2009, 16(5): 378-381 doi: 10.1109/LSP.2009.2016014 [43] Mendel J M. Tutorial on higher-order statistics (spectra) in signal processing and system theory: Theoretical results and some applications. Proceedings of the IEEE, 1991, 79(3): 278-305 doi: 10.1109/5.75086 [44] Nikias C L, Mendel J M. Signal processing with higher-order spectra. IEEE Signal Processing Magazine, 1993, 10(3): 10-37 doi: 10.1109/79.221324 [45] Hoang H G, Vo B T, Vo B N. A fast implementation of the generalized labeled multi-Bernoulli filter with joint prediction and update. In: Proceedings of the 18th International Conference on Information Fusion (Fusion). Washington, USA: IEEE, 2015. 999−1006 [46] Vo B N, Vo B T, Hoang H G. An efficient implementation of the generalized labeled multi-Bernoulli filter. IEEE Transactions on Signal Processing, 2017, 65(8): 1975-1987 doi: 10.1109/TSP.2016.2641392 [47] Tichavsky P, Muravchik C H, Nehorai A. Posterior Cramer-Rao bounds for discrete-time nonlinear filtering. IEEE Transactions on Signal Processing, 1998, 46(5): 1386-1396 doi: 10.1109/78.668800 [48] Schuhmacher D, Vo B T, Vo B N. A consistent metric for performance evaluation of multi-object filters. IEEE Transactions on Signal Processing, 2008, 56(8): 3447-3457 doi: 10.1109/TSP.2008.920469 [49] Kuhn H W. The Hungarian method for the assignment problem. Naval Research Logistics Quarterly, 1955, 2(1-2): 83-97 doi: 10.1002/nav.3800020109 期刊类型引用(0)
其他类型引用(3)
-
下载:
计量
- 文章访问数: 603
- HTML全文浏览量: 106
- PDF下载量: 132
- 被引次数: 3
