CN110703187B - 一种基于水下机动平台的声信标测向及定位方法 - Google Patents

Abstract

本发明公开了一种基于水下机动平台的声信标测向及定位方法，步骤1：信号测向：基于迭代更新的思想使波束方向逐步接近声信标方向，并通过多项式求根计算更新波束方向；步骤2：声信标方向跟踪：建立一阶角度跟踪模型，并通过粒子滤波算法利用历史信息实现声信标测向结果的平滑及测向精度的进一步提高；步骤3：声信标定位：使用新型角交汇定位模型，建立解算方程，并采用拟牛顿迭代法进行非线性方程组求解，实现定位解算。经过以上声信标方向根据后，可以充分挖掘了历史信息，提升了测向精度；且减少了野值点，降低了野值点对声信标定位的影响。

Description

一种基于水下机动平台的声信标测向及定位方法

技术领域

本发明属于声信标测向及定位的技术领域；具体涉及一种基于水下机动平台的声信标测向及定位方法。

背景技术

基于水下机动平台的声信标测向及定位，指的是利用水下机动平台搭载水听器阵列对声信标进行信号方位估计，进而根据水下机动平台在不同位置的方位估计结果解算声信标位置的技术。水下机动平台声信标测向与定位具有被动实施、定位精度高、鲁棒性好等优点。因此，广泛应用于黑匣子等声信标搜探、水下通信设备侦察等场合，具有重要的民用及军用价值。

现有的声信标测向方法其波束角度是离散的。当声信标方向与离散的波束角度不一致时，会造成测向精度的下降。此外，在测向时没有充分利用历史信息，测向精度有进一步提高的空间。利用目声信标测向信息进行定位，现有的方法通常假定声信标与水下机动平台在同一平面内。在实际情况中，声信标位于海底，而水下机动平台距离海底有一定距离。两者的失配会造成较大的解算误差。因此，无论对于声信标测向还是定位，其精度均有较大的改进空间。本专利涉及一种基于水下机动平台的声信标测向及定位方法，可以有效提高声信标测向精度及定位精度。

发明内容

本发明涉及一种基于水下机动平台的声信标测向及定位方法，水下机动平台在不同位置接收声信标所发射的水声信号，测得信号方位分别为θ₁,θ₂,...,θ_n；通过历史信息对信号方位进行跟踪，提高声信标的测向精度；利用不同位置的测向结果结合水下机动平台自身惯性导航信息解算声信标位置，实现定位功能。

本发明通过以下技术方案实现：

一种基于水下机动平台的声信标测向及定位方法，所述方法包括以下步骤：

步骤1：信号测向：基于迭代更新的思想使波束方向逐步接近声信标方向，并通过多项式求根计算更新波束方向；

步骤2：声信标方向跟踪：建立一阶角度跟踪模型，并通过粒子滤波算法利用历史信息实现声信标测向结果的平滑及测向精度的进一步提高；

步骤3：声信标定位：使用新型角交汇定位模型，建立解算方程，并采用拟牛顿迭代法进行非线性方程组求解，实现定位解算。

所述步骤1中的测向波束为单频声波信号，所述信号模型如下：

X＝a(θ_n)s+N (1)

X为N×L维接收信号矩阵，N为接收水听器个数，L为快拍数；s为1×L维发射信号矢量，发射信号为CW脉冲信号，则其发射信号矢量为：

a(θ_n)为N×1维阵列导向矢量，本发明采用均匀线阵，则其导向矢量为：

式中：ω_θ＝(dω₀/c)sinθ_n为空间频率，c为水下声速，d为阵元间距，θ_n为信号方向，[*]^T表示矩阵转置；N为N×L维噪声矩阵。

进一步的，所述步骤1的通过迭代逐步更新搜索角度与信号方向的差值Δ，更新值由信号后验概率分布的最大期望值给出，包括以下步骤：

步骤1.1：初始化方差矩阵R、噪声功率

搜索角度与信号方向差Δ，

方差矩阵为对角线矩阵，对角线元素初值可设为非负的随机数，噪声功率初值由环境噪声估算，搜索角度与信号方向差为信号搜索间隔的一半；

步骤1.2：求取信号后验概率分布的均值与方差；

其中，信号后验概率分布的方差为：

式中：

信号后验概率分布的均值为：

步骤1.3：更新方差矩阵R、噪声功率

搜索角度与信号方向差δ；

其中，方差矩阵R的更新值为：

R＝diag([γ₁,γ₂,...,γ_M]) (8)

式中：||*||₂表示矩阵的2-范数；*_(i.)表示矩阵中第i行的行向量；*_(i,i)表示矩阵中第i行第i列的值；

噪声功率

的更新值为：

式中：||*||_F表示矩阵的F-范数；Tr[*]表示矩阵的迹；

搜索角度与信号方向差δ的更新值为：

Δ＝-U^-1V (10)

式中：Tr[*]表示实部；*_(.j)表示矩阵中第j列的列向量。

表示共轭；

步骤1.4：重复第二步及第三步，直至Δ的变化量小于某一预设门限；

步骤1.5：计算信号方位。收敛之后，则与信号方向最接近的搜索角度为：

进而，信号方位估计结果为：

进一步的，所述步骤2的一阶角度跟踪模型具体为：

根据声信标方向跟踪，所建立的状态转移方程为：

Zⁿ⁺¹＝HZⁿ+R (14)

其中：

为状态矩阵，θ_n'为声信标方向变化率；

为状态转移矩阵；R为状态转移噪声矩阵；

相应的观测方程为：

Yⁿ＝θ_n+Q (15)

其中：Q为观测噪声。

进一步的，所述步骤2建立以上跟踪模型后，采用粒子滤波原理进行声信标方向跟踪的具体过程为:

步骤2.1：粒子初始化：令状态矩阵初值为

对状态矩阵进行粒子采样，令每个粒子的权值相等，即：

其中：K为采样粒子个数；

步骤2.2：预测：基于前一的状态值Z^i-1对当前时刻的状态Zⁱ⁺¹进行预测采样，即：Zⁱ(k)～P(Zⁱ|Z^i-1(k))；

步骤2.3：更新权值：采用下式计算每个粒子的权值并归一化：

w_i*(k)＝w_i-1(k)P(Y_i|X_i(k)) (16)

式中：P(*|*)表示条件概率分布；

步骤2.4：状态估计：当前时刻的状态估计值为所有粒子的加权平均，即：

步骤2.5：重采样：令

G(k)～U(0,1)，式中：U(0,1)表示0至1区间内的均匀分布，

则经过重采样后的第k个新粒子为：X_i(L(k))。重采样后将粒子的权值平均化，即：

步骤2.6：每当有新的观测样本，重复步骤2.2至步骤2.5，则经过跟踪后的声信标方向为：

进一步的，所述步骤2角交汇定位模型所建立三维定位解算方程为：

式中：[x_s,y_s]为待求解的声信标水平坐标；z_s为声信标深度，由压力传感器测量；[x₁,y₁,z₁]及[x₂,y₂,z₂]为水下机动平台在不同测量点的三维坐标，由惯导系统测量；[x_a,y_a,z_a]为水听器阵列的单位方向矢量，由水下机动平台的姿态测量仪测量；

及

为声信标与水下机动平台的距离。

进一步的，所述步骤3中新型角交汇定位模型，即公式(18)的具体求解为：

步骤3.1：确定未知数初值：

步骤3.2：计算偏微分矩阵初值：

步骤3.3：计算函数值初值：

步骤3.4：更新未知数：

步骤3.5：更新函数值：

步骤3.6：计算残差：rⁱ＝xⁱ⁺¹-xⁱ，yⁱ＝F(xⁱ⁺¹)-F(xⁱ)；

步骤3.7：更新偏微分矩阵：

步骤3.8：重复步骤3.4至步骤3.7，直至两次迭代的间隔小于某一门限β；

即：(xⁱ⁺¹-xⁱ)·(xⁱ⁺¹-xⁱ)^T＜β时，停止迭代；停止迭代时未知数的数值即为声信标水平坐标的定位结果。

本发明的有益效果是：

经过以上声信标方向根据后，可以带来两个显著优点；一方面，充分挖掘了历史信息，提升了测向精度；另一方面，减少了野值点，降低了野值点对声信标定位的影响。

附图说明

图1本发明工作态势图。

图2本发明流程图。

图3本发明水池模拟试验信号测向结果，图3-(a)试验结果图，图3-(b)结果误差图。

图4本发明水池模拟试验声信标方位跟踪结果，图4-(a)试验结果图，图4-(b)结果误差图。

图5本发明水池模拟试验声信标定位结果，图5-(a)试验结果图，图5-(b)结果误差图。

具体实施方式

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

实施例1

图1所示为基于水下机动平台的声信标测向及定位工作态势图。水下机动平台在不同位置接收声信标所发射的水声信号，测得信号方位分别为θ₁,θ₂,...,θ_n；通过历史信息对信号方位进行跟踪，提高声信标的测向精度；利用不同位置的测向结果结合水下机动平台自身惯性导航信息解算声信标位置，实现定位功能。

一种基于水下机动平台的声信标测向及定位方法，所述方法包括以下步骤：

步骤1：信号测向：测向波束角度的离散化，采用期望-最大化算法，基于迭代更新的思想使波束方向逐步接近声信标方向，并通过多项式求根计算更新波束方向；经过若干次迭代，当波束方向与声信标方向一致时，达到较高的测向精度；

步骤2：声信标方向跟踪：考虑到水下机动平台的运动速度有限，建立一阶角度跟踪模型，并通过粒子滤波算法利用历史信息实现声信标测向结果的平滑及测向精度的进一步提高；

步骤3：声信标定位：针对声信标及水下机动平台非共面问题，本发明提出了一种新型角交汇定位模型，建立解算方程，并采用拟牛顿迭代法进行非线性方程组求解，实现定位解算。

所述步骤1中的测向波束为单频声波信号，所述信号模型如下：

X＝a(θ_n)s+N (1)

X为N×L维接收信号矩阵，N为接收水听器个数，L为快拍数；s为1×L维发射信号矢量，发射信号为CW脉冲信号，则其发射信号矢量为：

a(θ_n)为N×1维阵列导向矢量，本发明采用均匀线阵，则其导向矢量为：

式中：ω_θ＝(dω₀/c)sinθ_n为空间频率，c为水下声速，d为阵元间距，θ_n为信号方向，[*]^T表示矩阵转置；N为N×L维噪声矩阵。

所述步骤1中的测向波束信号方向在角度搜索空间具有离散性，采用离散表示的方法来估计方位；离散表示的模型如下：

X＝AS+N (2)

A为N×M维字典矩阵，M为搜索角度数；S为M×L维系数矩阵；由于搜索角度是离散化的，当信号方向与搜索角度不一致时，会造成信号方向估计精度的下降；

令β＝{β₁,β₂,…,β_M}为离散的搜索角度，当

时，假定β_n是距离信号方向最近的搜索角度；利用泰勒展开对导向矢量进行近似表示，有：

a(θ_n)≈a(β_n)+b(β_n)(θ_n-β_n) (3)

b(β_n)是a(β_n)对β_n的导数，进而，(2)式可表示为：

B＝[b(β₁),b(β₂),…,b(β_M)]；Δ＝diag[δ]；δ＝[δ₁,δ₂,...,δ_M]；

diag[*]表示对角阵。

进一步的，所述步骤1的通过迭代逐步更新搜索角度与信号方向的差值Δ，更新值由信号后验概率分布的最大期望值给出，包括以下步骤：

步骤1.1：初始化方差矩阵R、噪声功率

搜索角度与信号方向差Δ，

方差矩阵为对角线矩阵，对角线元素初值可设为非负的随机数，噪声功率初值由环境噪声估算，搜索角度与信号方向差为信号搜索间隔的一半；

步骤1.2：求取信号后验概率分布的均值与方差；

其中，信号后验概率分布的方差为：

式中：

信号后验概率分布的均值为：

步骤1.3：更新方差矩阵R、噪声功率

搜索角度与信号方向差δ；

其中，方差矩阵R的更新值为：

R＝diag([γ₁,γ₂,...,γ_M]) (8)

式中：||*||₂表示矩阵的2-范数；*_(i.)表示矩阵中第i行的行向量；*_(i,i)表示矩阵中第i行第i列的值；

噪声功率

的更新值为：

式中：||*||_F表示矩阵的F-范数；Tr[*]表示矩阵的迹；

搜索角度与信号方向差δ的更新值为：

Δ＝-U^-1V (10)

式中：Tr[*]表示实部；*_(.j)表示矩阵中第j列的列向量。

表示共轭；

步骤1.4：重复第二步及第三步，直至Δ的变化量小于某一预设门限；

步骤1.5：计算信号方位。收敛之后，则与信号方向最接近的搜索角度为：

进而，信号方位估计结果为：

进一步的，所述步骤2的一阶角度跟踪模型具体为：

根据声信标方向跟踪，所建立的状态转移方程为：

Zⁿ⁺¹＝HZⁿ+R (14)

其中：

为状态矩阵，θ_n'为声信标方向变化率；

为状态转移矩阵；R为状态转移噪声矩阵；

相应的观测方程为：

Yⁿ＝θ_n+Q (15)

其中：Q为观测噪声。

进一步的，所述步骤2建立以上跟踪模型后，采用粒子滤波原理进行声信标方向跟踪的具体过程为:

步骤2.1：粒子初始化：令状态矩阵初值为

对状态矩阵进行粒子采样，令每个粒子的权值相等，即：

其中：K为采样粒子个数；

步骤2.2：预测：基于前一的状态值Z^i-1对当前时刻的状态Zⁱ⁺¹进行预测采样，即：Zⁱ(k)～P(Zⁱ|Z^i-1(k))；

步骤2.3：更新权值：采用下式计算每个粒子的权值并归一化：

w_i*(k)＝w_i-1(k)P(Y_i|X_i(k)) (16)

式中：P(*|*)表示条件概率分布；

步骤2.4：状态估计：当前时刻的状态估计值为所有粒子的加权平均，即：

步骤2.5：重采样：令

G(k)～U(0,1)，式中：U(0,1)表示0至1区间内的均匀分布，

则经过重采样后的第k个新粒子为：X_i(L(k))；重采样后将粒子的权值平均化，即：

步骤2.6：每当有新的观测样本，重复步骤2.2至步骤2.5，则经过跟踪后的声信标方向为：

进一步的，所述步骤2角交汇定位模型所建立三维定位解算方程为：

式中：[x_s,y_s]为待求解的声信标水平坐标；z_s为声信标深度，由压力传感器测量；[x₁,y₁,z₁]及[x₂,y₂,z₂]为水下机动平台在不同测量点的三维坐标，由惯导系统测量；[x_a,y_a,z_a]为水听器阵列的单位方向矢量，由水下机动平台的姿态测量仪测量；

及

为声信标与水下机动平台的距离。

进一步的，所述步骤3中新型角交汇定位模型，即公式(18)的具体求解为：

步骤3.1：确定未知数初值：

步骤3.2：计算偏微分矩阵初值：

步骤3.3：计算函数值初值：

步骤3.4：更新未知数：

步骤3.5：更新函数值：

步骤3.6：计算残差：rⁱ＝xⁱ⁺¹-xⁱ，yⁱ＝F(xⁱ⁺¹)-F(xⁱ)；

步骤3.7：更新偏微分矩阵：

步骤3.8：重复步骤3.4至步骤3.7，直至两次迭代的间隔小于某一门限β；

即：(xⁱ⁺¹-xⁱ)·(xⁱ⁺¹-xⁱ)^T＜β时，停止迭代；停止迭代时未知数的数值即为声信标水平坐标的定位结果。

综上所述，本专利所述基于水下机动平台的未知周期声信标高精度定位方法流程图如图2所示。

实施例2

水池模拟试验的参数如下：水下机动平台航迹为圆形轨迹，速度为2m/s，半径为500m，深度为100m；声信标位于圆形轨迹圆心位置，深度为200m并保持静止不动。声信标发射CW脉冲信号，周期1s，脉宽15ms，频率30kHz。水中声速为1500m/s，声信号接收信噪比为0dB。其它输入参数的估计误差为：声速估计误差1.5m/s，深度估计误差1m，惯导误差0.2％。

首先采用本发明涉及的方法进行信号测向。信号测向结果如图3所示，图中的平均测向精度为0.98°。可见，采用本发明涉及的信号测向方法可以实现高精度的信号到达角度测量。

其次采用本发明涉及的方法进行声信标方向跟踪。声信标方向跟踪的结果如图4所示，图中的平均跟踪精度为0.34°。可见，采用本发明涉及的声信标方向跟踪方法可以进一步提高声信标角度估计精度。

最后采用本发明涉及的声信标定位方法进行声信标定位解算。声信标定位解算的结果如图5所示，图中稳定后的平均定位误差为4.58m。可见，采用本发明涉及的声信标定位方法可以获得高精度的声信标位置信息，达到了发明目的。

Claims (3)

1.一种基于水下机动平台的声信标测向及定位方法，其特征在于，所述方法包括以下步骤：

步骤1：信号测向：基于迭代更新的思想使波束方向逐步接近声信标方向，并通过多项式求根计算更新波束方向；

步骤2：声信标方向跟踪：建立一阶角度跟踪模型，并通过粒子滤波算法利用历史信息实现声信标测向结果的平滑及测向精度的进一步提高；

步骤3：声信标定位：使用新型角交汇定位模型，建立解算方程，并采用拟牛顿迭代法进行非线性方程组求解，实现定位解算；

所述步骤2的一阶角度跟踪模型具体为：

根据声信标方向跟踪，所建立的状态转移方程为：

Zⁿ⁺¹＝HZⁿ+R (14)

其中：

为状态矩阵，θ_n'为声信标方向变化率；

为状态转移矩阵；R为状态转移噪声矩阵；

相应的观测方程为：

Yⁿ＝θ_n+Q (15)

其中：Q为观测噪声；

所述步骤3角交汇定位模型所建立三维定位解算方程为：

式中：[x_s,y_s]为待求解的声信标水平坐标；z_s为声信标深度，由压力传感器测量；[x₁,y₁,z₁]及[x₂,y₂,z₂]为水下机动平台在不同测量点的三维坐标，由惯导系统测量；[x_a,y_a,z_a]为水听器阵列的单位方向矢量，由水下机动平台的姿态测量仪测量；

及

为声信标与水下机动平台的距离；

所述步骤3中新型角交汇定位模型，即公式(18)的具体求解为：

步骤3.1：确定未知数初值：

步骤3.2：计算偏微分矩阵初值：

步骤3.3：计算函数值初值：

步骤3.4：更新未知数：

步骤3.5：更新函数值：

步骤3.6：计算残差：rⁱ＝xⁱ⁺¹-xⁱ，yⁱ＝F(xⁱ⁺¹)-F(xⁱ)；

步骤3.7：更新偏微分矩阵：

步骤3.8：重复步骤3.4至步骤3.7，直至两次迭代的间隔小于某一门限β；

即：(xⁱ⁺¹-xⁱ)·(xⁱ⁺¹-xⁱ)^T＜β时，停止迭代；停止迭代时未知数的数值即为声信标水平坐标的定位结果。

2.根据权利要求1所述的测向及定位方法，其特征在于，所述步骤1的通过迭代逐步更新搜索角度与信号方向的差值Δ，更新值由信号后验概率分布的最大期望值给出，包括以下步骤：

步骤1.1：初始化方差矩阵R、噪声功率

搜索角度与信号方向差Δ，

方差矩阵为对角线矩阵，对角线元素初值可设为非负的随机数，噪声功率初值由环境噪声估算，搜索角度与信号方向差为信号搜索间隔的一半；

步骤1.2：求取信号后验概率分布的均值与方差；

其中，信号后验概率分布的方差为：

式中：

信号后验概率分布的均值为：

步骤1.3：更新方差矩阵R、噪声功率

搜索角度与信号方向差Δ；

其中，方差矩阵R的更新值为：

R＝diag([γ₁,γ₂,...,γ_M]) (8)

式中：||*||₂表示矩阵的2-范数；*_(i.)表示矩阵中第i行的行向量；*_(i,i)表示矩阵中第i行第i列的值；

噪声功率

的更新值为：

式中：||*||_F表示矩阵的F-范数；Tr[*]表示矩阵的迹；

搜索角度与信号方向差Δ的更新值为：

Δ＝-U^-1V (10)

式中：Tr[*]表示实部；*_(.j)表示矩阵中第j列的列向量，

表示共轭；

步骤1.4：重复第二步及第三步，直至Δ的变化量小于某一预设门限；

步骤1.5：计算信号方位， 收敛之后，则与信号方向最接近的搜索角度为：

进而，信号方位估计结果为：

3.根据权利要求2所述的测向及定位方法，其特征在于，所述步骤2建立以上跟踪模型后，采用粒子滤波原理进行声信标方向跟踪的具体过程为:

步骤2.1：粒子初始化：令状态矩阵初值为

对状态矩阵进行粒子采样，令每个粒子的权值相等，即：

其中：K为采样粒子个数；

步骤2.2：预测：基于前一的状态值Z^i-1对当前时刻的状态Zⁱ⁺¹进行预测采样，即：Zⁱ(k)～P(Zⁱ|Z^i-1(k))；

步骤2.3：更新权值：采用下式计算每个粒子的权值并归一化：

w_i*(k)＝w_i-1(k)P(Y_i|X_i(k)) (16)

式中：P(*|*)表示条件概率分布；

步骤2.4：状态估计：当前时刻的状态估计值为所有粒子的加权平均，即：

步骤2.5：重采样：令

G(k)～U(0,1)，式中：U(0,1)表示0至1区间内的均匀分布，

则经过重采样后的第k个新粒子为：X_i(L(k))， 重采样后将粒子的权值平均化，即：

步骤2.6：每当有新的观测样本，重复步骤2.2至步骤2.5，则经过跟踪后的声信标方向为：

Images