CN103323815B - 一种基于等效声速的水下声学定位方法 - Google Patents

Abstract

本发明涉及一种基于等效声速的水下声学定位方法，其特征在于，包括以下步骤：步骤1：初始化，设定水下目标源的初始位置，在水下随意放置n个传感器接收器；步骤2：到达时间差TDoA的测量；步骤3：计算有效声音速度；步骤4：利用卡尔曼滤波器对目标源进行定位，减小背景噪声的影响，提高水下定位的精度。与现有技术相比，本发明具有同步三维定位、实时速度计算、定位低误差、高精度等优点。

Description

一种基于等效声速的水下声学定位方法

技术领域

本发明涉及水下声学定位领域，尤其是涉及一种基于等效声速的水下声学定位方法。

背景技术

精确的水下定位对于水下传感器网络的诸多应用非常重要，如水下精确导航、海洋资源开采、追踪海洋生物以及监视环境安全等等。为了从水下传感器网络中的水下传感器节点和水下自主航行器(AUV)得到更加有意义的数据，需要有更加准确的水下定位系统方案。

然而，水下环境的定位不同于陆地的定位，现有的水下航行器的定位方法大致可以分为三类：基于水下传感器网络的声学定位、基于地形信息的匹配定位以及最新的即时定位和构图的方法。其中最为常用的是水下声学定位方法。水下声学定位方法通过水下传感器网络实现，水下传感器网络使用声波通信，传播介质是水，这与无线电波在空中传播有显著的不同，它具有高时延，时延不断动态变化的特点。水声信号传播面临着衰减大、通信信道带宽低、多径效应严重、误码率高等问题。

目前，为水下传感器网络设计的声学定位方法大多数是基于距离的，基于距离的方法可以得到更加精确的结果。基于距离的定位方法通过测量水下传感器网络中节点的时间或者角度，转化为相应节点之间的距离值，通过多边形运算，推算出未知节点的位置。大部分水下基于距离的定位方法采用的是基于时间的测量方法，也就是通过测量信号到达时间(ToA)或者成对的信号到达时间差，再把时间乘以水下声音传播的速度，得到相应的距离信息。

然而，现有的定位系统有它们自身的缺陷，不能够为在重要的水下任务中追踪水下航行器提供足够精确的定位。首先，通过单次的信号到达时间差得出的定位存在很大的误差，这种误差主要是由于环境的噪声造成的。其次，在现有的基于距离的定位系统中，定位方法都假定水下声音传播的速度是恒定的，为1500m/s，同时，这些方法还假定在水下三维环境中，声音是沿着直线传播的。但是，这些假设都有一定的局限性，在对定位精度要求比较高的水下任务中是不能够成立的。

水下环境中，声音传播速度随着水下温度，压力以及盐度不断变化，产生声音传播的弯折现象。因此，即使有非常精准的水下同步方法，保证水下传感器节点在水下能够同步，并且传感器的测量误差非常小，声音传播模型的错误估计还是会给定位带来很大的误差。

发明内容

本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种同步三维定位、实时速度计算、定位低误差、高精度的基于等效声速的水下声学定位方法。

本发明的目的可以通过以下技术方案来实现：

一种基于等效声速的水下声学定位方法，其特征在于，包括以下步骤：

步骤1：初始化，设定水下目标源的初始位置p_m(x_m，y_m，z_m)，在水下随意放置n个传感器接收器即接收节点，记录其位置分别为p₁(x₁，y₁，z₁)，p₂(x₂，y₂，z₂)，…，p_n(x_n，y_n，z_n)；

步骤2：到达时间差TDoA的测量；

步骤3：计算有效声音速度；

步骤4：利用卡尔曼滤波器对目标源进行定位，减小背景噪声的影响，提高水下定位的精度。

所述的到达时间差的计算采用广义交叉相关方法(GCC)，已知目标源发出声波s(t)，定义第i个接收点接收到的声波信号为：r_i(t)＝h_i(t)*s(t)+n_i(t)，第j个接收点接收到的声波信号为：r_j(t)＝h_j(t)*s(t-τ)+n_j(t)，其中，n_i(t)和n_j(t)分别表示第i个接收点和第j个接收点接收到得声波的加性噪声，假设n_i(t)和n_j(t)是不相关的，τ表示两个接收节点收到信号的时间延迟，r_i(t)和r_j(t)的互相关函数为：

${\hat{R}}_{r_i,r_j}\left(\mathrm{\τ}\right)=\frac{1}{T-\mathrm{\τ}}\underset{\mathrm{\τ}}{\overset{T}{\∫}}{r}_1\left(t\right)r_2(t-\mathrm{\τ})\mathrm{dt}$

其中，T表示观测周期；

首先，对r_i(t)通过滤波器H1(f)和r_j(t)通过滤波器H2(f)进行预滤波，得到y_i(t)和y_j(t)，然后计算两者的互相关函数的最大值：

令 $\frac{\∂{\hat{R}}_{r_i,r_j}\left(\mathrm{\τ}\right)}{\∂\mathrm{\τ}}=0,$

此时τ的值即为函数的最大值；

然后，当滤波器H1(f)＝H2(f)＝1时，求出互相关函数峰值出现时的τ值，即时间延迟τ的估计值： $\widehat{\mathrm{\τ}}=\arg \underset{\mathrm{\τ}}{\max }{R}_{y_i{y}_j}\left(\mathrm{\τ}\right).$

所述的计算有效声音速度具体为：根据上次定位的两个接收点到目标源的距离差除以声波到达两接收点的时间差求得此次有效声音速度，在时刻t(k-1)，定位系统定位到的目标源位置为p_m(k-1)，第i个接收点的位置为p_i(k-1)，第1个接收点的位置为p₁(k-1)，定义第i个接收点到目标源的距离与第1个接收点到目标源的距离差为：

ΔR_i1(k-1)＝||p_i(k-1)-p_m(k-1)||-||p_l(k-1)-p_m(k-1)||，

则在时刻t(k)，也就是系统对目标源进行第k次定位，利用第k-1次定位的ΔR_i1(k-1)计算出第k次定位的有效声音速度

${\hat{c}}_i\left(k\right)=\frac{\mathrm{\Δ}{R}_{i1}(k-1)}{{\mathrm{\τ}}_i(k-1)},$ 其中

τ_i(k-1)为在时刻t(k-1)，第i个接收节点和第一个接收节点收到信号的信号到达时间差。

利用卡尔曼滤波器对目标源进行定位为：将多次信号到达时间差(TDoA)和有效声音速度融合到一次定位过程中。

所述的将多次信号到达时间差(TDoA)和有效声音速度融合到一次定位过程中，具体计算如下：

$z\left(k\right)=\frac{1}{2}\left[\begin{array}{c}{K}_2^2-R_{21}^2\left(k\right)-2R_1\left(k\right)R_{21}\left(k\right)\\ K_3^2-R_{31}^2\left(k\right)-2R_1\left(k\right)R_{31}\left(k\right)\\ \·\\ K_i^2-R_{i1}^2\left(k\right)-2R_1\left(k\right)R_{i1}\left(k\right)\\ \·\\ K_n^2-R_{n1}^2\left(k\right)-2R_1\left(k\right)R_{n1}\left(k\right)\end{array}\right]$

其中z(k)为测量向量，是第i个接收节点的三维位置，R_i1(k)是第i个接收节点和目标源的距离与第一个接收节点和目标源的距离的差值， $R_{i1}\left(k\right)=\mathrm{\τ}\left(k\right){\hat{c}}_i\left(k\right),$ i＝1...n。

与现有技术相比，本发明引入有效声音速度计算方法，实现对水下声速的实时估计，优化了水下声音传播模型，同时利用卡尔曼滤波器进行信息融合，减小了背景噪声的影响，提高了水下定位的精度。

附图说明

图1为本发明的系统结构示意图。

图2为本发明的信号时延估计的框图。

图3为本发明的定位过程框图。

图4为本发明在实例中的x-y平面的定位结果。

图5为本发明在实例中的深度方向的定位结果。

图6为本发明与经典定位方法在实例中的定位结果比较图。

具体实施方式

下面结合附图和具体实施例对本发明进行详细说明。

实施例1

一种基于等效声速的水下声学定位方法，包括以下步骤：

步骤1：初始化，设定水下目标源的初始位置p_m(x_m，y_m，z_m)，在水下随意放置n个传感器接收器即接收节点，记录其位置分别为p₁(x₁，y₁，z₁)，p₂(x₂，y₂，z₂)，…，p_n(x_n，y_n，z_n)；

步骤2：到达时间差(TDoA)的测量，到达时间差的计算采用广义交叉相关(GCC)方法，已知目标源发出声波s(t)，定义第i个接收点接收到的声波信号为：r_i(t)＝h_i(t)×s(t)+n_i(t)，第j个接收点接收到的声波信号为：r_j(t)＝h_j(t)×s(t-τ)+n_j(t)，其中，n_i(t)和n_j(t)分别表示第i个接收点和第j个接收点接收到得声波的加性噪声，假设n_i(t)和n_j(t)是不相关的，τ表示两个接收节点收到信号的时间延迟，r_i(t)和r_j(t)的互相关函数为：

${\hat{R}}_{r_i,r_j}\left(\mathrm{\τ}\right)=\frac{1}{T-\mathrm{\τ}}\underset{\mathrm{\τ}}{\overset{T}{\∫}}{r}_1\left(t\right)r_2(t-\mathrm{\τ})\mathrm{dt}$

其中，T表示观测周期；

首先，对r_i(t)通过滤波器H1(f)和r_j(t)通过滤波器H2(f)进行预滤波，得到y_i(t)和y_j(t)，

然后，当滤波器H1(f)＝H2(f)＝1时，求出互相关函数峰值出现时的τ值，即时间延迟τ的估计值： $\widehat{\mathrm{\τ}}=\arg \underset{\mathrm{\τ}}{\max }{R}_{y_i{y}_j}\left(\mathrm{\τ}\right).$

步骤3：计算有效声音速度，根据上次定位的两个接收点到目标源的距离差除以声波到达两接收点的时间差求得此次有效声音速度，在时刻t(k-1)，定位系统定位到的目标源位置为p_m(k-1)，第i个接收点的位置为p_i(k-1)，第1个接收点的位置为p₁(k-1)，定义第i个接收点到目标源的距离与第1个接收点到目标源的距离差为：

ΔR_i1(k-1)＝||p_i(k-1)-p_m(k-1)||-||p₁(k-1)-p_m(k-1)||，

则在时刻t(k)，也就是系统对目标源进行第k次定位，利用第k-1次定位的ΔR_i1(k-1)计算出第k次定位的有效声音速度：

步骤4：运用卡尔曼滤波器进一步提高定位精度。

在一次定位运算中，仅采用当前时刻的TDoA的结果进行位置估计，有较大的误差，主要来源于传感器环境周围的测量噪声。本发明提出用卡尔曼滤波的方法，将多次信号到达时间差(TDoA)和有效声音速度融合到一次定位过程中，达到更加精确的定位结果。

定义系统的测量向量：

$z\left(k\right)=\frac{1}{2}\left[\begin{array}{c}{K}_2^2-R_{21}^2\left(k\right)-2R_1\left(k\right)R_{21}\left(k\right)\\ K_3^2-R_{31}^2\left(k\right)-2R_1\left(k\right)R_{31}\left(k\right)\\ \·\\ K_i^2-R_{i1}^2\left(k\right)-2R_1\left(k\right)R_{i1}\left(k\right)\\ \·\\ K_n^2-R_{n1}^2\left(k\right)-2R_1\left(k\right)R_{n1}\left(k\right)\end{array}\right]$

其中，p_i(x_i，y_i，z_i)是第i个接收节点的三维位置，i＝1...n是第i个接收节点和AUV的距离与第一个接收节点和AUV距离的差值。通过把等效声速的概念应用到定位方法中，可以得出对AUV位置的估计。根据多边定位方法，p₁＝[000]，定义向量h＝[p₂p₃…p_n]，第k次对水下AUV位置估计的结果为p_m(k)＝(h^Th)^-1h^Tz(k)。因此，卡尔曼滤波的观测方程是：

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

其中， $H\left(k\right)={\left[\begin{array}{cccccc}{x}_2& y_2& z_2& 0& 0& 0\\ x_3& y_3& z_3& 0& 0& 0\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ x_n& y_n& z_n& 0& 0& 0\end{array}\right]}_{n\×6,}$ $S\left(k\right)={\left[\begin{array}{c}{x}_m\left(k\right)\\ y_m\left(k\right)\\ z_m\left(k\right)\\ v_x\left(k\right)\\ v_y\left(k\right)\\ v_z\left(k\right)\end{array}\right]}_{6\×1}$

v(k)是测量噪声向量，一般情况下，它是高斯分布的，均值为0。在每一个离散的定位时间点T_k，定位系统都执行一次基于卡尔曼滤波的定位方法。

本发明中，每隔T分钟估计一次水下AUV的位置p_m(x_m，y_m，z_m)和速度v_m(v_x，v_y，v_z)。用S(k)＝[p_m(k)v_m(k)]来表示要被估测的量，卡尔曼滤波的状态方程：

S(k+1)＝AS(k)+U(k)+w(k)

其中，矩阵A是：

$A=\left[\begin{array}{cccccc}1& 0& 0& T& 0& 0\\ 0& 1& 0& 0& T& 0\\ 0& 0& 1& 0& 0& T\\ 0& 0& 0& 1& 0& 0\\ 0& 0& 0& 0& 1& 0\\ 0& 0& 0& 0& 0& 1\end{array}\right]$

U(k)是输入的部分，w(k)是过程噪声，通常可以被认为是高斯分布的，均值为0，协方差矩阵为Q(k)。在没有关于AUV航行更多的信息的情况下，可以设定系统的输入，U(k)＝0。每隔T分钟，通过定位方法更新一次对AUV位置的估计。

定义为卡尔曼滤波器对水下AUV运动参数向量S的最优估计，P为协方差矩阵，表示卡尔曼滤波估计的不确定性。卡尔曼滤波器方程如下：

1)一步状态估计：

S_k＝AS_k-1

2)误差协方差矩阵的估计：

P(k)_(-)＝P(k-1)₍₊₎+Q(k)

3)最佳增益：

K(k)＝P(k)_(-)A^T[AP(k_)(-)A^T+R]^-1

4)更新的状态估计

S_k(+)＝S_k(-)+K(k)(z(k)-AS_k(-))

5)更新后的误差协方差矩阵

${P\left(k\right)}_{(+)}^{-1}={P\left(k\right)}_{(-)}^{-1}+A^T{R}^{-1}A$

将实时计算得到的有效声音速度运用到对目标源的定位中去，每隔一定时间对目标源进行一次定位更新，每一次定位信息可以作为下一次定位使用。

实施例2

如图1所示，以水下航行器AUV为例，一种基于等效声速的水下声学定位方法包括以下步骤：

步骤1：进行传感器接收节点的布置和参数统计。设定实施例的初始状态：水下航行器AUV最初的位置是p_m＝[1000，1000，10]^T，四个水下接收节点的位置分别为：p₁＝[0，0，0]^T，p₂＝[2000，2000，5]^T，p₃＝[1000，0，100]^T，p₄＝[0，1000，100]^T，在大部分的航行时间内，水下航行器都是沿着直线运动的，本实施例AUV的速度取v_m＝＝[-20，10，1]^Tm/s，负方向表示深度的增加。

步骤2：到达时间差(TDoA)的测量。通过图2所示的过程，相乘，求平均，最大化输出，即能求出两个信号到达的时间差。通过把两个接收节点收到的信号从时域转换到频域，加上合适的加权，就可以通过DFT反变换得到广义交叉相关性函数的峰值在时域的位置。估计的信号时间延迟就是τ的值，即互相关函数峰值出现的地方。

步骤3：计算有效声音速度。在仿真环境下，根据水下节点的坐标和TDoA的值计算得到声音在水下传播的等效声速。

步骤4：运用卡尔曼滤波器进行定位，并提高定位精度。通过把这些信号到达时间差融合到卡尔曼滤波器的观测方程中，定位方法可以给出AUV运动参数的估计。在卡尔曼滤波器中，过程噪声协方差矩阵被设定为Q＝10I₆，其中I₆为6×6的单位矩阵，初始的误差协方差矩阵被定为P(0)＝100I₆。假定对信号到达时间差的估计的误差是高斯分布的，均值为0，方差为σ＝0.01m/s²。

通过水下仿真来评估本发明技术方案的定位方法的表现。当一个水下航行器以匀速沿着直线运动时，运用本发明，对AUV的运动参数进行估计。通过每隔5分钟进行一次定位，一共运行30次定位方法。

图4是定位方法对AUV在水平面的定位结果，这里没有AUV的深度信息。图5是AUV在深度方向的定位结果。圆点连接的实线表示AUV在每个时刻真实的深度信息。用三角连接起来的实线表示采用本发明得出的对AUV位置深度信息的估测值。注意，在定位过程中本发明与传统方法不同的是，AUV的三维位置是在定位过程中同时被估测出来的，而不是通过投影的方法先估测二维信息，然后用深度计估测深度信息。可以看出，本发明的定位方法得出的结果和AUV真实的位置非常接近。

图6是与经典定位方法的比较，这部分通过实验仿真，比较本发明的定位方法和经典的基于最小二乘方法的定位方法。图中，方块是本发明的定位误差，圆点虚线是用传统的基于最小二乘方法的定位误差，可得出本发明显著地减小了水下定位的误差。在传统的定位方法中，定位结果只考虑当前采集到的信号到达时间差的信息，本发明中，将现有的测量结果与之前信息进行融合，通过卡尔曼滤波后，可以做到更加准确的位置估计，达到更加精确的定位结果。

Claims (2)

1.一种基于等效声速的水下声学定位方法，其特征在于，包括以下步骤：

步骤1：初始化，设定水下目标源的初始位置p_m(x_m,y_m,z_m)，在水下随意放置n个传感器接收器即接收节点，记录其位置分别为p₁(x₁,y₁,z₁),p₂(x₂,y₂,z₂),…,p_n(x_n,y_n,z_n)；

步骤2：到达时间差TDoA的测量；

步骤3：计算有效声音速度；

步骤4：利用卡尔曼滤波器对目标源进行定位，减小背景噪声的影响，提高水下定位的精度；

所述的计算有效声音速度具体为：根据上次定位的两个接收点到目标源的距离差除以声波到达两接收点的时间差求得此次有效声音速度，在时刻t(k-1)，定位系统定位到的目标源位置为p_m(k-1)，第i个接收点的位置为p_i(k-1)，第1个接收点的位置为p₁(k-1)，定义第i个接收点到目标源的距离与第1个接收点到目标源的距离差为：

R_i1(k-1)＝||p_i(k-1)-p_m(k-1)||-||p₁(k-1)-p_m(k-1)||，

则在时刻t(k)，也就是系统对目标源进行第k次定位，利用第k-1次定位的R_i1(k-1)计算出第k次定位的有效声音速度

${\hat{c}}_i\left(k\right)=\frac{R_{i1}(k-1)}{{\mathrm{\τ}}_i(k-1)},$ 其中

τ_i(k-1)为在时刻t(k-1)，第i个接收节点和第一个接收节点收到信号的信号到达时间差；

利用卡尔曼滤波器对目标源进行定位为：将多次信号到达时间差TDoA和有效声音速度融合到一次定位过程中；

所述的将多次信号到达时间差TDoA和有效声音速度融合到一次定位过程中，具体计算如下：

$z\left(k\right)=\frac{1}{2}\left[\begin{array}{c}{K}_2^2-R_{21}^2\left(k\right)-2R_1\left(k\right)R_{21}\left(k\right)\\ K_3^2-R_{31}^2\left(k\right)-2R_1\left(k\right)R_{31}\left(k\right)\\ .\\ K_i^2-R_{i1}^2\left(k\right)-2R_1\left(k\right)R_{i1}\left(k\right)\\ .\\ K_n^2-R_{n1}^2\left(k\right)-2R_1\left(k\right)R_{n1}\left(k\right)\end{array}\right]$

其中z(k)为测量向量，是第i个接收节点的三维位置，R_i1(k-1)是第i个接收节点和目标源的距离与第一个接收节点和目标源的距离的差值， $R_{i1}(k-1)={\mathrm{\τ}}_i(k-1){\hat{c}}_i\left(k\right),i=1...n.$

2.根据权利要求1所述的一种基于等效声速的水下声学定位方法，其特征在于，所述的到达时间差的计算采用广义交叉相关方法GCC，已知目标源发出声波s(t)，定义第i个接收点接收到的声波信号为：r_i(t)＝h_i(t)*s(t)+n_i(t)，第j个接收点接收到的声波信号为：r_j(t)＝h_j(t)*s(t-τ)+n_j(t)，其中，n_i(t)和n_j(t)分别表示第i个接收点和第j个接收点接收到得声波的加性噪声，假设n_i(t)和n_j(t)是不相关的，τ表示两个接收节点收到信号的时间延迟，r_i(t)和r_j(t)的互相关函数为：

${\hat{R}}_{r_i,r_j}\left(\mathrm{\τ}\right)=\frac{1}{T-\mathrm{\τ}}\underset{\mathrm{\τ}}{\overset{T}{\∫}}{r}_1\left(t\right)r_2(t-\mathrm{\τ})dt$

其中，T表示观测周期；

首先，对r_i(t)通过滤波器H1(f)和r_j(t)通过滤波器H2(f)进行预滤波，得到y_i(t)和y_j(t)，然后计算两者的互相关函数的最大值：

令 $\frac{\∂{\hat{R}}_{r_i,r_j}\left(\mathrm{\τ}\right)}{\∂\mathrm{\τ}}=0,$

此时τ的值即为函数的最大值；

然后，当滤波器H1(f)＝H2(f)＝1时，求出互相关函数峰值出现时的τ值，即时间延迟τ的估计值： $\widehat{\mathrm{\τ}}=\arg \underset{\mathrm{\τ}}{max}{R}_{y_i{y}_j}\left(\mathrm{\τ}\right).$