近年来，国内的大地学者提出将陆地-岛礁-海面基准向海底基准传递，逐步建设健全海底大地基准网^[1-4]，以期为海洋安全、海洋经济开发、海洋科学研究、海洋环境监测等提供综合定位、导航、授时服务^[5]。1985年美国斯克利普斯海洋研究所Spiess教授^[6]首次提出采用全球导航卫星/声学定位(Global Navigation Satellite System-Acoustic，GNSS-A)相结合的方式建立海底大地基准，GNSS-A技术是GNSS动态定位和水声测距的一种组合技术，集卫星定位、水声定位和水声通讯于一体，可实现由卫星向水下或海底目标的实时高精度定位^[7-9]。21世纪以来，美国、加拿大、日本、俄罗斯等国家基于GNSS-A技术，通过周期性测量和标校实现长期维持高精度的海洋大地测量基准^[10-13]。我国的海底大地基准建设起步较晚，但已取得了比较大的突破，2019年7月在水深3 000 m海域进行了海底大地基准网布设及观测实验，通过GNSS-A技术对5个海底基准点开展圆走航加十字型的对称观测，获取了高精度的定位结果^[14]。

在水下定位技术中，安装于船底的声学换能器可实时测量与目标之间的传播时间，乘以声速就是两点之间的距离。该种模式结构简单、操作便利，大大提高了海底目标定位的精度和效率^{[15, 16]}。为了提高水下定位精度，国内外研究学者致力于削弱各类误差。Xu等^[17]提出了基于对称航迹的水下差分定位技术，使用差分算子削弱系统误差，可以提高GNSS-A技术的定位精度。一些学者将未建模的系统误差设成未知参数或者建模后估计，以消除对定位的影响^{[18, 19]}。Zhao等^[20]从随机模型出发，建立了基于声线入射角的水下定位随机模型，实验证明其优于传统的等权随机模型定位结果。邝英才等^[21]和Zhao等^[22]提出了一种GNSS-A整体解算方法，不同类型观测值纳入统一模型，提升了海底控制点的定位精度和稳定性。这些研究为我国海底控制网的建设打下了理论基础，也为本文提供了参考和借鉴。

基于GNSS-A技术在垂直方向上存在固有的几何缺陷，使得海底大地基阵网垂直方向的定位精度受以声速为主的误差影响明显^[23]。此外，受涌浪影响，船体姿态对水面GNSS定位精度的影响可由运动传感器实时改正，水下测距的精度主要受声速误差和噪声的影响。水平观测结构对称的测线可以减弱系统误差对观测精度的影响，同时抵消部分粗差，从而提高海底目标在水平方向的精度^{[24, 25]}；利用实测声速剖面开展声线跟踪可以进一步削弱声速误差的影响，最终提高海底目标的定位精度^{[26, 27]}。

本文基于GNSS-A绝对标校技术，研究提高海底目标垂直定位精度的解算方案。针对双程传播时延问题，给出一种基于动态观测的海底目标定位算法；针对声学测距易受声速误差影响的特性，结合实测声速剖面开展声线跟踪并构造迭代。通过现场实验，验证本文所提出的海底目标定位算法的精度与可靠性。

1 GNSS-A海底目标定位原理 1.1 GNSS-A技术原理

GNSS-A技术是综合利用GNSS、声波传播时间和声速剖面测量等方法计算海底目标的位置，主要分为两个部分：一是通过建立严密的站心地平坐标系，确定声学换能器在大地坐标系下的实时位置；二是通过声线追踪并构造迭代计算声学换能器与海底应答器之间的距离^[28]。如图 1所示，测量船上安装有GNSS接收机(可接收GNSS卫星信号)、运动传感器和安装于船底的声学换能器。测量船绕海底应答器航行时，声学换能器向海底发射声学信号，换能器接收到海底应答器的反馈信号并记录声波的传播时间，结合实测的声速剖面数据，计算声学换能器与海底应答器之间的距离^{[29, 30]}。GNSS接收机可实时获取测量船的位置信息，运动传感器可同步测量船体姿态信息。在船体坐标系中，GNSS天线与传感器之间的位置关系在安装前经过精确标定以及运动姿态传感器数据的补偿就可确定声学换能器的位置。

1.2 声学换能器实时动态定位

GNSS接收机可实时获取测量船的位置，确定船底换能器和GNSS接收机的相对位置关系，通过坐标改正就可以计算声学换能器的位置。由于GNSS接收机得到的是天线相位中心在WGS-84坐标系下的坐标，需将GNSS的坐标改化到声学换能器中心。船体坐标系与以GNSS接收机相位中心为原点的站心地平坐标系的位置关系如图 2(红线代表GNSS天线相位中心与船体坐标系中心的连线)所示^[31]。

计算声学换能器在以GNSS接收机相位中心为原点的站心地平坐标系的坐标改正数的计算方法如公式(1)所示。

$ \begin{array}{l} \ \ \ \ \ \ \ \ \left(\begin{array}{l} \triangle X \\ \triangle Y \\ \triangle Z \end{array}\right)_{\text {站心地平坐标系 }}=H_{\text {heading }} \cdot H_{\text {pitch }} \cdot H_{\text {roll }} \cdot\\ \left(\begin{array}{l} x \\ y \\ z \end{array}\right)_{船体坐标系}, \end{array} $

(1)

其中，H_heading、H_pitch、H_roll分别为

$ \begin{aligned} &\ \ \ \ \ \ \ \ \ \ H_{\text {heading}}=\\ &\left(\begin{array}{ccc} \cos (\text {heading}) & -\sin (\text {heading}) & 0 \\ \sin (\text {heading}) & \cos (\text {heading}) & 0 \\ 0 & 0 & 1 \end{array}\right), \end{aligned} $

(2)

$ \begin{aligned} &\ \ \ \ \ \ \ \ \ \ H_{\text {pitch}}=\\ &\left(\begin{array}{ccc} \cos (\text {pitch}) & 0 & -\sin (\text {pitch}) \\ 0 & 1 & 0 \\ \sin (\text {pitch}) & 0 & \cos (\text {pitch}) \end{array}\right), \end{aligned} $

(3)

$ \begin{aligned} &\ \ \ \ \ \ \ \ \ \ H_{\text {roll}}= \\ &\left(\begin{array}{ccc} 1 & 0 & 0 \\ 0 & \cos (\text {roll}) & -\sin (\text {roll}) \\ 0 & \sin (\text {roll}) & \cos (\text {roll}) \end{array}\right), \end{aligned} $

(4)

式中，heading为偏航角，船右舷为正、左舷为负；pitch为纵摇角，船头向上为正、向下为负；roll为横摇角，船右舷为正、左舷为负；(x    y    z)_{船体坐标系}^T是GNSS接收机与声学换能器在船体坐标系下的基线向量。

2 GNSS-A海底目标定位算法 2.1 基于动态观测的海底目标定位算法

基于海面浮标方式的海底声学定位，其声学换能器在一次测量的发射与接收过程中，因浮标的移动速度通常较小，一次测量中双程传播时延造成的位置变动可忽略，因此可将发射与接收当作单次观测；而在船载走航动态观测中，忽略双程传播时延，则会引入新的误差，即测量船在换能器发射信号到接收信号时处在不同位置，本文针对这种影响，将声学换能器动态观测中的发射与接收信号的过程分别考虑^[32]。

如图 3所示，海底应答器坐标为P_T，设在第k次观测中，声学换能器向海底应答器发射信号时位置为P_S^k，海底应答器接收到声信号并作出响应，此时声学换能器至海底应答器的距离为ρ_S^k，声学换能器接收到海底应答器返回信号时的位置为P_R^k，此时海底应答器至声学换能器的距离为ρ_R^k，则观测方程为

$ \rho_S^k+\rho_R^k=f\left(P_S^k, P_T^k\right)+f\left(P_R^k, P_T^k\right)+\varepsilon_u^k+\varepsilon^k, $

(5)

$ \begin{array}{l} \ \ \ \ \ \ \ \ f\left(P_S^k, P_T^k\right)= \\ \sqrt{\left(x_t^k-x_s^k\right)^2+\left(y_t^k-y_s^k\right)^2+\left(z_t^k-z_s^k\right)^2}, \end{array} $

(6)

$ \begin{array}{l} \ \ \ \ \ \ \ \ f\left(P_R^k, P_T^k\right)= \\ \sqrt{\left(x_t^k-x_r^k\right)^2+\left(y_t^k-y_r^k\right)^2+\left(z_t^k-z_r^k\right)^2}, \end{array} $

(7)

式中，f(P_S^k, P_T^k)，f(P_R^k, P_T^k)分别为换能器在发射和接收信号时与海底应答器的距离，ε_u^k为与声学换能器的测距误差，ε^k为其他噪声引起的测距随机误差^[33]。

利用声速剖面仪实时获取的声速剖面结合公式(8)计算Harmonic平均声速:

$ C_H=\left(z-z_0\right)\left[\int_{z_0}^z \frac{\mathrm{d} z}{C(z)}\right]^{-1}。$

(8)

根据泰勒级数公式在P_T0处展开，式(5)可线性化为

$ \begin{array}{l} \ \ \ \ \ \ \ \ \rho_{S_0}^k+\rho_{R_0}^k-f_S^0-f_R^0=a_s^k \mathrm{~d} X+a_r^k \mathrm{~d} X+\varepsilon_u^k+\\ \varepsilon^k, \end{array} $

(9)

式中，ρ^k_S0和ρ^k_R0分别代表声学换能器发射和接收声信号时至海底应答器的欧式距离，且

$ \left\{\begin{array}{l} \rho_{S_0}^k=C_H \cdot t_s^k \\ \rho_{R_0}^k=C_H \cdot t_r^k \end{array}, \right. $

(10)

式中，C_H为平均声速，t_s^k和t_r^k分别为换能器至应答器和应答器至换能器的单程传播时间。a_s^k和a_r^k为雅各比矩阵，系数为

$ \left\{\begin{array}{l} a_s^k=\left[\begin{array}{lll} \frac{\partial f_s^k}{\partial x_0} & \frac{\partial f_s^k}{\partial y_0} & \frac{\partial f_s^k}{\partial z_0} \end{array}\right] \\ a_r^k=\left[\begin{array}{lll} \frac{\partial f_r^k}{\partial x_0} & \frac{\partial f_r^k}{\partial y_0} & \frac{\partial f_r^k}{\partial z_0} \end{array}\right] \end{array}\right., $

(11)

$ \begin{aligned} &\left\{\begin{array}{c} \frac{\partial f_i^k}{\partial x_0}=\frac{x_{t_0}^k-x_i^k}{f_i^0} \\ \frac{\partial f_i^k}{\partial y_0}=\frac{y_{t_0}^k-y_i^k}{f_i^0} \\ \frac{\partial f_i^k}{\partial z_0}=\frac{z_{t_0}^k-z_i^k}{f_i^0} \\ f_i^0=\sqrt{\left(x_{t_0}^k-x_i^k\right)^2+\left(y_{t_0}^k-y_i^k\right)^2+\left(z_{t_0}^k-z_i^k\right)^2} \end{array}\right., \\ &i=s, r, \end{aligned} $

(12)

式中，(x^k_t0, y^k_t0, z^k_t0)为海底应答器的初始坐标值。根据式(9)得距离的残差v: v=ε_u^k+ε^k，利用最小二乘v^Tv=min解算海底应答器的位置。

2.2 基于声线跟踪的距离迭代算法

水平对称结构的观测可削弱系统误差在水平方向的影响，对于改善水平方向的定位精度具有较好的效果。但在垂直方向上，由于声速误差等系统误差的影响，精度改善并不明显^{[34, 35]}。为了研究该问题，本文给出基于实测声速剖面的距离迭代计算方法，通过提高声学测距精度，改善海底应答器垂直方向定位精度。

首先根据水深、换能器吃水深度和实测的声速剖面将换能器至应答器的距离划分为n层，利用计算得到的目标初始深度和测量船与海底应答器的平距计算初始入射角θ，基于层内常梯度声线跟踪方法计算水平、垂向分量和传播时间^{[36, 37]}：

$ \left\{\begin{array}{l} y_T=\sum\limits_{i=1}^n y_i, z_T=\sum\limits_{i=1}^n z_i \\ y_i=\frac{\sqrt{1-\left(p c_i\right)^2}-\sqrt{1-p^2\left(c_i+g_i \Delta z_i\right)^2}}{p g_i} \\ t_i=\frac{1}{\left|g_i\right|} \ln \left[\frac{c_i+g_i \Delta z_i}{c_i} \frac{1+\sqrt{1-\left(p c_i\right)^2}}{1+\sqrt{1-p^2\left(c_i+g_i \Delta z_i\right)^2}}\right] \end{array}\right., $

(13)

式中，y_i和t_i分别为经历水层厚度Δz_i的水平位移量和传播时间，p为snell常数；c_i、θ_i分别为波束经历第i层上界面的声速和入射角，c_i+1、θ_i+1分别为波束经历第i层下界面的声速和入射角，声速在层内以常梯度g_i=(c_i+1-c_i)/(z_i+1-z_i)变化，入射角可以根据以下公式计算：

$ \theta=-\arctan \left(y_T / z_T\right) 。$

(14)

利用解算后的坐标根据(14)式重新计算入射角，再结合观测时间开展声线跟踪更新应答器的坐标。重复上述过程，直到相邻两次确定的坐标小于给定的限差ε。

$ \left\{\begin{array}{l} y_T^{i+1}-y_T^i \leqslant \varepsilon_y \\ \left|z_T^{i+1}-z_T^i\right| \leqslant \varepsilon_z \end{array}\right., $

(15)

式中，(y, z)_Tⁱ和(y, z)_Tⁱ⁺¹分别为第i次和第i+1次信标点的定位结果，ε_y和ε_z为迭代终止坐标阈值。

本文提出的基于GNSS-A的海底目标定位算法如图 4所示，过程大致分3步。

第一步：根据GNSS接收机和声学换能器基阵之间的向量关系与实时获取的姿态数据，结合式(1)计算声学换能器在以GNSS接收机相位中心为原点的站心地平坐标系下的坐标改正数，进而求得声学换能器的实时位置P_S(X_S, Y_S, Z_S)和P_R(X_R, Y_R, Z_R)。

第二步：根据声学换能器测得的与海底应答器的传播时间t_S和t_R，结合声速剖面计算的平均声速C_H，计算声学换能器与海底应答器之间的距离ρ_S和ρ_R；采用换能器的平面坐标平均值作为海底应答器的平面坐标初始值(X_T0, Y_T0)、水深和换能器吃水深度的差值作为海底应答器初始垂直坐标Z_T0，采用最小二乘(Least Squares, LS)法并结合式(9)进行平差，计算海底应答器的坐标改正数(Δx, Δy, Δz)和残差向量V₁，得到平面坐标较准确的海底应答器坐标(X_T, Y_T, Z_T1)。

第三步：根据解算得到的海底应答器坐标(X_T, Y_T, Z_T1)和声速剖面数据开展声线跟踪，利用式(14)、式(15)迭代出更为准确的传播时间t和入射角θ，计算观测时间t_S和t_R与传播时间t的时间延迟ΔT，并依据$ \sum\limits_{i=1}^n \Delta T^2=\min$计算距离的残差向量V₂，进而求得Z_T。

3 验证实验

为验证本文定位算法的精度和可靠性，海洋地质九号船于2021年4月3日在中国南海进行了GNSS-A海底目标定位实验。实验涉及的仪器有GNSS星站差分定位仪Veripos LD5，水平和垂向定位精度分别为10 cm和20 cm；运动传感器Kongsberg MRU5，角速率动态精度和加速度精度分别为0.02°和0.01 m/s²；超短基线定位系统Kongsberg HiPap 102P，其声学换能器水下测距精度为0.02 m；移动船载声速测量系统AML MVP300，最大测量深度为3 350 m，声速测量精度为0.05 m/s，其搭载的四电极电导率传感器精度为0.001 S/m，精密热敏电阻温度传感器精度为0.005℃，应变片式压力传感器精度为测量水深的0.05%。实验开始前将声学信标绑定到可回收的装置中投入海底，并以此位置为圆心，分别设计两个不同半径的圆形对称实验测线。图 5所示的是实验区域及设计的测线和部分仪器。

实验区水深约1 750 m，为确保测量实验的精度和几何结构，设计以约0.6倍水深和1.2倍水深为半径的圆形测线，设计半径分别为1 000 m和2 000 m，可以保证声波入射角在合理范围内。

实验期间测量船以4节船速匀速行驶，顺序是测量船先沿大圆测线航行后沿小圆测线航行，图 6中红色圆圈是海底信标的位置，黑色圆圈是测量船在实验期间的实际航迹图，图 7为实验中声学换能器发射声波至海底信标的旅行时间的统计，大圆平均每次测距时间约为1.68 s(单程)，小圆平均每次测距时间约为1.33 s(单程)。移动船载声速测量系统获取的实验时间期间两次测量的声速剖面数据如图 8所示，两个声速剖面的平均声速约为1 492 m/s。

采用观测数据的平均值作为海底信标的初始坐标，结合观测数据代入观测式(9)中解算海底信标的位置，图 9为利用LS法计算的声学换能器与海底信标的距离残差。声学换能器历元间测量是独立同精度观测，因此对残差进行统计，并采用以下单位权中误差对测量结果精度评定^[38]：

$ \sigma=\sqrt{\frac{V^{\mathrm{T}} V}{n-t}}, $

(16)

式中，V为测量残差向量，t为必要观测的个数，本文t=3。

如图 9所示，红色框内是大圆测线解算的残差，黄色框内是小圆测线解算的残差，残差均匀分布在零轴附近，且小圆的残差小于大圆的残差。残差最大为0.98 m，最小为-1.09 m，在水深1 750 m的情况下，单位权中误差为0.29 m，说明观测结果具有较好的测量精度。对于海底信标的位置解算结果，采用北东天(North East Up，NEU)3个方向坐标分量和相对点位中误差对平差结果精度评定。平差参数中误差是评价观测数据内符合精度的指标，其值越小代表内符合精度越好，采用LS法解算的海底信标点的相对点位中误差为0.91 cm，北方向、东方向和垂直方向的误差分量分别为0.58 cm、0.60 cm和0.35 cm。

实验中垂直方向上的定位精度主要受声速误差的影响，因垂直方向上观测结构不对称，无法抵消垂直方向上的误差，采用声线跟踪可以减弱声速误差对观测距离的影响，提高海底信标在垂直方向上的定位精度。利用LS法解算的结果结合常梯度声线跟踪进行观测距离的迭代计算，图 10给出了初始入射角和迭代入射角的对比。

声线跟踪采用初始入射角进行迭代计算，利用迭代后的位置信息更新入射角，直至迭代出更准确的海底信标的位置。通过初始入射角和迭代后入射角的对比可以看出，大圆入射角为43-47°，小圆入射角为26-31°，迭代后入射角均略大于初始入射角，证明声线跟踪有效追踪声线在水中的实际传播路径，提高了观测距离的精度。图 11给出了本文方法和LS法解算观测距离的残差对比结果，表 1统计了两种方法的精度信息。

从对比结果可以看出，采用本文方法解算结果的残差分布更密集，残差的极值也均小于LS法的解算结果；在观测距离为2 600 m的情况下，单位权中误差从0.29 m减小到了0.20 m，说明本文方法可以有效减弱声速误差对测距精度的影响，提高观测的精度，进而提高海底信标在垂直方向的定位精度。

4 结论

针对声速误差对海底目标定位精度的影响，本文研究了GNSS-A海底目标高精度定位的算法，通过南海实测的深海实验数据验证了本文算法的有效性。实验结果表明：基于动态观测的海底目标定位算法能综合利用发射和接收过程中的观测数据，解决双程传播时延引起的位移误差问题，提高了海底目标的水平定位精度，获得了点位中误差优于1 cm的平差结果；基于声线跟踪的距离迭代算法可有效减弱声速误差对测距精度的影响，距离残差的单位权中误差从0.29 m减小到了0.20 m。本文算法可为提高海底目标垂直方向的定位精度提供参考。