CN105676181A - 基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法 - Google Patents

Abstract

本发明公开一种基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，首先引入两两传感器之间能量比的对数作为观测值，得到不同时刻的观测序列；然后结合目标的运动状态，将非线性的观测方程在误差可控范围内化为近似线性观测方程，构造线性状态-空间模型；然后导出状态-空间模型的扩展卡尔曼滤波序贯迭代求解算法；利用线性最小二乘(ER-LS)定位算法，对进入水下传感器网络范围内的运动目标定位，作为扩展卡尔曼滤波的初始值；最后，通过序贯滤波得到性能改进的水下目标运动轨迹。

Description

基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法

技术领域

本发明属于水下传感器网络框架下的滤波跟踪技术，具体涉及一种特定状态-空间模型下的扩展卡尔曼滤波跟踪方法。

背景技术

海洋覆盖了地球70％以上的面积，蕴藏了巨大的物质资源。由于海水的特殊性，海水对无线电波和可见光吸收太快，衰减迅速，从而不能实现能量的远距离传输，制约了这些技术在海洋探测和水下目标跟踪方面的应用。声波是唯一可以在水下进行远程传播的载体，不同的传播距离有着不同的应用场景，如传播距离为米到百米级可以利用声波进行高分辨力成像，传播距离为公里到百公里级可以利用声波进行目标探测和定位，传播距离为百公里到千公里级可以利用声波进行信息传输。而且，声波探测跟踪技术在声呐探潜、海洋生物跟踪等军事和民用领域具有广泛的应用。

无线传感器网络由多个通过无缆通信的传感器节点构成，为了实现其低功耗、快速布放、自组织性以及容错性的优势，节点通常由单个传感器构成，因而节点无法进行阵信号处理。传感器网络用于目标跟踪的方法有中心式和分布式之分，中心式指的是所有传感器将获得的观测值发送给指定的唯一中心节点处理，这种结构需要可观的宽带进行传感器节点之间的信息传输，这会对节点电池寿命提出挑战，中心式的结构弊端在于其鲁棒性和可靠性。为了减少中心式结构的通信负担，分布式结构应运而生，分布式结构每个节点先进行信号处理，压缩所需传输的数据量，再将处理结果发送给中心处理节点。通常选取靠近目标即信噪比的节点为中心处理节点，随着跟踪过程的进行，中心处理节点是不断变迁的，将各时刻用于搜集传感器信息的节点视为一个子集，同时将其它节点设置为节能睡眠模式，分布式结构较之中心式结构鲁棒性更好。

经典的声源目标定位方法主要有三类：基于波达方向定位(DOA)，基于时延定位(TDOA)和基于信号能量定位(RSSI)。RSSI对传输带宽以及传感器精度要求较低，更适合无线传感器网络。RSSI包括最大似然定位和最小二乘定位等，最大似然方法精度高、耗能大、计算量大，由于单个传感器节点体积小、能量有限，节能是首要考虑因素。引入能量比的概念，单个节点只需运行计算量很小的能量检测器，就可以得到有效压缩后的能量观测值数据，因而计算量可有效降低，采用最小二乘方法可实现目标位置估计。

发明内容

本发明的目的是针对水下传感器网络对运动目标的跟踪问题，提出一种特定的线性状态-空间模型以及与之适应的扩展卡尔曼滤波跟踪方法，可以有效地通过两两传感器之间能量比得到的观测序列，对运动状态进行时间更新和测量更新，通过序贯滤波得到性能改进的跟踪结果。

由于海洋环境的动态变化性以及目标的运动特性，使得海洋声信号在实际应用中常常涉及到非线性，构建的状态-空间模型通常是非线性的。而本发明提出的方法是对所有有效节点的能量观测值数据两两求比值再作对数变换，得到的能量差在对数尺度上可化为近似线性状态方程，序贯贝叶斯滤波在线性情况下的最佳解是扩展卡尔曼滤波(EKF)，能够很好的适应于这种非线性不严重的水下运动目标定位情形。

本发明的具体技术方案如下：

一种基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，包括步骤：

1)在水下布放由M个传感器节点组成的网络，每个节点由单个传感器组成，接收目标辐射出的声信号，对获得的声信号按一定频率进行采样，序贯计算各个节点固定时间段内声信号的能量，获得能量测量值；

2)根据目标在水下的运动规律建立运动状态方程；

3)根据步骤1)中各传感器节点获得的测量值，构建原始观测方程；以目标运动规律为依据对原始测量方程做线性化处理，得到线性的等效观测方程；

4)配置滤波迭代方程，得到状态估计的时间更新方程和测量更新方程、估计误差协方差的时间更新方程和量测更新方程；

5)根据步骤4)滤波迭代方程进行滤波迭代后,得到各时刻最优状态。

进一步的，在步骤1)中，在待测浅海海域海底平面布放的由M个传感器节点组成的网络，各个节点单独接收目标辐射出的声信号大小，根据目标运动状态选取适当的坐标系，假定海底平坦，每个传感器节点的位置由二维笛卡尔坐标系来表示：

(x_i,y_i),i＝1,2,…,M(1)

其中，x_i和y_i分别表示第i个传感器在二维平面内的横坐标和纵坐标。

当运动目标进入监控区域时，若节点全部可用，则可以得到M个接收信号，对获得的声信号按一定频率进行采样，序贯计算各个节点固定时间段内声信号的能量。

分布的各个传感器节点将接收到来自运动目标的信号，各个传感器节获得的目标辐射出的声能信号测量值为：

$E_i\left(k\right)={\mathrm{\γ}}_i\·\frac{S\left(k\right)}{\left|\right|L_i\left(k\right)-L_T\left(k\right)|{|}^{\mathrm{\α}/2}}+n_i\left(k\right)---\left(2\right)$

式中，E_i(k)表示第k时刻第i个传感器接收的目标辐射的能量大小；γ_i表示第i个传感器的接收增益；S(k)表示距离声源1m处获得的声能强度；L_T(k)和L_i(k)是3×1向量，分别表示声源目标和第i个传感器的位置，||L_T(k)-L_i(k)||是声源目标和第i个传感器的欧氏距离；α是衰减因子，与海底海面反射的影响有关；n_i(k)是方差为σ_i ²的零均值高斯白噪声，n_i(k)～N(0,σ²)，σ²为噪声方差。

进一步的，在步骤2)中，水下运动目标的速度一般都比较小，因此在较短的时间间隔内，此处特指相邻观测时间间隔内，其位移和深度等状态参量都可以近似为一阶线性变化附加二阶扰动项，因此单一目标的运动状态方程可以近似认为是线性的。假设运动目标在距离海平面等深处的二维平面内做匀速直线运动。设运动目标声源状态向量X(k)写作：

X(k)＝[x_T(k)v_x(k)y_T(k)v_y(k)]^T(3)

公式(3)中四个分量依次是运动目标的x轴坐标，沿x轴速度分量，y轴坐标，沿y轴速度分量。目标的运动状态方程写作：

X(k+1)＝Φ·X(k)+Γ·w(k)(4)

状态转移矩阵： $\mathrm{\Φ}=\left[\begin{array}{cccc}1& T& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& T\\ 0& 0& 0& 1\end{array}\right]$

T是相邻观测时间间隔。

过程噪声矩阵： $\mathrm{\Γ}=\left[\begin{array}{cc}{T}^2/2& 0\\ T& 0\\ 0& T^2/2\\ 0& T\end{array}\right]$

w(k)是概率分布为N(0,Q(k))的零均值高斯白噪声，

过程噪声的协方差矩阵： $Q\left(k\right)=\left[\begin{array}{cc}qT& 0\\ 0& qT\end{array}\right]$

q是过程噪声强度，T是相邻采样点间的时间间隔。

进一步的，在步骤3)中，观测向量为

Z_i(k)＝log{E_i(k)/E_m(k)}(5)

$Z\left(k\right)=\left[\begin{array}{c}log\left\{E_1\right(k)/E_m(k\left)\right\}\\ log\left\{E_2\right(k)/E_m(k\left)\right\}\\ ...\\ log\left\{E_i\right(k)/E_m(k\left)\right\}\\ ...\end{array}\right]=\left[\begin{array}{c}...\\ ...\\ Z_i\left(k\right)\\ ...\\ ...\end{array}\right]---\left(6\right)$

E_m(k)表示第k时间间隔内所有E_i(k)中的最大值，表示在第k时间间隔内第m个传感器接收能量最大；Z_i(k)表示第k时间间隔内第i个传感器与第m个传感器的能量比取对数；Z(k)表示第k时刻各观测向量Z_i(k)构成的列向量。

因为目标等速度运动，则等效观测方程简化为两维，观测向量以目标运动规律为依据对原始测量方程做线性化处理，得到线性的等效观测方程可以写作

ΔZ(k)＝H_x·Δx_T(k)+H_y·Δy_T(k)+v(k)(7)

其中，ΔZ(k)是观测向量在第k时间间隔内增量：ΔZ_i(k)＝Z_i(k)-Z_i(k-1)是一个标量；Z_i(k-1)表示第k-1时间间隔内第i个传感器与第m个传感器的能量比取对数；H_x和H_y分别表示x轴和y轴相关的观测矩阵；Δx_T(k)和Δy_T(k)分别表示第k时间间隔内x轴和y轴上的位移量；v(k)是零均值高斯白噪声，符合高斯概率密度分布N(0,R(k))；R(k)是观测噪声的协方差矩阵；

进一步的，在步骤3)中，观测模型x轴和y轴相关的观测矩阵H_x，H_y中对应的元素使用的初始估计值由基于能量比的最小二乘定位算法产生的。基于能量比的最小二乘定位算法本质是求解一个矩阵方程。

进一步的，在步骤4)中，目标运动只在x轴上做匀速直线运动，y轴位移分量和速度分量都为零，状态向量和状态运动方程如下，

X(k)＝[x_T(k)v_x(k)y_T(k)v_y(k)]^T(8)

X(k-1)＝[x_T(k-1)v_x(k-1)y_T(k-1)v_y(k-1)]^T(9)

ΔX(k)＝[Δx_T(k)000]^T(10)

状态估计的时间更新方程：

ΔX(k+1|k)＝ΔX(k|k)(11)

状态估计的测量更新方程：

ΔX(k+1|k+1)＝ΔX(k+1|k)+K(k+1)·[ΔZ(k+1)-H_x·ΔX(k+1|k)](12)

误差协方差的时间更新方程：

ΔP(k+1|k)＝ΔP(k|k)(13)

误差协方差的测量更新方程：

$\begin{array}{c}\mathrm{\ΔP}(k+1|k+1)=[I-K(k+1)\·H_x]\·\mathrm{\ΔP}(k+1|k)\\ \·{[I-K(k+1)\·H_x]}^T+K(k+1)\·R\·K{(k+1)}^T\end{array}---\left(14\right)$

卡尔曼增益：

K(k+1)＝ΔP(k+1|k)·H_x ^T·[H_x·ΔP(k+1|k)·H_x+R]^-1(15)

R表示过程噪声；ΔP表示Δx_T(k)的误差协方差矩阵。

本发明的有益效果是：本发明适用于目标相邻观测时刻内的位移远远小于目标与传感器节点距离的跟踪场景，根据运动状态构造线性状态-空间模型，通过序贯滤波得到相对精度改善的水下传感器网络框架下的基于能量比的扩展卡尔曼滤波跟踪方法。如果将非线性方程作为观测方程的话，就得求雅可比矩阵及其幂，计算量大，且只能在非线性程度轻微时候使用。本发明比传统的求雅可比矩阵的扩展卡尔曼滤波在计算量上有优势。由于本发明构造的状态-空间模型的应用有限制，需要假设在每一个相邻时间间隔中，运动目标移动的距离的绝对值远远小于它和参与跟踪的传感器节点之间的距离，因此对跟踪场景有限制，一旦突破这一限制，扩展卡尔曼滤波的优点无法凸显。

附图说明

图1为水下分布式传感器网络跟踪模型；

图2为M＝20个节点水下传感器网络基于能量比最小二乘定位算法结果与实际路径的比较；

图3为M＝20个节点水下传感器网络基于能量比最小二乘定位算法结果的均方误差；

图4为M＝20个节点水下传感器网络应用扩展卡尔曼滤波跟踪算法结果与实际路径的比较；

图5为M＝20个节点水下传感器网络应用扩展卡尔曼滤波跟踪算法结果的均方误差。

具体实施方式

下面结合附图和具体实例对本发明做进一步的描述，以验证本发明的有效性。图1为水下分布式传感器网络跟踪模型，由运动目标、头节点、激活节点和睡眠节点组成，在算法的各个定位周期中，自适应地激活一定范围内的一定数量的传感器节点，依靠多个激活的传感器协同合作来定位跟踪运动目标。具体的实施过程如下：

步骤1：设定目标运动状态模型，观测模型和海洋环境参数。海深均匀为28m，所用的海洋声速梯度是在舟山浅海实验所得，传感器网络有M＝20个节点，在深度28m海底处6000*6000范围内随机布放。运动目标声源在深度15m处做匀速直线运动，目标速度1m/s，假设运动目标发出的声源信号是一个1000Hz以下的宽带信号，幅值能量均值为5000，方差可调的信号。所有传感器节点的背景噪声设为方差为1均值为0的标准正态分布。因为目标在水深15m水平面内运动，可以用二维直角坐标系表示其位置，目标初始位置为(0，0)，每一次测量时间间隔T为1sec，且运动目标由初始位置沿x轴匀速运动。以上设置可以看出，满足在每一个相邻时间间隔中，运动目标移动的距离的绝对值远远小于它和参与跟踪的传感器之间的距离这一条件，符合专利所述的状态-空间模型的应用条件。

初始状态

X(0)＝[0100]^T

初始估计状态

$\hat{X}\left(0\right)={\left[\begin{array}{cccc}0& 0& 0& 0\end{array}\right]}^T$

初始状态误差协方差矩阵

$P\left(0\right)=\left[\begin{array}{cc}100& 0\\ 0& 100\end{array}\right]$

测量噪声协方差矩阵

$R\left(k\right)=\left[\begin{array}{cc}20& 0\\ 0& 20\end{array}\right]$

过程噪声协方差矩阵

$Q\left(k\right)=\left[\begin{array}{cc}10& 0\\ 0& 10\end{array}\right]$

步骤2：采用分布式定位跟踪方法，共20个传感器节点，用基于能量比的最小二乘定位算法产生初始估计值。在算法的各个定位周期中，自适应地激活一定范围内的一定数量的传感器，依靠多个激活的传感器协同合作来定位跟踪运动目标。由于在跟踪的全过程中，每个跟踪时刻不参与跟踪的节点被看作睡眠节点，因此公式(30)中常常是求伪逆的过程。图2显示最小二乘定位算法的结果，有部分时刻单个时刻定位误差达到2米多，由于各时刻是单独估计的，所以定位结果按时间轴看起来不平滑，没有连贯性。图3显示最小二乘定位算法部分时刻的目标轨迹均方误差最大可达5米多。

步骤3：运动状态时间更新和测量更新。根据状态转移矩阵和状态噪声对初始估计进行时间更新。根据传感器网络各节点获得的每一时刻的全新测量值来更新初始状态变量估计值，根据观测方程中的观测矩阵以及观测噪声，状态方程中的过程噪声来计算卡尔曼增益，由此更新初始误差协方差矩阵估计值。观测矩阵在目标匀速直线运动情况下近似线性，矩阵中各个常量使用的初始估计值由基于能量比的最小二乘定位算法产生。图4给出扩展卡尔曼滤波方法得到的目标轨迹曲线。由于采用序贯滤波，利用了目标运动的演化知识，因而得到的是平滑的曲线。图5为扩展卡尔曼滤波方法的跟踪误差，各时刻的跟踪误差都小于1，比最小二乘估计结果有了显著改善。可见，使用了构造的线性状态-空间模型基础上的扩展卡尔曼滤波，跟踪结果得到了平滑，跟踪精度得以提高。

以上所述仅为本发明的较佳实施举例，并不用于限制本发明，凡在本发明精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，其特征在于，包括步骤：

1)在水下布放由M个传感器节点组成的网络，各传感器节点用以接收目标辐射出的声信号，并对获得的声信号按频率进行采样，然后序贯计算各个传感器节点对应时间段内声信号的能量测量值；

2)根据目标在水下的运动规律建立运动状态方程；

3)根据步骤1)中各传感器节点获得的能量测量值，计算M个传感器节点中每两个节点间声能的比值并取对数，作为观测量，构建原始观测方程；以目标运动规律为依据对原始观测方程做线性化处理，得到线性的等效观测方程；

4)配置滤波迭代方程，得到状态估计的时间更新方程和测量更新方程、估计误差协方差的时间更新方程和量测更新方程；

5)根据步骤4)滤波迭代方程进行滤波迭代后，得到各时刻最优状态。

2.根据权利要求1所述的基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，其特征在于，所述的步骤1)中，各个传感器节点获得的目标辐射出的声信号为：

$E_i\left(k\right)={\mathrm{\γ}}_i\·\frac{S\left(k\right)}{\left|\right|L_i\left(k\right)-L_T\left(k\right)|{|}^{\mathrm{\α}/2}}+n_i\left(k\right)$

式中，E_i(k)表示第k时刻第i个传感器接收的目标辐射的能量大小；γ_i表示第i个传感器的接收增益；S(k)表示距离声源1m处获得的声能强度；L_T(k)和L_i(k)是3×1向量，分别表示声源目标和第i个传感器的位置，||L_T(k)-L_i(k)||是声源目标和第i个传感器的欧氏距离；α是衰减因子，与海底海面反射的影响有关；n_i(k)是方差为σ_i ²的零均值高斯白噪声，n_i(k)～N(0,σ²)，σ²为噪声方差。

3.根据权利要求1所述的基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，其特征在于，所述的步骤2)中，目标的运动状态是匀速直线运动。

4.根据权利要求1所述的基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，其特征在于，所述的步骤3)中，目标等速度运动，则等效观测方程简化为两维，表示为：

ΔZ(k)＝H_x·Δx_T(k)+H_y·Δy_T(k)+v(k)

其中

$\mathrm{\Δ}Z\left(k\right)=\left[\begin{array}{c}{\mathrm{\ΔZ}}_1\left(k\right)\\ {\mathrm{\ΔZ}}_2\left(k\right)\\ ...\\ {\mathrm{\ΔZ}}_i\left(k\right)\\ ...\end{array}\right]=\left[\begin{array}{c}{Z}_1\left(k\right)-Z_1(k-1)\\ Z_2\left(k\right)-Z_2(k-1)\\ ...\\ Z_i\left(k\right)-Z_i(k-1)\\ ...\end{array}\right]+v\left(k\right)$

Z_i(k)为能量比测量值，即为步骤3)中的观测量，表示第k时间间隔内第i个传感器的能量E_i(k)与第m个传感器的能量E_m(k)比值取对数：

Z_i(k)＝log{E_i(k)/E_m(k)}

式中，Z_i(k-1)表示第k-1时间间隔内第i个传感器与第m个传感器的能量比取对数；ΔZ_i(k)＝Z_i(k)-Z_i(k-1)表示Z_i(k)与Z_i(k-1)的差值，是一个标量；ΔZ(k)表示第k时间间隔内所有传感器节点各自ΔZ_i(k)构成的向量；H_x和H_y分别表示x轴和y轴相关的观测矩阵；Δx_T(k)和Δy_T(k)分别表示第k时间间隔内目标二维坐标中x轴和y轴上的位移量；v(k)为零均值高斯白噪声，符合高斯概率密度分布N(0,R(k))，R(k)为协方差矩阵。

5.根据权利要求4所述的基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，其特征在于，所述观测矩阵中使用的初始估计值由基于能量比的最小二乘定位算法产生。

6.根据权利要求1所述的基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，其特征在于，所述的步骤3)中，要满足观测方程线性化的条件是：目标单位时间内位移远远小于传感器距离目标的距离，

x_T(k)＝x_T(k-1)+Δx_T(k),y_T(k)＝y_T(k-1)+Δy_T(k)

Δx_T(k)≈0,Δy_T(k)≈0

式中，x_T(k)和y_T(k)分别表示第k时刻目标所在位置。

7.根据权利要求1所述的基于分布式传感器能量比的水下运动目标扩展卡尔曼滤波跟踪方法，其特征在于，在所述的步骤4)中，建立坐标系，令目标运动只在x轴上存在位移，y轴位移分量为零，则状态估计的时间更新方程：

Δx_T(k+1|k)＝Δx_T(k|k)

状态估计的测量更新方程：

Δx_T(k+1|k+1)＝Δx_T(k+1|k)+K(k+1)·[ΔZ(k+1)-H_x·Δx_T(k+1|k)]

误差协方差的时间更新方程：

ΔP(k+1|k)＝ΔP(k|k)

误差协方差的测量更新方程：

ΔP(k+1|k+1)＝[I-K(k+1)·H_x]·ΔP(k+1|k)

·[I-K(k+1)·H_x]^T+K(k+1)·R·K(k+1)^T

卡尔曼增益：

K(k+1)＝ΔP(k+1|k)·H_x ^T·[H_x·ΔP(k+1|k)·H_x+R]^-1

I表示单位矩阵；

R表示过程噪声；

ΔP表示Δx_T(k)的误差协方差矩阵。