三维坐标转换的通用整体最小二乘算法

方 兴,曾文宪,刘经南,2,姚宜斌

1.武汉大学测绘学院,湖北 武汉 430079;2.武汉大学卫星导航定位技术研究中心,湖北 武汉 430079

三维坐标转换的通用整体最小二乘算法

方 兴1,曾文宪1,刘经南1,2,姚宜斌1

1.武汉大学测绘学院,湖北 武汉 430079;2.武汉大学卫星导航定位技术研究中心,湖北 武汉 430079

三维坐标转换模型属于非线性EIV(errors-in-variables)模型,现有整体最小二乘算法均设定了某些特殊的假设条件,如仅适用于小角度或者属于非统计意义上的数值解,并且不能用于结构性的系数矩阵等,算法适用性受到极大限制。本文提出三维坐标转换模型的通用加权整体最小二乘算法,该算法适用于任意旋转角度以及一般性权矩阵的三维坐标转换,并且将结构性系数矩阵,同时包含随机和非随机元素的系数矩阵等情况纳入统一的算法。实例计算表明,本文提出的算法具有通用性,适用于实际应用中的各类三维坐标转换模型。

整体最小二乘估计;三维坐标转换;变量含误差模型;非线性算法

1 引 言

三维坐标转换模型是大地测量、摄影测量、地图投影、计算机视觉等专业领域常用的变量含误差(errors-in-variables,EIV)模型之一,模型描述了源坐标系统[X Y Z]至目标坐标系统[x y z]的函数映射关系[1]

式中,n表示坐标点个数;1n表示n×1的元素为1的向量;

其中,M表示坐标旋转矩阵;α1、α2和α3表示坐标系绕x、y和z轴旋转角度参数;M1、M2和M3表示相应的旋转矩阵;[ΔxΔyΔz]T表示平移参数向量;μ表示尺度参数;式(1)为七参数三维坐标转换模型。

坐标转换模型的整体最小二乘估计(total least squares method,TLS)具有渐进无偏性[2]。可以看到,模型(1)的系数矩阵为[A 1n],包含随机量与非随机量,并且向量化后的系数矩阵中的随机量呈结构性特征[3];其次,旋转参数[α1α2α3]与系数矩阵呈复杂的非线性关系。因此,三维坐标转换模型无法直接采用TLS算法求解。文献[4]将式(1)转换为线性化的Gauss-Helmert模型(Gauss-Helmert model,GHM)求解,由于算法的迭代更新项存在问题[5],得到的并不是最优的TLS解。文献[6]提出了一种仅适用于小旋转角度的TLS解。文献[7]基于Gauss-Markov模型和EIV模型的正交回归解的等价性提出了加权TLS算法,但对权阵有严格的限制条件。文献[1,8—9]采用迭代GHM方法推导了一般权阵条件下的加权TLS算法,迭代GHM方法尽管被认为是解决非线性问题的有效方法,但到目前为止,该方法并没有从理论上予以证明[10]。近年来,国内大地测量领域对三维坐标转换模型的TLS算法展开了研究[11-13]。文献[14]讨论了小角度、特殊权矩阵条件下的基于奇异值分解的TLS算法。文献[15]在单位权阵条件下,将旋转矩阵全部9个元素视为待估参数,构造了附有限制条件的间接平差模型,采用奇异值分解方法,得到了任意旋转角度的TLS解法。由于基于奇异值分解的TLS算法本质上属于数值解法,无法得到统计意义上的估计值,在测量中的应用受到了限制[16-17]。文献[18]通过将式(1)转换为附有限制条件的平差模型导出了特殊权矩阵条件下的加权TLS算法。

可以看到,现有三维坐标转换模型的TLS算法均受到特殊条件限制,实际应用范围有限,如摄影测量中物方坐标系到像方坐标系的转换是大角度三维坐标系统之间的转换,不同大地基准下点的坐标测量技术手段可能不同,从而坐标点的精度不等且相关,因此,有必要研究一般情况下的TLS算法。本文通过提取系数矩阵的随机量作为待估计参数,研究了三维坐标转换模型的通用加权整体最小二乘算法,该算法适用于任意旋转角度以及一般性的权矩阵,将结构性系数矩阵,同时包含随机和非随机元素的系数矩阵等情况,纳入统一的坐标转换模型算法,能够很好地满足现代大地测量领域对高精度三维坐标转换模型估计的需要。

2 三维坐标转换的加权整体最小二乘算法

为采用统一的方法求解结构性系数矩阵的EIV模型,将模型(1)中系数矩阵的随机量与坐标转换参数一并作为待估计参数,则三维坐标转换模型可描述为

式中,¯A表示将系数矩阵的随机量(即源坐标)作为模型的待估计参数所构成的矩阵。

经过以上变换,坐标转换的EIV模型可表示为非线性GM模型形式

在整体最小二乘准则下,进一步将式(3)的求解转化为经典的非线性模型的最优化问题。令式(3)中观测值l的一般性权阵为W,不限定观测数据(源坐标和目标坐标值)的精度或者相关性,式(3)的整体最小二乘最优估计可表示为以下最优化问题

导出式(4)的梯度即一阶偏导项

式中,各项梯度公式分别为

导出三维坐标转换模型的梯度公式(5)—式(10)后,根据非线性最优化理论,可以采用多种方法求解模型参数的最佳估值。如根据梯度公式,在参数的近似值处采用梯度下降法,即在反梯度方向搜索模型的最优估计值,或者采用高斯牛顿法求解等。这些方法均为非线性估计的基本方法,具体请参阅文献[19]。采用传统的大地测量求解非线性Gauss-Markov模型的算法同样可以求解模型(3)[20]。运用各种迭代算法求解坐标转换参数时,参数初始值可选取。以拟牛顿法为例,列出三维坐标转换模型式(2)的整体最小二乘算法的计算步骤如下:

(2)迭代计算拟牛顿法中的海森(Hessian)矩阵的近似表达式量,直到前后两次参数估计之差Δβi小于设定的正微小量

(3)根据以下Broyden-Fletcher-Goldfarb-Shanno(BFGS)算法公式,得到三维坐标转换模型的整体最小二乘估计结果

3 实 例

笔者选用文献[7]中大角度三维坐标转换模型实例数据,采用本文提出的算法进行计算。坐标转换数据来源于两个表面模型坐标转换的实测观测值,共包含4个公共点,相应的源坐标值和目标坐标值见表1,模型估计结果见表2。

表1 公共点三维坐标数据[7]Tab.1 Public points 3D coordinates

表2 三维坐标转换参数估计结果Tab.2 Results of the 3D datum transformation

当坐标数据等权即权阵为单位阵(W＝I24)的情况下,表2的第2列和第3列分别为转换参数的最小二乘估计值和整体最小二乘估计值,第4列为设定坐标数据非等权情况下(W＝I6⊗diag (1,4,6.25,16))的整体最小二乘估计值,表中的估计结果与文献[1,7]一致。对比第2列和第3列数据可以看到,整体最小二乘估计结果由于同时顾及了目标坐标和源坐标的随机误差,残差平方和远小于最小二乘的残差平方和,与文献[21—23]实例中坐标转换模型的最小二乘和整体最小二乘残差平方和对比结论一致。同时,表中旋转参数的最小二乘和整体最小二乘估计结果相等,与文献[7]证明的当观测数据的权阵为一定特殊形式时,三维坐标转换模型旋转参数的最小二乘和整体最小二乘估计结果相等的结论一致。第3列和第4列数据表明,权阵对整体最小二乘估计结果影响显著,若没有采用正确的权阵,例如没有考虑坐标数据不等精度、相关性,或者权阵没有考虑到系数矩阵的结构性特点等,都会导致整体最小二乘的估计结果产生偏差。

4 结 论

三维坐标转换模型是大地测量、摄影测量等众多专业领域常用的基本EIV模型之一,现代大地测量中高精度、三维动态大地坐标系统的建立对三维坐标转换模型的估计算法提出了更高的要求。现有TLS算法均受到一定限制条件的制约,如仅适用于小角度坐标转换模型、属于非统计意义上的数值解以及需要对坐标转换模型的结构性系数矩阵的权阵进行特殊处理等。本文通过将模型中的源坐标与坐标转换参数一并作为待估计参数,推导了三维坐标转换模型的TLS通用算法,该算法适用于任意旋转角的坐标转换以及一般性的权矩阵,并可直接应用于坐标转换模型结构性系数矩阵的情况。本文通过大角度三维坐标转换的实例对该算法进行了说明和验证,结果表明新算法突破了现有算法的限制条件,能有效地应用于坐标转换模型的一般情况,其通用性有利于应用中的编程实现。

[1] FANG Xing.Weighted Total Least Squares Solution for Application in Geodesy[D].Hanover:Leibniz University of Hanover,2011.

[2] HUFFEL S V,VANDEWALLE J.The Total Least-squares Problem:Computational Aspects and Analysis[M].Philadelphia:Society for Industrial and Applied Mathematics,1991.

[3] FANG Xing.A Structured and Constrained Total Leastsquares Solution with Cross-covariances[J].Studia Geophysica et Geodaetica,2014,58(1):1-16.

[4] BLEICH P,ILLNER M.Strenge Lösung der Räumlichen Koordinatentransformation Durch Iiterative Berechnung[J].Allgemeine Vermessungs-Nachrichten,1989,96(4):133-144.

[5] POPE A.Some Pitfalls to be Avoided in the Iterative Adjustment of Nonlinear Problems[C]∥Proceedings of the 38th Annual Meeting of American Society Photogrammetry.Washington DC:[s.n.],1972:449-473.

[6] ACAR A,ÖYLÜDEMIR MT,AKYILM AZ O,et al.Deformation Analysis with Total Least Squares[J].Natural Hazards and Earth System Sciences,2006,6: 663-669.

[7] FELUS F,BURTCH R.On Symmetrical Three-dimensional Datum Conversion[J].GPS Solutions,2009,13 (1):65-74.

[8] AKYILMAZ O.Solution of the Heteroscedastic Datum Transformation Problem[R].Munich:International Association of Geodesy,2012.

[9] LU Jue,CHEN Yi,FANG Xing,et al.Performing 3-D Similarity Transformation Using the Weighted Total Least-squares Method[R].Munich:International Association of Geodesy,2012.

[10] 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.

[11] LI Bofeng,SHEN Yunzhong,LI Weixiao.The Seamless Model for Three-dimensional Datum Transformation[J].Science China:Earth Science,2012,55(12):2099-2108.

[12] LI Bofeng,SHEN Yunzhong,ZHANG Xingfu,et al.Seamless Multivariate Affine Error-in-variables Transformation and Its Application to Map Rectification[J].International Journal of Information Science,2013,27(8):1572-1592.[13] CHEN Yi,LU Jue.Performing 3D Similarity Transformation by Robust Total Least Squares[J].Acta Geodaetica et Cartographica Sinica,2012,41(5):715-722.(陈义,陆珏.以三维坐标转换为例解算稳健总体最小二乘方法[J].测绘学报,2012,41(5):715-722.)

[14] LU Jue,CHEN Yi,ZHENG Bo.Applying Total Least Squares to the Three-dimensional Datum Transformation [J].Journal of Geodesy and Geodynamics,2008,28(5): 77-81.(陆珏,陈义,郑波.总体最小二乘方法在三维坐标转换中的应用[J].大地测量与地球动力学,2008,28(5): 77-81.)

[15] XU Chaoqian,YAO Yibin,XIONG Siting,et al.3 D Rectangular Coordinate Transformation Adapted to Arbitrary Rotation Angle Based on Total Least-SquaresRegression[J].Journal of Geomatics,2010,35(5):46-48.(许超钤,姚宜斌,熊思婷,等.三维任意旋转角度坐标转换的整体最小二乘回归解法[J].测绘信息与工程, 2010,35(5):46-48.)

[16] XU Peiliang,LIU Jingnan,SHI Chuang.Total Least Squares Adjustment in Partial Errors-in-variables Models:Algorithm and Statistical Analysis[J].Journal of Geodesy,2012,86(8): 661-675.

[17] FANG Xing.Weighted Total Least Squares:Necessary and Sufficient Conditions,Fixed and Random Parameters[J].Journal of Geodesy,2013,87(8):733-749.

[18] GE Xuming,WU Jicang.Iterative Method of Weight Constraint Total Least Squares for Three Dimensional Datum Transformation[J].Science of Wuhan University:Geomatics and Information,2012,37(2):178-182.(葛旭明,伍吉仓.三维基准转换的约束加权混合整体最小二乘的迭代解法[J].武汉大学学报:信息科学版,2012,37(2):178-182.)

[19] NOCEDAL J,WRIGHT S.Numerical Optimization[M].Berlin:Springer,2006.

[20] LENZMANN L,LENZMANN E.Zur Lösung des Nichtlinearen Gauss-Markov-Modells[J].Zeitschrift für Geodäsie, Geoinformation und Landmanagement,2007,132:108-120.

[21] KONG Jian,YAO Yibin,WU Han.Iterative Method for Total Least-squares[J].Science of Wuhan University: Geomatics and Information,2010,35(6):711-714.(孔建,姚宜斌,吴寒.整体最小二乘的迭代解法[J].武汉大学学报:信息科学版,2010,35(6):711-714.)

[22] SCHAFFRIN B,FELUS Y A.On the Multivariate Total Least-squares Approach to Empirical Coordinate Transformations:Three Algorithms[J].Journal of Geodesy, 2008,82(6):373-383.

[23] CAI Jianqing,GRAFAREND EW.Systematical Analysis of the Transformation between Gauss-Krueger-coordinate/ DHDN and UTM-coordinate/ETRS89 in Baden-Württemberg with Different Estimation Methods[C]∥Geodetic Reference Frame International Association of Geodesy Symposia.Munich:International Association of Geodesy,2009,134: 205-211.

(责任编辑:陈品馨)

A General Total Least Squares Algorithm for Three-dimensional Coordinate Transformations

FANG Xing1,ZENG Wenxian1,LIU Jingnan1,2,YAO Yibin1
1.School of Geodesy and Geomatics,Wuhan University,Wuhan 430079,China;2.Research Center of GNSS, Wuhan University,Wuhan 430079,China

The model of three-dimensional coordinate transformation is a nonlinear errors-in-variables model.The methods proposed in published literatures are always restricted in practice for their special assumptions,such as size of rotation angles,structured design matrix and special weight matrix.A general weighted TLS algorithm of a three-dimensional coordinate transformation is investigated in this paper.The new algorithm can be applied in transformations with an arbitrary rotation angles and any applicable weights of the observations,as well as the structured design matrix or the design matrix with both random and fixed elements.The example given in this paper illustrates that this algorithm is general and suits to all kinds of three-dimensional coordinate transformations in practice.

total least-squares method;three-dimensional coordinate transformations;errors-in-variables model;nonlinear program

FANG Xing(1981—),male,PhD,majors in the theory and method of surveying data processing.E-mail:xfang＠sgg.whu.edu.cn

ZENG Wenxian

P207

A

1001-1595(2014)11-1139-05

国家自然科学基金(41474006;41404005;41231174;41174012);中央高校基本科研基金(2042014kf053)

2014-01-13

方兴(1981—),男,博士,现主要从事测量数据处理理论与应用研究。

曾文宪

E-mail:wxzeng＠sgg.whu.edu.cn

FANG Xing,ZENG Wenxian,LIU Jingnan,et al.A General Total Least Squares Algorithm for Three-dimensional Coordinate Transformations[J].Acta Geodaetica et Cartographica Sinica,2014,43(11):1139-1143.(方兴,曾文宪,刘经南,等.三维坐标转换的通用整体最小二乘算法[J].测绘学报,2014,43(11):1139-1143.)

10.13485/j.cnki.11-2089.2014.0193

修回日期:2014-08-29