水下定位与导航是海洋探测的重要内容，在海洋工程建设、海洋科学研究、海洋权益维护等方面应用广泛^[1-3]。电磁波在水体中快速衰减的性质使全球导航卫星系统(global navigation satellite system, GNSS)难以直接用于水下定位，而声波在海水中良好的传播性质使水下声学定位系统得到了广泛研究和应用，如长基线系统(long baseline system, LBL)、短基线系统(short baseline system，SBL)、超短基线系统(ultra short baseline system，USBL)和GNSS智能浮标系统(GNSS intelligent buoys system，GIB)^[4-5]。其中，超短基线系统通过测角和测距实现对目标的跟踪定位，精度相对较低但操作维护简便^[6]；长基线系统、短基线系统以及GNSS智能浮标系统都是采用距离交会定位原理，多余观测量的存在使得定位精度明显提高^[7-9]。

对于采用距离交会原理的水下声学定位系统，误差来源主要有声学时延的测量误差、仪器的安装校准偏差、声速相关误差等^[10-11]。对于前两类误差，目前可通过精细校准、模型改进等方法有效消除或减弱。对于声速相关误差，又可细分为两类。其一是声波折射效应引起的误差，一般考虑结合实测声速剖面采用声线跟踪计算的方式消除折射效应的影响^[12-14]。其二是由于非实时实地测量声速剖面引起的声速剖面代表性误差，但受制于目前的声速测量手段，实现一定范围内海水声速的时空连续监测非常困难^[15-16]。为此，一种处理思路是通过构建三维声速场获取更加精细的声速分布信息^[17]；另一种是通过构建模型或增加约束的方式来消除误差影响，代表性的方法包括测距误差拟合^[18]、水下差分定位^[19-20]、深度约束定位^[21]等。显然，无论实测声速剖面精确与否，或采用何种水下定位模型和误差模型，利用声线跟踪消除折射效应影响是需要解决的一个基础性问题。本文重点讨论水声测距定位的声线跟踪问题。

在声速剖面准确的情况下，声线跟踪是消除声波折射效应最直接的方法^[22-23]，尤以常梯度声线跟踪方法应用最为广泛^[24-25]，有常梯度层内弧长积分和常梯度层内弦长积分两种不同的积分方法，实际计算差别不大但推荐理论上更严密的弧长积分方法^[26]。常梯度声线跟踪计算精度较高但计算量较大。等效声速剖面法用一个积分面积相同的常梯度声速剖面代替复杂的实际声速剖面^[27]，可以有效提高计算的效率，但由于忽略了波束入射角变化对计算的影响，因此本质上是一个近似计算模型。对于采用距离交会原理的水下声学定位，无论采用何种声线跟踪方法，所面临的最大难题就是波束入射角未知的问题。现有方法中已考虑到波束入射角计算的是一种基于双曲线定位(单差分定位)的声线跟踪方法^[28]，而对于最具一般性的非差分定位声线跟踪方法研究较少。针对上述问题，本文重点研究面向高精度水下声学定位的常梯度声线跟踪定位算法，基于一般性的非差分定位模型提出了定位解算与声线跟踪的迭代计算方法。本文方法的核心是波束入射角的反算问题，可为相关的科学研究和实际解算提供参考。

1 水下声学定位方法

距离交会定位是水声定位最普遍的定位形式，本文以水下基准点的定位计算为例对相关方法进行阐述。如图 1所示，测量船可利用搭载的GNSS接收机获得高精度的位置信息，再结合船体姿态和位置参数，通过坐标转换获得声学换能器在T_i时刻的坐标(x_i, y_i, z_i)，观测历元i=1, 2, …, m，声学换能器向水下基准点发射声信号并接收返回信号，通过测量声信号由发射到接收的时间延迟2t_i，在声速c已知的情况下，可以得到声学换能器与基准点的距离ρ_i=ct_i。

水下声学定位的非差分定位模型可以表示为

式中，δd_i为换能器时间延迟测量相关的系统性误差，如收发换能器、水下基准点声学设备的响应误差；δp_i为位置参数校准相关的系统性误差，如仪器校准不准确产生的位置偏差；δc_i为声速相关的系统性误差，如声速剖面不准确产生的声速剖面代表性误差；ε_i为测量偶然误差；(x_t, y_t, z_t)为水下基准点的坐标；距离交会方程f_i(x_t, y_t, z_t)可以表示为

在忽略各类系统误差的情况下，非差定位的最小二乘解为

式中，估计参数X=[Δx  Δy  Δz]^T；系数矩阵

观测向量

由上述可知，非差定位本质上忽略了系统误差的影响，为保证水下定位精度，必须在测量前对声学仪器的响应误差、仪器的安装校准偏差进行控制，并采集测区作业期间较多的声速剖面以减少声速剖面代表性误差的影响。此外，由式(2)可知，只有当距离观测值ρ_i为换能器与水下基准点间的直线距离时，才能保证水下定位的精度。但海水为非均匀介质体，海水中的声速存在时空上的变化，因在垂直方向上变化明显，所以通常用包含了水深和声速信息的声速剖面来表示，此时，声速c显然不是一个常值，采用平均声速等近似声速无法完全消除折射效应的影响。

2 常梯度声线跟踪模型

由于声波在水体中的传播性质较好，因此声学设备被广泛应用于各类海洋探测中，如单波束测深仪、多波束测深系统、侧扫声呐以及各类水下声学定位系统。但海水声速的时空变化使得声波的传播路径并非直线，而是发生连续折射的曲线，此即声波的折射效应。对于海洋测量，一般是将几何意义上的直线距离作为观测量，需要利用实测声速剖面进行声线跟踪计算，将连续弯曲的声线改正为直线距离。

目前，常用的声线跟踪模型主要有常梯度声线跟踪法和等效声速剖面法。其中，常梯度声线跟踪法的计算精度最高但计算量较大；等效声速剖面法是一种近似模型，提高了计算效率却也损失了一定的精度。本文主要探讨面向高精度水下声学定位的声线跟踪定位算法，因此选用常梯度声线跟踪法开展研究。

常梯度声线跟踪是将实际声速剖面分成若干层，假设层内声速以常梯度变化，采用层追加的方法计算声线的传播时间、垂向距和侧向距。如图 2所示，若第k层内深度区间为(z_k, z_k+1)，声速区间为(c_k, c_k+1)，波束入射角为θ_k，k=0, 1, …，n，根据射线声学的理论，声线在声速常梯度层中的轨迹为圆弧，圆弧半径R_k为

式中，p为Snell常数；g_k为层内声速梯度，可分别表示为

第k层声线的垂向距、侧向距可以表示成

改正后的距离为

声波沿圆弧从s_k传播到s_k+1所需要的时间t_k可表示为

反推得由t_k计算波束入射角θ_k+1的公式

由上述推导过程可以看出，采用常梯度声线跟踪法的前提是已知声速梯度g、时间延迟t以及波束入射角θ₀。对于采用距离交会定位的水下声学定位系统，声速梯度g可由实测声速剖面计算得到，时间延迟t为系统的观测值，只有波束入射角θ₀无法直接获得，因而成为水下定位声线跟踪的难点。

3 常梯度声线跟踪水下定位算法

对于水下声学定位的声线跟踪，关键是如何计算得到准确的波束入射角，一种几何近似计算的方法是利用声学换能器在T_i时刻的坐标(x_i, y_i, z_i)和水下基准点的坐标(x_t, y_t, z_t)来推算近似的波束入射角，可以表示为

这种几何近似计算很难保证声线跟踪计算的精度。从声线跟踪的过程也可以看出，波束入射角与目标坐标是相关的，因此二者应该是同步迭代的关系，即利用声线跟踪修正距离观测值，从而提高定位计算的精度，再由定位计算结果推算更为准确的波束入射角，进而提高声线跟踪计算的精度。

3.1 搜索法确定波束入射角

结合常梯度声线跟踪法，提出一种顾及波束入射角的常梯度声线跟踪水下定位算法。首先说明波束入射角的确定方法。由式(10)可以推算得到各层的传播时间为

将式(5)代入式(10)，并保留波束入射角θ₀可得

各层内的声波传播时间之和应该等于时间延迟观测量，即

式中，波束入射角θ₀可视为未知量，各层内的声速c₀, c₁, …, c_n+1和声速梯度g₀, g₁, …, g_n可根据换能器和基准点水深从声速剖面获得。为了获得准确的波束入射角θ₀的值，最简单的方法是直接在一定的区间范围内对θ₀进行搜索，依次代入式(14)并将最接近时间延迟观测量的搜索值作为波束入射角的最优值。但这种方法会带来巨大的计算量，尤其是在声速剖面采样点较多、搜索步长较小以及观测量较多的情况下。

3.2 解算法确定波束入射角

为了提高计算效率，本文利用解算方程的方法求解未知的波束入射角。为进一步简化计算方程，将式(5)代入式(13)，并保留Snell常数p，可得

根据各层内的声波传播时间之和应该等于时间延迟观测量，可得

式(17)是Snell常数p作为未知数的非线性方程。这里采用牛顿迭代法来对这个超越方程进行解算，牛顿迭代法可以表示为

式中，迭代次数j=0, 1, …，r。对于式(17)的解算，可具体为

通过牛顿迭代法解算超越方程的方式来确定波束入射角，可以有效提高计算效率。这种方法唯一的要求就是要有一个相对准确的初值以保证迭代收敛，这里推荐使用式(12)、式(5)来确定未知参数Snell常数p的初值。当确定了Snell常数p的准确值后，就可以由式(5)推算得到准确的波束入射角。实际上，也可以将波束入射角作为未知参数进行方程解算，但会使解算方程更复杂，增加不必要的计算量。

3.3 水声定位计算流程

以上便是在常梯度声线跟踪模型下推算波束入射角的核心方法。结合流程图 3对整个顾及波束入射角的常梯度声线跟踪水下定位算法的流程说明如下：

(1) 计算目标概略坐标。由声速剖面计算平均声速c，计算观测距离的近似值ρ_i=ct_i，代入式(2)—(3)解算得到目标点的概略坐标(x_t, y_t, z_t)。

(2) 确定波束入射角。由声速剖面、目标水深z_t、时间延迟t，可根据式(14)—(15)采用搜索法确定波束入射角θ₀，或根据式(16)—(21)采用解算法确定波束入射角θ₀。

(3) 声线跟踪计算。由波束入射角θ₀、声速剖面、时间延迟t，根据式(4)—(9)计算修正后的距离ρ_i。

(4) 计算目标坐标。由修正后的距离ρ_i代入式(2)—(3)解算目标坐标(x_t, y_t, z_t)。

(5) 重复步骤(2)—(4)，直至目标坐标的变化量小于阈值，阈值可根据相邻两次迭代结果的差值满足精度要求进行设定。

4 试验与分析

2019年7月，在我国南海海域开展了海洋大地测量基准与海洋导航新技术深海综合试验，在试验区域内以圆形航迹对海底基准点进行了连续观测(图 4)。水声数据由哈尔滨工程大学的深海高精度水声定位系统采集得到，测区内测量的声速剖面如图 5所示。分别采用平均声速法、常梯度声线跟踪-搜索法、常梯度声线跟踪-解算法进行水下目标的定位计算。表 1中统计了各方法计算得到的目标坐标、单位权中误差σ₀、迭代次数及计算时间；图 6为常梯度声线跟踪-搜索法在不同搜索波束入射角下计算得到的时间延迟偏差；图 7为两种常梯度声线跟踪定位方法(搜索法、解算法)计算的波束入射角以及利用坐标由式(12)计算的近似入射角。试验采用的计算机配置为Intel Core i3-4150处理器、8.00GB内存、Windows 7旗舰版(64位操作系统)。

对试验结果的讨论分析如下：

(1) 从解算得到的目标坐标分析。由表 1可以看出，常梯度声线跟踪(搜索法、解算法)的解算结果完全相同，表明两种方法都能实现准确的声线跟踪。当采用不同的平均声速进行定位解算时，尽管目标坐标的变化在平面方向不超过2m，但在深度上的变化超过了百米，因此平均声速法难以消除折射效应的影响，其精度与平均声速的取值直接相关。

(2) 从测距偏差的标准差分析。各方法测距残差的标准差相差不大，其中常梯度声线跟踪的单位权中误差最小，而平均声速法的单位权中误差均大于常梯度声线跟踪，且与平均声速的取值相关。需要特别说明的是，不论采用平均声速法或声线跟踪法，都是对距离观测值进行修正，再利用修正后的距离观测值进行最小二乘解算，因此单位权中误差反映的是定位解算的内符合精度。本次试验采用的是对称性的圆形航迹，波束入射角的变化相对较小(图 7)，这种情况下折射效应对观测距离的影响接近系统性的常值。这就解释了为何平均声速法与声线跟踪法的单位权中误差均为分米级，但深度方向的坐标却相差达到了数十米甚至数百米。

(3) 对常梯度声线跟踪定位算法进一步分析。不同波束入射角下的时间延迟偏差如图 6所示，通过搜索可以得到一条近似的“V”形曲线，曲线最小值点对应着搜索的最优波束入射角，也可以看出波束入射角与时间延迟偏差的关系是非线性的。图 7为常梯度声线跟踪计算得到的波束入射角，同时也给出了利用坐标反算得到的近似入射角，二者之间的偏差在2°以上，因此采用近似入射角必然使跟踪计算产生偏差。而对于本例，直接采用近似入射角进行声线跟踪定位解算会发生迭代不收敛的问题。

(4) 从定位计算效率分析。由表 1统计结果可知，平均声速法尽管迭代次数较多，但由于算法简单，所用的计算时间非常短；常梯度声线跟踪-搜索法理论严密且简单易行，但计算效率较低；常梯度声线跟踪-解算法通过方程解算的方式来提高计算效率，计算效率较搜索法提高了15倍以上。由于常梯度声线跟踪的理论模型最为严密，因此适用于高精度的定位解算。

综上所述，本文提出的两种顾及波束入射角的常梯度声线跟踪水下定位算法都能有效消除折射效应的影响。其中搜索法简单易行但计算量较大；解算法利用方程解算的形式，在保证计算精度的同时有效提高了计算效率；两种方法都是基于水下非差分定位模型，具有更好的普适性和一般性，可直接或间接用于类似的声线跟踪计算问题中。

5 结论

由海水声速的时空变化引起的声波折射效应，是影响声学测距精度的一个普遍问题。在声速剖面可准确测量的情况下，声线跟踪是消除折射效应的常用方法。但基于距离交会定位原理的水下声学定位系统，通常不会像超短基线系统、多波束测深系统等对波束入射角进行直接测量，也就导致了常用的声线跟踪方法无法直接应用于水下声学的距离交会定位。本文的意义在于采用更具一般性的非差分定位模型，结合常梯度声线跟踪，提出了一套完整的声线跟踪和定位计算流程。分别采用搜索法和解算法解决了波束入射角精确反算的关键技术问题，不仅提高了声线跟踪计算的精度和效率，也更清晰地梳理了声线跟踪与定位计算间的关系，对于更好地理解声速误差的影响特点也有一定的参考价值。

致谢: 特别感谢国家重点研发计划“海洋大地测量基准与海洋导航新技术”项目组及深海试验全体科研人员的合作与支持。