文章快速检索  
  高级检索
GNSS-A水下定位的动态非线性Gauss-Helmert模型及其抗差总体卡尔曼滤波算法
邝英才1,2, 吕志平1,3, 李林阳1, 王方超1, 许国昌3     
1. 信息工程大学地理空间信息学院, 河南 郑州 450001;
2. 火箭军士官学校作战保障系, 山东 青州 262500;
3. 哈尔滨工业大学(深圳)空间科学与应用技术研究院, 广东 深圳 518055
摘要:GNSS-声学组合式观测是确定海底控制点位置的重要手段, 但会受到声速不确定性、海面平台定位偏差等误差因素的干扰, 而基于误差传播定律的常规方法对各类误差的处理策略使得海底点坐标解算不准确。针对这一问题, 本文将声速测距误差非时变项设为待解参数, 在水下观测方程的系数矩阵中讨论声速测距误差时变项与换能器位置误差的影响, 构建了GNSS-声学水下定位的动态非线性高斯-赫尔默特(Gauss-Helmert, GH)模型, 并推导了该模型的总体卡尔曼滤波解。在此基础上, 进一步考虑扩展后的观测信息受到粗差污染的情况, 给出了模型的抗差处理方法及解算步骤。最后分别通过仿真试验和胶州湾海域实测试验进行了验证, 试验结果表明, 在不同深度或不同换能器位置误差大小的无粗差设定下, 本文方法解算精度及稳定性较常规方法均更高; 当观测信息含有粗差时, 模型的抗差滤波算法能更准确地识别及定位异常信息, 其三维点位精度明显更优, 解算效果达到最佳。
关键词GNSS-声学技术    海底控制点    声速测距误差    非线性GH模型    抗差估计    
Dynamic nolinear Gauss-Helmert model and its robust total Kalman filter algorithm for GNSS-acoustic underwater positioning
KUANG Yingcai1,2, Lü Zhiping1,3, LI Linyang1, WANG Fangchao1, XU Guochang3     
1. Institute of Surveying and Mapping, Information Engineering University, Zhengzhou 450001, China;
2. Department of Combat Support, Rocket Force NCO College, Qingzhou 262500, China;
3. Institute of Space Science and Applied Technology, Harbin Institute of Technology, Shenzhen 518055, China
Abstract: The GNSS-acoustic combined observing is an important means to determine the position of seafloor control points, but it will be interfered by error factors such as the uncertainty in sound velocity and the positioning deviation of the sea surface platform. However, the processing strategy of general method based on the error propagation law for various errors makes the seafloor point coordinate solution inaccurate. To solve the above problems, this paper sets the time-invariant term of sound velocity ranging as the parameter to be solved, and discusses the influence of time-varying error of sound velocity ranging and transducer position error in the coefficient matrix of underwater observation equation. Thus, the dynamic nonlinear Gauss-Helmert (GH) model for GNSS-acoustic underwater positioning is constructed, and the total Kalman filter solution of this method is derived. On this basis, taking into account the gross errors polluting of the observations, the robust method and solution steps of the new model are given. Finally, simulation experiments and a testing experiment in the sea area near Jiaozhou Bay are used to verify the performance of the new model. The results show that under conditions with no gross errors and either different water depths or different transducer position errors, the accuracy and stability of the proposed method are both higher than those of the general method. When the observations are polluted by gross errors, the robust filter algorithm of the new model can accurately identify and locate the abnormal information. The precision of its 3D point deviation results can be obviously optimized, and the solution performance is superior to that of the general method.
Key words: GNSS-A technology    seafloor control point    sound velocity ranging error    nonlinear GH model    robust estimation    

海洋大地控制网是经略海洋的基础保障,作为重要的起算基准问题,海底控制网点绝对坐标的可靠获取已成为海洋测绘领域的研究热点[1-2]。目前,采用GNSS-声学(GNSS-acoustic,GNSS-A)联合定位技术对海底基准进行布设和定位的方案应用广泛。该方案利用海面浮标或测量船的瞬时位置组成虚拟基线网,水下由换能器与海底控制点(应答器)建立声学应答联系,以获得全球坐标框架下海底基准的三维位置信息[3]

围绕这一热点问题,不少学者在走航模式、定位模型、参数估计方法等方面进行研究[4-6],究其本质,是为了削弱这个测量平差问题中固有的系统误差和粗差[7-8]。声速时空变化引起的测距误差是影响GNSS-A技术定位精度的主要因素,对其根源更直观的描述是由于声速的不确定性,包括背景声速剖面测量造成的非时变类误差和海洋内波、潮汐、日温等现象造成的周期性、随机性时变类误差[9-11]。目前具有代表性的解决方法是进行历元间作差和增加待估参数[12-14]。前者采用差分模型进行解算以消除长周期项误差,但该模型没有考虑短周期项及随机项误差,且要求参与差分的历元具有相似的观测结构,因此并不具备普适性;后者将所有观测时段的系统误差作为统一的待估参数参与平差解算,实际上只估计出了不随时间变化的声速测距常数项误差,定位模型有待进一步完善。而如船底换能器位置误差等则通常被归入随机模型中,基于误差传播定律得到协方差矩阵,再采用最小二乘算法(least squares,LS)进行解算[12]

上述误差处理策略虽然易于实现,但实质上可以被认为是一种间接或者近似的处理手段。更严密的思路应是将声速测距误差时变项、换能器位置误差等也一并保留在函数模型中进行处理,但随之而来的问题是,此时水下观测模型的系数矩阵将不再是常系数矩阵。当系数矩阵存在误差时,构建更合理的误差变量模型(errors-in-variables,EIV),通过整体最小二乘算法(total least squares,TLS)求解渐近无偏解[15-16]是较最小二乘算法更优的方法[17-19]。但以这一思路为基础求解GNSS-A水下定位问题尚待更细致深入的研究。

基于以上分析,本文尝试构建一种函数关系描述更为精细的海底控制点定位模型并给出其稳健算法,以保证海底点坐标的严密解算。首先从理论推导出发,剖析了基于误差传播定律的海底点定位方法存在的问题;其次,通过观测量的扩展与误差项的整合,构建了海底控制点定位的非线性GH模型,并得到了该模型的总体卡尔曼滤波解;然后,顾及扩展后的观测量中可能存在的粗差影响,将抗差估计引入参数求解中,给出了新模型的抗差总体卡尔曼滤波解法步骤;最后,通过仿真和实测试验对新模型的定位效果进行了验证。

1 基于误差传播理论的海底点定位方法

将采集的GNSS观测信息进行预先定位解算可以得到海面平台位置,结合姿态角信息、声学时延信息等,理想状态下计算海底控制点三维坐标的观测方程可以用式(1)表示[12]

(1)

式中,t表示声学量测时延;Cavg表示该海域的加权平均声速,可利用声速剖面信息修正计算得到;(xTd, yTd, zTd)表示船底声学换能器位置;(xTp, yTp, zTp)表示海底控制点(应答器)位置;δρdelay是由硬件延迟或网络通信引起的测距误差的合并表示,充分考虑测量仪器的改正参数,并进行相应的时间配准处理后[20],此类误差可忽略不计;δρvbδρvt分别表示由声速时空变化引起的常数项及时变类测距误差;ε表示偶然误差及其他未模型化的随机误差。

实际应用中,GNSS-A水下定位过程受到多种时变误差e*的影响,此时,式(1)应被改写为

(2)

式中,约定δρdelay已得到有效改正;et表示传播时延观测量中的随机误差;(exTd, eyTd, ezTd)表示基准转换后换能器坐标中的随机误差,主要来源于海面平台定位过程。海面平台的差分定位模式易受作业距离限制,精密单点定位模式依赖连续精密产品支持,因此即使测量时选择了抗多径的GNSS天线以避免主要的系统误差,海面基准位置仍会存在偏差。经过误差损失量级较小的基准转换后[21],换能器坐标也包含定位误差;eρvt表示声速测距误差中的时变部分,包括海洋日变化现象引起的长周期项误差、内波引起的短周期项误差和区域随机性误差。

针对上述误差,一种常规的解决思路是,增加统一的声速测距误差参数,将换能器位置误差归入随机模型中,利用误差传播定律[22],得到逐观测量对应的方差(暂不考虑粗差),进而得到模型对应的协方差矩阵,再利用最小二乘算法进行求解。

为推导简洁明了,暂设每个历元仅有一个测量值。假设在i历元对式(2)进行线性化处理

(3)

式中,(xTd0, i, yTd0, i, zTd0, i)和(xTp0, yTp0, zTp0)分别表示换能器和海底控制点的近似坐标;D0i=表示换能器与控制点间的近似几何距离;δX=[δxTp δyTp δzTp δρv]T表示模型待估参数,包括海底控制点三维坐标及声速测距误差δρvdi=D0i+δρv0表示线性化后的方程常数项;eTdi=[-(xTd0, ixTp0)/D0i -(yTd0, iyTp0)/D0i -(zTd0, izTp0)/D0ie(xTdi, yTdi, zTdi)表示线性化后的换能器位置误差。若把换能器位置误差归入随机模型中讨论,基于最小二乘原理,可以求解得到

(4)

式中,ATp= (aTp1)T … (aTpn)T T表示系数矩阵,在此模型中为常量矩阵;P=σ02·Σ-1表示模型对应的权矩阵,σ02为单位权方差,Σ为协方差矩阵,其对角线元素由各历元的换能器位置方差和声学时延方差经误差传播定律得到;l=[d1t1·Cavg1dntn·Cavgn]T表示模型的自由项。为表述简洁,后文将这种方法简称为EPL方法(error propagation law)。

剖析上述推导过程可以发现,换能器坐标及其误差间具有密切联系,声速测距误差也同时包含常数项和可变类误差,而EPL方法在对式(2)进行线性化处理时,虽然增设了统一的声速测距误差参数,但它并不能完全反映出声速不确定性引起的各类测距误差影响;且该方法将换能器位置误差归入随机模型中进行间接处理的策略也不够严密。

2 GNSS-A水下定位的动态非线性GH模型及其抗差解法 2.1 动态非线性GH模型的构建及其解法

为构建更为精细的GNSS-A水下定位模型,应当充分利用已知量、变量及误差项间的固有函数关系,对常规模型进行适当拓展。若要在拓展后的函数模型中处理各类误差项,则需要在系数矩阵中含有误差的情形下进行讨论,此时这个平差问题转变为不同精度条件下EIV模型的解算问题[23]。解算EIV模型需要进行线性化,而水下声学观测方程本身也是非线性的,因而此问题被扩展为了一个双非线性问题[24]。GH模型是一种通用平差模型,其对函数关系的描述更简洁灵活[18]。通过构建GNSS-A水下定位非线性GH模型,可以在不增加函数关系限制的情况下,迭代求解得到更具一致性与渐进无偏性的结果[25-26]

不同精度条件下的EIV模型函数关系及对应的随机模型可表示为[27]

(5)
(6)

式中,EA表示系数矩阵A相应的随机误差矩阵;ξ表示模型的待估参数向量;eL表示观测量矩阵L相应的随机误差向量;vec(·)符号表示矩阵的拉直运算;QLQvec(EA)分别表示eL和vec(EA)对应的协因数矩阵。EIV模型中系数矩阵含有随机误差,不再为固定的常数矩阵,其系数可视为一类伪观测量,此时将观测向量和误差向量进行相应扩展

(7)

式(5)和式(6)可改写为

(8)

式(8)即为更具一般性的非线性GH模型。在i历元时,将海底控制点观测方程表示为非线性GH模型形式

(9)

由于此模型中,系数矩阵含有多类误差,会随着待估误差项的变化而变化,因此需要进行多次迭代以达到不断更新的目的。假设在第j次迭代循环中,将上述模型在(Ye(j), ξ(j))处按泰勒级数线性化展开并化简,可得到

(10)

将偏导数写为矩阵形式

(11)

式中,则可进一步简化为

(12)

式(12)即为采用线性化后的GH模型描述的海底点定位问题。该模型中,待估参数δξ包括海底点坐标和声速测距常数项误差,总误差向量e包括换能器坐标随机误差、声速测距时变项误差和传播时延随机误差。

考虑到卡尔曼滤波在处理系统性误差时较最小二乘算法更为有效,当系统拓展至动态条件下,可在模型中加入状态方程[28],并利用总体卡尔曼滤波(total Kalman filter, TKF)方法进行求解[29-30]。基于此,将上述模型扩展为动态EIV模型,增加的状态方程可以写为

(13)

式中,δξ表示由上一历元信息传递得到的本历元推估值;W表示本历元预测值的残差。

式(12)、式(13)构成了海底点定位的动态EIV模型,而后采用总体卡尔曼滤波算法求解该模型。基于拉格朗日必要条件可以构建如下目标函数

(14)

式中,表示总误差向量e的协因数阵。将目标函数分别对各变量求偏导,使偏导数为0,经过等式变形与化简最终可得

(15)
(16)

式(15)和式(16)即为利用新模型解算海底控制点定位问题的总体卡尔曼滤波解。由上述推导

过程可知,与系数阵AeA(j)BeA(j)有关,而A(j)eABeA(j)又受到总误差估值向量的影响,因此不同于常规解算方法,本文模型的TKF解算是一个逐渐逼近无偏估值的迭代过程。但需注意的是,本次解出的待估参数需要在上次迭代结果的基础上进行更新,而总误差估值向量的结果是不具有迭代传递性的。

2.2 残差协因数阵的推导

在采用非线性GH模型处理海底点定位问题的基础上,进一步考虑模型中两类观测量受到粗差干扰的情况,为此,本文引入了抗差估计原理[31]以提高模型的稳健性,其核心是构造等价权阵削弱粗差对观测信息的影响。为保证解算方法在观测空间和结构空间均有较强的抵抗能力,采用标准化残差构造权因子函数[32]

(17)

式中,表示残差协因数阵中的对应元素,单位权方差σ0可由中位数函数计算得到[33]。因此问题的关键是求得新模型的残差协因数阵。可将式(16)改写为

(18)

M(j)=AeA(j)·δξ(j)w(j),联系式w(j)=f(Ye(j), ξ(j))+BeA(j)·e(j),可发现w(j)中的随机变量仅为模型观测量Y,因此根据协因数阵传播定律[22]可以得到

(19)

进而可以写出M(j)矩阵协因数阵传播形式

(20)

再结合式(18)及其协因数阵传播形式,可得到新模型的残差协因数阵表达式为

(21)

利用式(17)和式(21)得到本次迭代的标准化残差后,结合如IGG3函数等抗差权函数[34]的变形式,可以构造出下次迭代所需替换的等价协方差元素

(22)

式中,k0k1为经验阈值。将等价协方差矩阵更新后,代入式(15)和式(16)中进行下一次迭代解算,直至满足设定的终止条件。

2.3 海底点定位动态非线性GH模型的抗差解算步骤

与基于最小二乘算法的常规方法不同,求解本文模型涉及多个量的更新和迭代。给定总误差向量e及其协方差阵Q的初值,可按以下解算步骤多次迭代进行。

(1) 根据由传播时延随机误差、换能器位置误差、声速测距时变项误差构成的总误差向量初值e0,构造矩阵AeA(0)BeA(0)w(0),并结合初值Q0,采用式(15)和式(16)求解模型的TKF算法结果,再将其视为抗差迭代过程的初始值。

(2) 假设为第j次抗差迭代,根据式(21)求解此时的残差协方差矩阵,根据式(17)计算声学时延观测量和伪观测量对应的标准化残差,再结合抗差权函数的变形式得到本次迭代的等价协方差矩阵Q(j)

(3) 根据本次迭代的先验信息,更新矩阵AeA(j)BeA(j)w(j),再结合本次迭代更新后的等价协方差矩阵Q(j),采用式(15)和式(16)再次进行模型的抗差求解,得到本次迭代的参数估值向量和总误差估值向量

(4) 根据设定的终止条件(本文试验设定条件为,其中阈值常数ε取为10-8)进行判断。若不满足条件,则回到步骤(2)开始新一轮的迭代;若满足条件,则停止迭代,此时便得到了模型的抗差解。

根据上述描述,可以用图 1表示本文提出的海底控制点定位动态非线性GH模型的抗差TKF算法流程。

图 1 海底点定位动态非线性GH模型的抗差TKF算法 Fig. 1 Algorithm of robust TKF algorithm for dynamic nonlinear G H model of seafloor datum point positioning

3 试验分析

为了验证本文提出的非线性GH模型及其抗差TKF解法在海底控制点定位中的应用效果,首先模拟了GNSS-A技术的海上测量环境,设计了多组仿真试验,然后对比了不同条件下海底应答器位置解算精度,最后将本文方法应用于胶州湾海域实测数据处理中。

3.1 仿真试验与结果

模拟环境中,测量船在一定深度的海面进行水下数据采集,海底布设了一个应答器作为控制点。为使试验方案更加严密,分别模拟了深度小于100 m的浅海环境和深度为1000 m和3000 m的深海环境。假设测量船沿平行线航迹走航,在浅海区域单条航线长度为200 m,平行线间隔为50 m。在深海区域,深度为1000 m时,单条航线长度为4000 m,平行线间隔为1000 m;深度为3000 m时,单条航线长度为6000 m,平行线间隔为1500 m。以模拟海底点为坐标原点,模拟海底水平面为XOY面,基于右手定则建立当地坐标系。坐标系的Z轴垂直于海底指向天顶方向,X轴指向测量船走航轨迹的长边方向。航迹与应答器相对平面位置的示意图如图 2所示。

图 2 仿真走航轨迹与应答器相对位置 Fig. 2 The relative positions of the trajectories and the transponder

假设测量船速度为4节,水面模拟波高为2 m,周期为15 s的余弦波动。试验中测量船沿预设轨迹进行了多圈走航,以保证有效观测量的采集。测量船静止状态下,换能器吃水为3.5 m,其坐标由船载天线相位中心坐标经精确的姿态改正得到。声学信号量测系统误差按照文献[12]中的方法给定,测距随机误差大小服从均值为0,方差为5 cm的正态分布。采用分层等梯度声线跟踪算法[35]仿真声学信号传播时间作为水下观测量。声速剖面采用Munk经验模型构建[36],如图 3所示。

图 3 Munk理想声速剖面 Fig. 3 The sound velocity profile constructed by the Munk empirical model

首先在无粗差条件下,讨论常规解法和本文方法在海底控制点定位中的应用效果。根据采集深度不同及换能器坐标随机误差大小不同这两种限定条件,设计了如下两个仿真试验。

试验1是为了考察当换能器位置误差一定时,不同深度条件下两种方法的海底点解算效果。不同深度设置时,水下声信号的传播时间各异,声速误差的变化情况也不相同。试验模拟采集深度分别为50、100、1000和3000 m,在换能器坐标水平方向和垂直方向上分别加入方差为0.3、0.5 m的随机误差,采用的计算方案如下。

方案1:采用本文第1部分所述EPL方法解算。

方案2:基于严密的EIV模型理论,构建本文提出的非线性GH模型进行TKF算法解算。

采用海底模拟点的三维点位偏差结果作为定位精度指标,计算公式为

(23)

式中,(xTp, yTp, zTp)为试验计算结果,(xTp*, yTp*, zTp*)为海底模拟点真值坐标。每种方法均模拟进行200次,三维定位精度对比如图 4所示。另分别对每种方法200次解算结果的三维精度序列统计最大值(max)、最小值(min)、均方根(root mean square,RMS)和标准差(standard deviation,STD),见表 1,其中RMS和STD的计算公式为

图 4 不同深度条件下两种方案200次模拟三维误差结果 Fig. 4 The positioning results of the two schemes after 200 simulations under different depth conditions

表 1 不同深度下两种方案解算精度统计结果 Tab. 1 The statistical positioning results of the two schemes at different depths  m
深度 方案 max min RMS STD
50 方案1 0.268 0.045 0.113 0.047
方案2 0.034 0.025 0.029 0.002
100 方案1 0.420 0.054 0.177 0.073
方案2 0.040 0.031 0.035 0.002
1000 方案1 0.233 0.055 0.118 0.038
方案2 0.044 0.031 0.037 0.003
3000 方案1 0.261 0.056 0.125 0.047
方案2 0.041 0.034 0.039 0.002

(24)
(25)

分析图 4表 1可以发现,在试验设计尽可能规避了测量模式等影响因素的情况下,浅水或深海环境中的换能器位置误差、声学测距误差仍对定位结果造成了一定影响。方案1采用常规方法进行解算,在不同深度情况下均只能取得分米量级结果;方案2采用本文提出的非线性GH模型进行TKF解算,其定位结果均可达厘米量级,较方案1更优。由图 4中可以看出,方案1经200次模拟后的结果序列较方案2波动更大,精度及稳定性均较差。这一点从序列统计结果中也能体现,方案2解算结果序列的RMS值分别为方案1的25.6%、19.8%、31.1%和30.9%;STD值分别为方案1的4.0%、2.6%、6.6%和3.4%,在不同深度条件的模拟环境下,本文提出的新模型及TKF解算结果均非常稳定,较常规算法的定位效果提升显著。分析出现以上结果的原因是,常规算法虽然增加了声速系统误差参数,但其估计值并不能完全代表声速测距的常数项与时变项误差的影响;而本文提出的非线性GH模型考虑到海底点定位问题的本质是函数模型的系数矩阵包含误差,因此扩展了观测量,将换能器位置误差和难以参数化的声速测距时变项误差作为误差向量,采用迭代逼近的方式进行TKF参数估计,较常规算法更为严密与合理。

试验2是为了考察当采集深度一定时,不同换能器位置误差大小条件下两种方法的海底点解算效果。考虑到近岸差分结果较深远海域海面动态定位结果更稳定,试验环境模拟为深度3000 m的深海区域。在换能器坐标水平方向和垂直方向上加入随机误差,其方差大小分别为(0.3 m,0.5 m)、(0.5 m,0.8 m)、(0.8 m,1.0 m)和(1.0 m,3.0 m),依然采用试验1中的两种方案进行定位计算。每种方法模拟进行200次,三维定位精度对比如图 5所示。解算结果的三维精度序列统计见表 2

图 5 不同换能器误差大小条件下两种方案200次模拟三维误差结果 Fig. 5 The positioning results of the two schemes after 200 simulations with different transducer position errors

表 2 不同换能器误差大小下两种方案解算精度统计结果 Tab. 2 The statistical positioning results of the two schemes with different transducer position errors  m
换能器位置误差 方案 max min RMS STD
(0.3, 0.5) 方案1 0.270 0.054 0.122 0.046
方案2 0.043 0.033 0.038 0.002
(0.5, 0.8) 方案1 0.454 0.057 0.178 0.076
方案2 0.043 0.033 0.038 0.002
(0.8, 1.0) 方案1 0.491 0.062 0.216 0.097
方案2 0.043 0.033 0.038 0.002
(1.0, 3.0) 方案1 1.482 0.066 0.497 0.232
方案2 0.043 0.033 0.038 0.002

图 5表 2可以看出,对比方案2,方案1解算精度及稳定性均较差,当水平方向和垂直方向上换能器位置误差的方差大小增大至(1.0 m,3.0 m)时,方案1解算结果最大三维精度接近1.5 m,RMS值接近0.5 m。在不同的换能器位置误差设置下方案2均是更优的解算方案,结果序列的RMS值分别是方案1的31.0%、21.3%、17.4%和7.6%;STD值分别为方案1的3.9%、2.4%、1.8%和0.7%。当海面动态定位精度不高导致换能器坐标存在较大误差时,常规方法仅仅将换能器位置误差归入随机模型中进行近似处理,这样的策略不利于构建精细的定位模型,难以得到高精度的海底点定位结果;此时本文提出的新模型更具优势,其构造的总误差向量兼顾两类主要误差的影响,依靠更严密的TKF迭代算法进行参数估计,从而能保证可靠的解算效果。

其次,在以上试验结果基础上,对换能器坐标或声学时延观测量增加随机粗差,讨论在含有粗差的条件下,EPL方法和提出的新方法在海底控制点定位中的应用效果。设定水深为100 m,换能器坐标在水平方向和垂直方向上的随机误差方差分别为0.3 m和0.5 m。当模拟粗差指定出现时,其出现的历元和具体位置均是随机的,即换能器坐标和时延观测量不一定同时含有粗差,各方向上的坐标参数也不一定同时含有粗差,其大小为随机误差的10倍。时延观测量先验方差按照入射角余弦函数进行设置[37]。在粗差位置完全随机的条件下,分别设计了以下6种试验方案。

方案1:未加入粗差时,基于第1部分所述EPL方法解算。

方案2:未加入粗差时,构建非线性GH模型进行TKF算法解算。

方案3:加入粗差后,基于第1部分所述EPL方法解算。

方案4:加入粗差后,构建非线性GH模型进行TKF算法解算。

方案5:加入粗差后,基于第1部分所述EPL方法,进行抗差最小二乘解算。

方案6:加入粗差后,构建非线性GH模型,将TKF迭代结果作为初值,采用迭代抗差TKF算法解算。

方案5和方案6在构造等价协方差矩阵时,均采用如式(22)所示的抗差协方差函数,取经验阈值k0=1.5,k1=7.0。每种试验方案各仿真300次的三维偏差结果如图 6所示,各项指标的统计结果见表 3

图 6 各方案解算精度对比 Fig. 6 3D point deviations of each scheme

表 3 各方案解算精度统计结果 Tab. 3 The statistical results of each scheme  m
方案 max min RMS STD
方案1 0.423 0.008 0.169 0.080
方案2 0.041 0.030 0.035 0.002
方案3 2.841 0.270 1.223 0.487
方案4 1.007 0.398 0.629 0.094
方案5 1.609 0.141 0.715 0.282
方案6 0.290 0.012 0.077 0.041

图 6表 3结果,可以得出如下结论。

(1) 对比加入粗差前后的方案结果可知,方案1和方案2解算精度明显更优,说明当水下声学观测信息或换能器坐标信息受到粗差污染时,海底点定位精度将大幅度降低,因此有必要削弱其影响。而对比方案3和方案4结果可知,当存在粗差时,本文提出的海底点定位新模型的解算效果仍优于常规方法,这与前文得到的结论是一致的。分析原因是,扩展观测量类别和分类描述主要误差影响后,新模型对异常信息有更高的敏感度,采用TKF算法进行解算时,包含异常信息的总误差向量随着迭代不断更新,一定程度上削弱了异常信息影响,但由于TKF算法稳健性不佳,粗差干扰无法真正消除,因此方案4结果虽然相对同样存在粗差的方案3更优,但仍存在不小的偏差。

(2) 对比方案3、方案4结果与引入粗差处理策略的方案5、方案6结果可知,方案5结果的RMS值较方案3提升了41.6%,STD值提升了42.2%;方案6结果的RMS值较方案4提升了87.7%,STD值提升了56.0%,说明在含粗差条件下,将抗差估计引入常规方法和本文提出的非线性GH模型中可以提高算法的稳健性,是十分必要的。

(3) 对比方案5和方案6结果可知,在含粗差条件下,采用本文提出的新模型及其抗差解法进行海底点求解,较基于抗差估计的常规算法有更佳的定位效果,其三维点位偏差结果序列RMS值可提高至0.077 m,STD值可提高至0.041 m,分析原因是,由于模拟粗差在频次与位置上存在随机性,当观测信息受粗差污染特别是换能器坐标存在粗差时,常规方法将线性化后的换能器误差直接归入随机模型,并且忽视了声速测距误差的时变部分,不严密的模型误差处理策略使其无法准确识别和定位异常观测信息,进而难以构造正确的等价权削弱其影响。而方案6将系数矩阵中的换能器位置误差、声速测距时变误差等归为统一的误差项,采用非线性GH模型严密处理这几类误差因素,其抗差TKF算法兼顾结构空间和观测空间的抗干扰能力,通过独立构造每类观测信息的等价协方差矩阵,真正实现了对观测粗差及各类异常信息干扰的有效抵制。

3.2 实测数据处理应用

为进一步验证本文提出的模型及算法在实测数据处理中的应用效果,采用了胶州湾海域试验中的后处理数据进行解算。试验在海底布设了一个应答器,测量船搭载AeroAntenna公司AT1675系列GPS天线、POS-MV定位定姿系统、Valeport Mini系列声速剖面仪、HEUUSBL6000声学基准换能器及RBRduet潮位仪等传感器进行了约40 min的有效观测。其中声学基准换能器含有8个阵元,均匀分布在以换能器中心为圆心的圆上。阵元与船体坐标系三轴间的位置关系相对固定,因此其坐标均可由换能器中心坐标转换得到。对角线阵元的距离为0.26 m,同次采集可以得到8个水下声学时延观测值。图 7为本次试验中采集到的实测声速剖面。在部分运动状态相对稳定的观测历元,基于采集到的声学时延信息、水位深度信息、声速剖面信息等,建立了高精度的分层等梯度声线跟踪模型,经过深度改正、表层入射角精化等步骤,利用声线跟踪技术迭代计算出波束终点位置[35],在后续试验中,将此位置作为可靠的外部检核值。为保证后续解算的统一,以此坐标结果为原点,建立局部海底坐标系,约定坐标系的X轴指向子午线北,Z轴指向法线正向,Y轴指向正东方向,图 8为此坐标系下测量船走航轨迹与应答器位置的相对示意。

图 7 实测声速剖面 Fig. 7 The measured sound velocity profile

图 8 测量船航迹与应答器相对位置 Fig. 8 The trajectory of the survey ship and the relative position of the transponder

顾及仿真试验结论,分别采用基于最小二乘算法的EPL方法、引入抗差估计的EPL方法(EPL_R)、基于TKF解法的非线性GH模型和基于抗差TKF解法的非线性GH模型(RTKF)进行本次试验的数据解算。抗差协方差函数均采用式(22)形式,取经验阈值k0=1.5,k1=4.0。对于基于最小二乘算法的两种方案,将不同方案的解算结果与前述声线跟踪结果作差,并转换到本地坐标系下,计算其水平方向及三维点位偏差结果,作为这两种方案解算结果的外符合精度ΔEN、Δ3D;另计算迭代抗差后的单位权中误差和应答器坐标参数对应的权倒数,利用公式min=算各坐标平差参数对应的中误差,进而得到水平方向及三维点位的验后中误差mENm3D,作为这两种方案解算结果的内符合精度,见表 4;对于基于总体卡尔曼滤波算法的两种方案,将各历元结果与声线跟踪结果作差,得到逐历元的结果偏差。按式(24)计算偏差序列在水平方向和三维方向上的RMSEN、RMS3D,作为这两种方案解算结果的外符合精度;按式(25)计算偏差序列在水平方向和三维方向上的STDEN、STD3D,作为这两种方案的内符合精度,见表 5图 9为TKF算法和RTKF算法在各方向上逐历元的解算结果对比。

表 4 两种最小二乘算法解算结果外符合精度和内符合精度统计 Tab. 4 The statistical results of the outer and inner precision of each LS scheme  m
方法 ΔEN Δ3D mEN m3D
EPL 1.051 4.207 0.306 1.310
EPL_R 0.857 1.135 0.216 0.913

表 5 两种TKF算法解算结果外符合精度和内符合精度统计 Tab. 5 The statistical results of the outer and inner precision of each TKF scheme  m
方法 RMSEN RMS3D STDEN STD3D
TKF 0.271 0.401 0.146 0.226
RTKF 0.101 0.145 0.059 0.086

图 9 两种滤波算法逐历元解算结果对比 Fig. 9 The comparison of solution results of two filtering algorithms in each epoch

分析表 4表 5图 9结果可知,常规EPL方法难以保证高精度的估计结果,与精度相对可靠的声线跟踪结果对比,其三维方向上点位偏差超过4 m,验后中误差超过1.3 m;引入抗差估计后,三维点位偏差由4.207 m提高至1.135 m,验后中误差由1.310 m减小至0.913 m,解算效果有所改善,但其估值结果仍存在较大偏差。说明常规方法解算效果不佳不仅是因为观测信息受到粗差污染,更本质的原因是其模型构建的缺陷。实测环境下,海面平台定位偏差、复杂声速时空变化等因素会带来观测时延粗差、换能器位置误差及多类声速测距误差,统一的声速系统误差参数估值不能完全代表各类误差对结果的真实影响。此外,将换能器位置误差归入随机模型中的处理策略,也不利于抗差算法对异常信息的识别与定位,因而不精细的函数模型及随机模型导致结果出现了偏差。而基于新模型的估计结果较常规方法有较大提升,主要是因为构建新模型时剖析了此定位问题中的真实函数关系,将声速测距误差中的非时变与时变误差项分类进行建模,通过扩展常规模型为动态误差变量模型,使得几类误差的处理更为严密。分析图 9结果可知,基于抗差TKF算法的方案定位精度及稳定性最优,三维点位偏差序列的RMS值和STD值分别为0.145 m和0.086 m,当观测信息受到粗差污染时,该方案的估值结果始终优于TKF算法,说明RTKF算法抵抗粗差的能力更强,其合理的模型构建,也有利于针对性地削弱异常信息影响。以上结果可以证明,与常规方法相比,本文提出的新模型及其抗差算法在海底控制点定位解算中有更佳的应用效果。

4 结语

利用GNSS-A技术定位海底点是目前获得海底基准信息的有效方案之一。GNSS动态定位误差引起的换能器位置偏差不可避免,声速变化的不确定性也会引起常数项和时变类测距误差,这些误差均对海底控制点定位过程有直接影响,当同时存在观测粗差时,海底点位精度将受到严重干扰。针对这一问题,本文剖析了常规方法在求解海底点坐标时存在的缺陷,提出并推导了GNSS-A水下定位的动态非线性GH模型及其抗差TKF算法,通过仿真和实测试验得到了以下结论。

(1) 在海底控制点定位中,声速时空变化引起的测距误差与GNSS动态定位引起的换能器位置偏差是影响定位结果的主要因素。常规方法在函数模型中增加声速测距误差参数,在随机模型中加入线性化后的换能器误差,这种不严密的近似处理使得结果并不精确。在不同深度条件的试验中,该方法只能取得分米级结果,且随着换能器位置误差的增大,其结果的RMS值和STD值均显著变大。

(2) 本文提出构建更合理的动态非线性GH模型进行海底点定位,将观测量进行了扩展,设计了包含换能器位置误差和声速测距误差时变项的误差向量,通过构建动态误差变量模型并采用TKF算法迭代更新结果,在不同水深的仿真试验环境中,其定位精度及稳定性较常规方法有显著提升,当换能器位置误差增大时,也能取得更可靠的解算结果。

(3) 当观测信息受到粗差污染时,常规策略解算结果严重失真,引入抗差估计后,在一定程度上能抵抗粗差影响,但因其模型构建思路的缺陷,导致异常信息不能被准确识别与定位,从而不能更好地削弱模型误差和粗差的干扰。而本文提出的新模型及其抗差解法对函数关系的描述更为严密清晰,通过从构造空间和观测空间合理构造各类观测信息的等价协方差矩阵,能更好地提升稳健性,达到最佳的定位效果。在仿真试验中,其三维点位偏差结果序列的RMS值可达0.077 m,STD值可达0.041 m,在胶州湾海域实测试验中,其外符合精度可达0.145 m,内符合精度可达0.086 m。仿真试验和实测试验结果表明,本文方法应用于海底点位解算中是可行且可靠的。


参考文献
[1]
YANG Yuanxi, XU Tianhe, XUE Shuqiang. Progresses and prospects of marine geodetic datum and avigation in China[J]. Journal of Geodesy and Geoinformation Science, 2018(1): 16-24.
[2]
刘经南, 陈冠旭, 赵建虎, 等. 海洋时空基准网的进展与趋势[J]. 武汉大学学报(信息科学版), 2019, 44(1): 20-40.
LIU Jingnan, CHEN Guanxu, ZHAO Jianhu, et al. Development and trends of marine space-time frame network[J]. Geomatics and Information Science of Wuhan University, 2019, 44(1): 17-37.
[3]
SPIESS F N. Analysis of a possible sea floor strain measurement system[J]. Marine Geodesy, 1985, 9(4): 385-398. DOI:10.1080/15210608509379536
[4]
曾安敏, 杨元喜, 明锋, 等. 海底大地基准点圆走航模式定位模型及分析[J]. 测绘学报, 2021, 50(7): 939-952.
ZENG Anmin, YANG Yuanxi, MING Feng, et al. Positioning model and analysis of the sailing circle mode of seafloor geodetic datum points[J]. Acta Geodaetica et Cartographica Sinica, 2021, 50(7): 939-952. DOI:10.11947/j.AGCS.2021.20200529
[5]
孙文舟, 殷晓冬, 曾安敏, 等. 附加深度差和水平距离约束的深海控制点差分定位算法[J]. 测绘学报, 2019, 48(9): 1190-1196.
SUN Wenzhou, YIN Xiaodong, ZENG Anmin, et al. Differential positioning algorithm for deep-sea control points on constraint of depth difference and horizontal distance constraint[J]. Acta Geodaetica et Cartographica Sinica, 2019, 48(9): 1190-1196. DOI:10.11947/j.AGCS.2019.20180514
[6]
邝英才, 吕志平, 王方超, 等. GNSS/声学联合定位的自适应滤波算法[J]. 测绘学报, 2020, 49(7): 854-864.
KUANG Yingcai, LV Zhiping, WANG Fangchao, et al. The adaptive filtering algorithm of GNSS/acoustic joint positioning[J]. Acta Geodaetica et Cartographica Sinica, 2020, 49(7): 854-864. DOI:10.11947/j.AGCS.2020.20190393
[7]
YAMADA T, ANDO M, TADOKORO K, et al. Error evaluation in acoustic positioning of a single transponder for seafloor crustal deformation measurements[J]. Earth, Planets and Space, 2002, 54(9): 871-881. DOI:10.1186/BF03352435
[8]
FUJITA M, ISHIKAWA T, MOCHIZUKI M, et al. GPS/ acoustic seafloor geodetic observation: method of data analysis and its application[J]. Earth, Planets and Space, 2006, 58(3): 265-275. DOI:10.1186/BF03351923
[9]
SPIESS F N, CHADWELL C D, HILDEBRAND J A, et al. Precise GPS/acoustic positioning of seafloor reference points for tectonic studies[J]. Physics of the Earth and Planetary Interiors, 1998, 108(2): 101-112. DOI:10.1016/S0031-9201(98)00089-2
[10]
OSADA Y, FUJIMOTO H, MIURA S, et al. Estimation and correction for the effect of sound velocity variation on GPS/acoustic seafloor positioning: an experiment off Hawaii Island[J]. Earth, Planets and Space, 2003, 55(10): e17-e20. DOI:10.1186/BF03352464
[11]
SUN Wenzhou, YIN Xiaodong, ZENG Anmin. The relationship between propagation time and sound velocity profile for positioning seafloor reference points[J]. Marine Geodesy, 2019, 42(2): 186-200. DOI:10.1080/01490419.2019.1575938
[12]
XU Peiliang, ANDO M, TADOKORO K. Precise, three-dimensional seafloor geodetic deformation measurements using difference techniques[J]. Earth, Planets and Space, 2005, 57(9): 795-808. DOI:10.1186/BF03351859
[13]
ZHAO Jianhu, ZOU Yajing, ZHANG Hongmei, et al. A new method for absolute datum transfer in seafloor control network measurement[J]. Journal of Marine Science and Technology, 2016, 21(2): 216-226. DOI:10.1007/s00773-015-0344-z
[14]
SUN Wenzhou, YIN Xiaodong, BAO Jingyang, et al. Semi- parametric adjustment model methods for positioning of seafloor control point[J]. Journal of Geodesy and Geoinformation Science, 2020, 3(1): 85-92.
[15]
ADCOCK R J. Note on the method of least squares[J]. The Analyst, 1877, 4(6): 183. DOI:10.2307/2635777
[16]
FISHER G W, HUFFEL S V, VANDEWALLE J. The total least squares problem: computational aspects and analysis[J]. Mathematics of Computation, 1992, 59(200): 724. DOI:10.2307/2153088
[17]
SCHAFFRIN B, WIESER A. On weighted total least-squares adjustment for linear regression[J]. Journal of Geodesy, 2008, 82(7): 415-421. DOI:10.1007/s00190-007-0190-9
[18]
NEITZEL F. Generalization of total least-squares on example of unweighted and weighted 2D similarity transformation[J]. Journal of Geodesy, 2010, 84(12): 751-762. DOI:10.1007/s00190-010-0408-0
[19]
MAHBOUB V. On weighted total least-squares for geodetic transformations[J]. Journal of Geodesy, 2012, 86(5): 359-367. DOI:10.1007/s00190-011-0524-5
[20]
KUANG Yingcai, LU Zhiping, LI Linyang, et al. Robust constrained Kalman filter algorithm considering time registration for GNSS/acoustic joint positioning[J]. Applied Ocean Research, 2021, 107: 102435. DOI:10.1016/j.apor.2020.102435
[21]
CHADWELL C D. Shipboard towers for global positioning system antennas[J]. Ocean Engineering, 2003, 30(12): 1467-1487. DOI:10.1016/S0029-8018(02)00141-5
[22]
崔希璋, 於宗俦, 陶本藻, 等. 广义测量平差[M]. 2版. 武汉: 武汉大学出版社, 2009.
CUI Xizhang, YU Zongchou, TAO Benzao, et al. Generalized surveying adjustment[M]. 2nd ed. Wuhan: Wuhan University Press, 2009.
[23]
FANG X. Weighted total least squares solutions for applications in geodesy[D]. Hannover: Leibniz University, 2011.
[24]
胡川, 陈义. 非线性整体最小平差迭代算法[J]. 测绘学报, 2014, 43(7): 668-674.
HU Chuan, CHEN Yi. An iterative algorithm for nonlinear total least squares adjustment[J]. Acta Geodaetica et Cartographica Sinica, 2014, 43(7): 668-674.
[25]
CHANG Guobin. On least-squares solution to 3D similarity transformation problem under Gauss-Helmert model[J]. Journal of Geodesy, 2015, 89(6): 573-576. DOI:10.1007/s00190-015-0799-z
[26]
方兴, 曾文宪, 刘经南, 等. 基于非线性高斯-赫尔默特模型的混合整体最小二乘估计[J]. 测绘学报, 2016, 45(3): 291-296.
FANG Xing, ZENG Wenxian, LIU Jingnan, et al. Mixed LS-TLS estimation based on nonlinear Gauss-Helmert model[J]. Acta Geodaetica et Cartographica Sinica, 2016, 45(3): 291-296. DOI:10.11947/j.AGCS.2016.20150157
[27]
AMIRI-SIMKOOEI A, JAZAERI S. Weighted total least squares formulated by standard least squares theory[J]. Journal of Geodetic Science, 2012, 2(2): 113-124. DOI:10.2478/v10156-011-0036-5
[28]
SCHAFFRIN B, IZ H B. Towards total Kalman filtering for mobile mapping[C]//Proceedings of the 5th International Symposium on Mobile Mapping Technology. Padua, Italy: ISPRS, 2007, 270-275.
[29]
MAHBOUB V, SAADATSERESHT M, ARDALAN A A. A general weighted total Kalman filter algorithm with numerical evaluation[J]. Studia Geophysica et Geodaetica, 2017, 61(1): 19-34. DOI:10.1007/s11200-016-0815-7
[30]
余航, 王坚, 王乐洋, 等. 动态EIV模型及其总体卡尔曼滤波方法[J]. 测绘学报, 2018, 47(4): 480-489.
YU Hang, WANG Jian, WANG Leyang, et al. Total Kalman filter method of dynamic EIV model[J]. Acta Geodaetica et Cartographica Sinica, 2018, 47(4): 480-489. DOI:10.11947/j.AGCS.2018.20170098
[31]
周江文. 经典误差理论与抗差估计[J]. 测绘学报, 1989, 18(2): 115-120.
ZHOU Jiangwen. Classical theory of errors and robust estimation[J]. Acta Geodaetica et Cartographica Sinica, 1989, 18(2): 115-120.
[32]
王彬, 李建成, 高井祥, 等. 抗差加权整体最小二乘模型的牛顿-高斯算法[J]. 测绘学报, 2015, 44(6): 602-608.
WANG Bin, LI Jiancheng, GAO Jingxiang, et al. Newton- Gauss algorithm of robust weighted total least squares model[J]. Acta Geodaetica et Cartographica Sinica, 2015, 44(6): 602-608. DOI:10.11947/j.AGCS.2015.20130704
[33]
ROUSSEEUW P, LEROY A. Robust regression and outlier detection[M]. New York: John Wiley and Sons, 1987.
[34]
YANG Y X. Robust estimation for dependent observations[J]. Manuscripta Geodaetica, 1994, 19: 10-17.
[35]
CHADWELL C D, SWEENEY A D. Acoustic ray-trace equations for seafloor geodesy[J]. Marine Geodesy, 2010, 33(2/3): 164-186.
[36]
MUNK W H. Sound channel in an exponentially stratified ocean, with application to SOFAR[J]. The Journal of the Acoustical Society of America, 1974, 55(2): 220-226. DOI:10.1121/1.1914492
[37]
ZHAO Shuang, WANG Zhenjie, HE Kaifei, et al. Investigation on underwater positioning stochastic model based on acoustic ray incidence angle[J]. Applied Ocean Research, 2018, 77: 69-77.
http://dx.doi.org/10.11947/j.AGCS.2023.20210467
中国科学技术协会主管、中国测绘地理信息学会主办。
0

文章信息

邝英才,吕志平,李林阳,王方超,许国昌
KUANG Yingcai, Lü Zhiping, LI Linyang, WANG Fangchao, XU Guochang
GNSS-A水下定位的动态非线性Gauss-Helmert模型及其抗差总体卡尔曼滤波算法
Dynamic nolinear Gauss-Helmert model and its robust total Kalman filter algorithm for GNSS-acoustic underwater positioning
测绘学报,2023,52(4):559-570
Acta Geodaetica et Cartographica Sinica, 2023, 52(4): 559-570
http://dx.doi.org/10.11947/j.AGCS.2023.20210467

文章历史

收稿日期:2021-08-19
修回日期:2022-04-28

相关文章

工作空间