海底控制点定位的半参数平差模型法

2019-02-13孙文舟殷晓冬暴景阳曾安敏

测绘学报 2019年1期

孙文舟,殷晓冬,暴景阳,曾安敏

1. 海军大连舰艇学院军事海洋与测绘系，辽宁大连 116018； 2. 武汉大学测绘学院，湖北武汉 430079； 3. 地理信息工程国家重点实验室，陕西西安 710054

高精度确定海底控制点的三维坐标是建立海底大地控制网的关键环节[1]。船载GNSS结合水下声学测距是实施水下声标(控制点)交会定位的有效手段[2-10]，影响这种方法定位精度最主要的因素是声速在时间和空间上复杂变化所引起的测距误差。采用圆走航的方式可以保证水面观测点相对于海底控制点对称分布[11]，对提高控制点水平坐标的精确度有重要帮助，然而垂直方向的解依然较差且不稳定。

为了提高海底控制点垂直解的精确度，国内外学者进行了不同的尝试。文献[12]提出了利用观测历元之间差分的方法，这种方法可以有效地消除声速误差长周期项的影响，但是短周期项依然存在。文献[13]利用传播时间建立了声速误差模型，提出了在500 m水深情况下，不依赖声速剖面的水下静态目标三维坐标快速求解方法。文献[14]认为声速误差项为常数，通过南海的试验验证了该方法相比于传统方法的有效性。近年来，压力传感器的测深精度不断提高，又有学者提出了深度约束下的海底控制点坐标确定方法[15-17]，进一步提高了三维坐标解的稳健性和精度。

上述方法都是对声速变化引起的测距误差中可以参数化的部分进行建模改进。根据文献[2,18]对水下声学测距误差的研究，由内波引起的测距误差形态复杂，难以用少量的参数表达，而引入参数过多又往往会导致法方程的病态甚至秩亏。因此，本文尝试利用半参数平差模型进行控制点坐标的解算，以期提高垂直解的精确度。

1 海底控制点坐标半参数平差模型解算方法

1.1 声速变化引起的测距误差分析

影响声速变化的主要因素是太阳辐射和海水运动[19]。声速剖面的变化包括日变化、季节变化以及附加之上的随机扰动，声速剖面的季节变化主要与太阳直射点在南北回归线之间周期性的运动有关，而声速剖面日变化主要与太阳高度角的日变化有关，由外力扰动引起的海洋内波造成水体微团垂向运动，使得上下层海水混合引起温度的变化从而导致海水声速变化[20]。文献[2]对声学测距误差的研究表明，由声速变化引起的测距系统误差项包含两部分，长周期项与短周期项，其中长周期误差项与潮汐的影响有关，而短周期误差项则主要与海洋内波现象有关。之后文献[18]在夏威夷岛附近进行了相关试验，再一次验证了这个结论，并且认为声速主要受温度变化影响，盐度变化对其影响微乎其微[18]。

为进一步分析文献[2,18]的试验结论，研究声速对测距误差的影响规律，本文进行了如下试验，首先利用文献[21]的标准方程构建了背景声速剖面[21-22]

(1)

式中，C1是声道轴处的声速；z1是声道轴的深度；ε=7.4×10-3是扰动系数；B=1300为声道厚度尺度；z是水层深度。

利用表层温度的变化规律以及热传导方程构建声速时空变化的部分[23]

Δc=Δc0e-βzsinωt-βz-φ

(2)

式(2)可表示声速剖面在不同时间不同深度的声速日变化值，通过式(1)和式(2)可构造一个稳定变化的声速剖面。

在海底控制点坐标解算的过程中往往采用某一固定的声速剖面，用上述方法构建的声速剖面簇的平均值作为固定声速剖面，计算在水深为3000 m，声速采样间隔为5 m的条件下，声速变化所引起的测距误差。为了不引入入射角偏差而导致的误差，仿真试验采用垂直波束，设表层温度日变化的振幅为0.2℃，每小时采样一次，共采样24次，表层海水声速变化与测距误差的计算结果如图1所示。

图1 表面声速变化与测距误差曲线Fig.1 Surface sound velocity variation and ranging error curve

温度的日变化只影响表层海水声速，图1(a)是对表层海水等间隔取样得到的海水声速平均值在时间上的变化，图1(b)是各采样时间所对应的测距误差的变化。从图中可以看出，测距误差的变化规律与声速的变化规律基本相似，周期相同但相位存在差异。

文献[2,18]的试验结果显示，测距误差长周期项的变化与潮汐的变化规律相近。这种现象是由于潮波的振动产生了日潮内波或半日潮内波，在潮内波的影响下海水产生垂向运动，从而导致了上层温度较高的海水与下层温度较低的海水进行热量交换，引起了海水温度的变化，温度变化是引起海水声速变化最主要的因素，其周期与潮内波的周期相同。测距误差短周期项的变化与海洋内波有关，与潮内波的作用机理相同，海洋内波造成了海水垂向运动，引起了海水温度的变化进而导致声速的变化。相比于声速长周期项变化，短周期项变化较为复杂，因此与其对应的测距误差短周期项随时间变化同样具有复杂的形式，难以进行参数化建模。

1.2 圆走航距离交会定位模型

水下控制点定位的观测方程通常表示为

ρi=fXi,Xo+δρdi+δρvi+εi

(3)

式中，ρi是第i时刻换能器到海底应答器的距离观测值；f(Xi,Xo)是两者之间的直线距离，Xi是第i时刻换能器的位置，Xo是海底应答器的位置；δρdi是由应答器电路延迟引起的系统性误差，是可计算的改正量；δρvi是声速变化及声线弯曲引起的系统性误差，为主要误差；εi是第i时刻测量的随机噪声。

若δρdi已通过外部设备进行改正，则式(3)线性化展开得[14]

(4)

(5)

若进行了n次观测，可将式(5)写成矩阵形式

V=BX+δρv-L

(6)

由上文分析可知，由声速变化引起的测距误差变化复杂，尤其是测距误差的短周期项，难以对其进行参数化的建模，若将每次观测的系统误差δρv都作为未知参数，则此时未知参数的个数为n+3，而观测方程的个数为n，因此式(6)无法得到未知参数的唯一解。本文采用半参数模型进行控制点三维坐标的解算。

1.3 控制点坐标解算的半参数模型

根据式(4)，可写出如下矩阵形式的观测方程

L=Bx+S+Δ

(7)

式中，S=[δpv1δpv2…δpvn]T，是由声速变化引起的测距系统误差的n维未知向量。式(7)即为海底控制点坐标解算的半参数模型[24-25]。

式(7)的误差方程是

(8)

由于未知参数多于方程的个数，因此无法直接利用最小二乘原理进行求解，需要选择新的平差准则[22]

VTPV+αSTRS=min

(9)

式中，α是光滑因子，对两参数V和S起平滑作用；R是一个给定的正定矩阵，本文中选定R=GTG，其中

同时，若令所有历元声速变化引起的测距误差和为0，则又可设置一个限制条件

eTS=0

(10)

利用式(8)—式(10)可构造拉格朗日函数

2γeTS

(11)

根据所构造的拉格朗日函数即可推导出式(8)—(9)的解，这里直接给出推导结果

(12)

将计算得到的坐标改正数加上所设定坐标初始值即可获得海底控制点的三维绝对坐标。

(13)

式中，δpvl是测距误差的长周期项；δpv0是长周期项测距误差的振幅；T是测距误差的长周期项的周期；φ是初始相位。对长周期项参数化处理后，短周期项仍然认为不可参数化，同样根据式(12)计算坐标改正数，加上坐标初始值获得控制点三维绝对坐标。

2 仿真试验分析

为了验证半参数平差模型解算海底控制点坐标的有效性，本文进行了如下仿真试验：水深设置为3000 m，应答器布放于海底，测量船以应答器为中心，半径为3000 m的圆航迹航行，速度为6节(约3 m/s)，每10 s实施一次水声测距。解算的坐标系采用笛卡儿坐标系，应答器为原点，x轴指向东，y轴指向北，z轴指向天顶。测距误差的仿真采用与文献[12]相同的方法，表示为

(14)

式中，误差共包含4项，第1项和第2项分别是声速变化引起的短周期项和长周期项，第3项是测区相关性误差，第4项为随机性误差。参照文献[2]在北太平洋的试验和文献[18]在夏威夷岛附近的试验，测距误差中的长周期项与潮汐变化规律近似，周期为12 h，振幅约20 cm，由内波引起的测距误差短周期项的周期从几十分钟到几个小时不等，振幅约为12 cm，因此令c1=12 cm、c2=20 cm、c3=2 cm，Tw=20 min，随机误差满足方差为5 cm的高斯分布，t0是初始时刻，t是任意时刻，x是t时刻测量船的三维坐标，x′是海底应答器的三维坐标。除此之外，设置潮汐的周期为12 h，振幅为5 m，海面波浪的周期为20 s，振幅为2 m，测量船水平方向的定位精度为5 cm，垂直方向定位精度为10 cm。

试验中将本文的方法与其他两种方法进行对比，方法1是直接对观测值利用最小二乘进行解算，方法2是采用文献[12]所提出的差分算法进行解算，本文方法为方法3。

第1次仿真试验的采样数为4320，采样时长刚好为测距误差长周期项的一个周期，试验中分别考虑存在内波和不存在内波两种情况，表1给出了不同方法解算所得控制点坐标偏差的统计参数。

表1 不同方法解算所得控制点坐标偏差的统计参数

对比水平方向和垂直方向坐标偏差的精确度(MSE)可以发现，声速误差对控制点垂直解的精确度影响较大，而对水平解的精确度影响小，这与圆走航测距定位消除水平方向误差的设计相吻合，也被国内外研究所证实。对比两个试验采用差分算法计算的垂直坐标偏差的精确度，可验证文献[12]的结论，即差分可以消除测距误差中长周期项的影响，但是短周期项的影响依旧存在。当存在内波时，水平方向依旧可以达到较高精确度，但是垂直方向解的精确度会明显变差。将3种方法进行横向比较，当不存在内波时，差分算法是3种方法中垂直解精确度最高的，高于半参数模型和最小二乘法，而在实际的海洋环境中，更普遍的情况是存在海洋内波，此时，差分算法垂直解的精确度降低，而最小二乘法与半参数模型解的精确度较高。

第2次试验的采样数为720，采样时长为测距误差长周期项周期的1/6，试验中同样分别考虑存在内波与不存在内波两种情况，不同方法解算控制点坐标偏差的统计参数见表2。

表2 不同方法解算所得控制点坐标偏差的统计参数

将试验2与试验1垂直坐标偏差的精确度进行对比可得，无论是否存在内波，当观测时长不是测距误差长周期项的整数倍时，最小二乘法垂直解的精确度都会降低。差分方法的结论与前一试验一致，整体精确度降低，这主要与观测数据的数量减少有关。半参数模型的垂直解则依旧可以得到较高的精确度。图2—图5中测距误差是指采用文献[12]的方法仿真得到的真实测距系统误差，而误差估值是指半参数模型在不同试验条件下对每次测距误差的估值，可以看出半参数模型可以较准确地描述测距误差的变化规律，对比图2和图4、图3和图5，可以发现在内波存在试验中，半参数模型对测距误差估值的结果明显发散，精度明显下降，表现为图3中测距误差与误差估值的中误差分别为0.17和0.22 cm；图5中测距误差与误差估值的中误差分别为0.10和0.18 cm。对比表2中半参数模型在存在和不存在内波两种条件下垂直坐标偏差的精确度，可以得出内波存在时垂直解误差的精确度由6.43增加到25.53，证明内波会影响半参数模型解算的精确度，但依旧小于最小二乘法的253.37和差分算法的389.20。从表2中还可以看出，两种条件下，虽然最小二乘法垂直解的精度高于半参数模型，但由于其准确度较低，所以其精确度远小于半参数模型。第2次试验对测距误差长周期项进行了参数化处理，结果显示，无内波时对振幅的估计为20.02 cm，而当内波存在时，对振幅的估计为15.92 cm，相对于所设定的值20 cm，估值的准确性明显降低，对应垂直方向坐标解的误差也增大，这个结果表明，当测距误差同时包含参数项和非参数项时，非参数项误差会影响对参数项估值的准确性。

图2 半参数模型对测距误差的估值(无内波)Fig.2 Ranging error estimated by semi-parametric model (no internal waves)

图3 半参数模型对测距误差的估值(有内波) Fig.3 Ranging error estimated by semi-parametric model (with internal waves)

图4 半参数模型对测距误差的估值(无内波)Fig.4 Ranging error estimated by semi-parametric model (no internal waves)

图5 半参数模型对测距误差的估值(有内波) Fig.5 Ranging error estimated by semi-parametric model (with internal waves)

在控制点三维坐标解算的过程中，最小二乘法和差分算法要想达到最优解，都需要较为苛刻的条件，最小二乘法要求观测的时长为测距误差长周期项的整数倍，而差分算法要求不存在短周期项的测距误差，前者显然实际工程应用价值较低，因为测量船连续圆走航12 h甚至24 h难度大，而后者所需要不存在内波现象的海洋环境条件更加难以满足。因此半参数模型的优势则明显突出出来，它不但可以取得较为精确的坐标解算值，而且对客观环境的依赖度相对较低。