-
二维基准转换在实际工程中使用频繁，不同平面坐标系统之间的转换就显得非常重要，而常用的二维基准转换方法为相似变换或正形投影变换[1]。为了求解转换参数，必须已知足够的公共点，建立高斯-马尔科夫(Gauss-Markov,G-M)平差模型[2]。实际上采用的公共点坐标本质上也是通过测量手段获取的，公共点在不同坐标系下的坐标均含有误差。若公共点两套坐标都含有误差，则平面坐标转换中误差方程的系数矩阵B不再是完全准确的，这种情况下经典的G-M模型将不再适用[3]。
为了解决坐标转换中公共点在两套坐标系下的坐标均含有误差这一问题，近年来，众多研究提出了各种不同的坐标转换模型或理论，无论是在二维平面还是三维空间坐标转换中都具有良好的改善效果。例如无缝三维基准转换模型l[4]、非线性坐标转换的稳健估计[5]、整体最小二乘法[6]、加权总体最小二乘法[7]、改进的加权总体最小二乘法[8]、约束总体最小二乘法以及约束加权混合整体最小二乘法[3, 9]等。以上所有模型都是考虑系数矩阵B含有误差eB，基于目标函数elTQl-1el+eBTQB-1eB=min，或同时考虑限制条件和定权方法而建立的，其模型复杂，求解过程较为繁琐。为了高效、简便地获取基准转换参数，本文在二维平面坐标系统下，基于相似变换原理，探讨了坐标参数化的最小二乘(coordinates parameterization least squares,CPLS)平面基准四参数转换方法。基本思想是：考虑公共点在两套平面坐标系下的坐标均含有误差，将公共点在 Ⅰ 套坐标系下的坐标(XⅠ，YⅠ)和转换参数作为未知参数，按照目标函数VⅡTPⅡVⅡ=min求解转换参数并评定精度，VⅡ为公共点在Ⅱ套坐标系中的坐标改正数。
-
平面基准转换的相似变换原理为：
(1) 其矩阵形式为：
(2) 式中，公共点i=1，2，…，n；(XⅡi,YⅡi)、(XⅠi,YⅠi)分别为第i个公共点在Ⅱ套和Ⅰ套坐标系下的坐标；[X0Y0]T为平移参数；θ为旋转角；m为尺度因子。
令Z1=mcosθ、Z2=msinθ，以[X0 Y0 Z1 Z2]T为参数，按照G-M模型建立误差方程：
(3) 式中，[δX0δY0δZ1δZ2]ⅡiT为参数改正数；(X0,YⅡi0)为第i个公共点在Ⅱ套坐标系下的近似坐标，可以通过参数近似值和公共点的Ⅰ套坐标推算得出；[vXⅡ1vYⅡ1…vXⅡnvYⅡn]T为各点在Ⅱ套坐标系下坐标的改正数。
传统的LS平面基准转换不考虑Ⅰ套坐标系下[XⅠiYⅠi]ⅡT的误差，只考虑Ⅱ套坐标系下[XⅡiYⅡi]T的误差，按照VTPⅡVⅡ=min的准则进行平差计算，其中VⅡ=[VXⅡVYⅡ]T。不考虑公共点在Ⅰ套坐标系下坐标的误差，显然不符合实际情况。为此衍生出了TLS等多种改进的求解模型和方法[2-9]。文献[2-9]中所提方法均以系数矩阵B含有误差eB为基础，对整体最小二乘坐标转换进行推导和分析。此类方法求解过程复杂，所需模型较多而且步骤繁琐。因此，本文探讨了更加简便的坐标参数化的平面基准转换方法，该方法的核心是将公共点在Ⅰ套坐标系下的坐标作为未知参数进行平差，从而解决公共点两套坐标均含有误差这一问题，使得转换参数更加准确。
-
对式(2) 进行变形得：
(4) 将尺度参数m与Ⅰ套坐标(XⅠ，YⅠ)相乘构造缩放后的Ⅰ套坐标(a,b)。本节内容根据公共点的Ⅰ套缩放坐标和Ⅱ套坐标求解平移参数和旋转角，尺度参数m在§1.3节中解算。
令ai=mXⅠi，bi=mYⅠi，c=cosθ，d=sinθ，则式(1) 可转化为：
(5) 设XⅡi、X0、Y0、ai、bi、c、d的近似值分别为XⅡi0、X00、Y00、ai0、bi0、c0、d0，相应的改正数为vXⅡi、δX0、δYo、δai、δbi、δc、δd，忽略微小量δai、δbi、δc、δd之间的乘积，则式(5) 可转化为：
(6) 式中，δai、δbi为第i个公共点Ⅰ套缩放坐标的改正数。
考虑到公共点在两套坐标系下的坐标都含有误差，则平面基准转换的观测方程为：
(7) 设vai、vbi为第i个公共点Ⅰ套缩放坐标的残差，vai=ai－ai0，vbi=bi－bi0，将vai、vbi和式(6) 代入到式(7) 中得误差方程为：
(8) 方程中观测值个数为4n，未知数个数为2n+4，误差方程的矩阵形式为：
(9) 式中，I为2n阶单位阵；A为2n阶对角阵；B为2n×4阶矩阵。
通过以上步骤，考虑公共点在两套坐标系下的坐标误差，建立了G-M平差模型，按照间接平差法进行平差计算，假设观测值的权为单位权，由式(9) 得最小二乘解为：
(10) -
根据式(10) 可得：
(11) 由于公共点Ⅰ套缩放坐标ai=mXⅠi，bi=mYⅠi，为了求解尺度参数m，由全微分公式得：
(12) 式中，m0为尺度近似值；δm为尺度改正数；(XⅠi0,YⅠi0)为公共点在Ⅰ套坐标系下的坐标。实例证明了在尺度接近1而且坐标转换范围为小区域的情况下，Ⅰ套坐标系下各点坐标改正数的微小差异对尺度的精度影响可以忽略不计。因此，假设各公共点在Ⅰ套坐标系下的坐标改正数相同，采用虚拟的坐标改正数(δXⅠ,δYⅠ)与尺度改正数δm构成误差方程(13) 的解向量。
式(12) 的矩阵形式为：
(13) 式(13) 为超定方程组，采用超定方程组的最小二乘解法可得：
(14) 从而得到尺度参数：
(15) 将式(10) 代入式(9) 可以得到CPLS平面基准转换的残差：
(16) 按照目标函数F=VTV=min，可得验后单位权中误差：
(17) -
本文按照平移参数X0=5.000 0 m、Y0=8.000 0 m、尺度参数m=1.0和旋转角θ=50.0°构造了6个公共点在两套坐标系下的坐标(${\hat{X}}$Ⅰ,${\hat{Y}}$Ⅰ)和(${\hat{X}}$Ⅱ，${\hat{Y}}$Ⅱ)，然后对所有的坐标加0~2 mm的随机误差，采用含有误差的坐标作为算例分析的样本数据，公共点在两套坐标系下的坐标如表 1所示。对表 1中含扰动误差的数据分别采用LS和CPLS两种方法解算参数并进行精度评定，对比结果如表 2所示。为了比较不同方法求解参数的可靠性，分别采用两种方法解算的转换参数，根据公共点在Ⅰ套坐标系下的坐标计算公共点在Ⅱ套坐标系下的坐标(XⅡLS，YⅡLS)和(XⅡCPLS，YⅡCPLS)，并与XⅡ、YⅡ比较计算坐标差ΔXⅡLS、ΔYⅡLS、ΔXⅡCPLS和ΔYⅡCPLS，如表 3所示。将ΔXⅡCPLS、ΔYⅡCPLS的绝对值与ΔXⅡLS、ΔYⅡLS的绝对值坐标差得ΔXⅡCPLS-LS、ΔYCPLS-LSⅡ，如表 4所示。
表 1 平面基准转换样本点/m
Table 1. Sample Points of Plane Datum Transformation/m
点号 无误差 加扰动误差 ${\hat{X}}$Ⅰ ${\hat{Y}}$Ⅰ ${\hat{X}}$Ⅱ ${\hat{Y}}$Ⅱ XⅠ YⅠ XⅡ YⅡ 1 1.000 00 8.000 00 11.771 14 12.376 26 1.000 88 8.000 27 11.772 68 12.376 43 2 1.000 00 9.000 00 12.537 19 13.019 04 1.001 77 9.001 94 12.538 21 13.020 22 3 1.000 00 10.000 00 13.303 23 13.661 83 1.001 06 10.001 81 13.304 65 13.663 61 4 3.000 00 8.000 00 13.056 72 10.844 17 3.001 42 8.000 09 13.058 37 10.845 02 5 3.000 00 9.000 00 13.822 76 11.486 96 3.000 51 9.000 20 13.823 73 11.487 82 6 3.000 00 10.000 00 14.588 81 12.129 74 3.001 99 10.000 83 14.590 30 12.130 81 表 2 平面基准转换不同方法解算的精度和参数估值
Table 2. Estimation of Precision and Parameters from Different Plane Datum Transformation Methods
方法 σ0 /mm X0 /m Y0 /m θ/ (°) 尺度参数 LS 0.60 4.999 68 8.000 45 49.995 938 1.000 054 99 CPLS 0.42 5.002 47 7.996 23 49.996 494 1.000 000 00 表 3 平面基准转换不同方法计算的Ⅱ套坐标及坐标差
Table 3. Coordinates and Difference in Coordinate System Ⅱ from Different Plane Datum Transformation Methods
点号 XⅡLS/m YⅡLS/m XⅡCPLS/m YⅡCPLS/m ΔXⅡLS/mm ΔYⅡLS/mm ΔXⅡCPLS/mm ΔYⅡCPLS/mm 1 11.771 66 12.376 92 11.771 66 12.376 15 1.02 -0.49 1.02 0.28 2 12.539 55 13.020 20 12.539 16 13.019 80 -1.34 0.02 -0.95 0.42 3 13.305 04 13.663 53 13.304 30 13.663 51 -0.39 0.08 0.35 0.09 4 13.057 61 10.844 31 13.058 38 10.844 32 0.75 0.71 -0.01 0.70 5 13.823 16 11.487 96 13.823 53 11.488 35 0.58 -0.14 0.20 -0.53 6 14.590 63 12.130 11 14.590 65 12.130 87 -0.33 0.70 -0.35 -0.07 表 4 平面基准转换不同方法计算的Ⅱ套坐标差之差/mm
Table 4. Difference of Coordinate Difference from Different Plane Datum Transformation Methods/mm
点号 |ΔXⅡLS| |ΔYⅡLS| |ΔXⅡCPLS| |ΔYⅡCPLS| ΔXⅡCPLS-LS ΔYⅡCPLS-LS 1 1.02 0.49 1.02 0.28 0.00 -0.21 2 1.34 0.02 0.95 0.42 -0.39 0.39 3 0.39 0.08 0.35 0.09 -0.04 0.02 4 0.75 0.71 0.01 0.70 -0.74 -0.01 5 0.58 0.14 0.20 0.53 -0.37 0.39 6 0.33 0.70 0.35 0.07 0.02 -0.63 从表 2中可以看出，两种方法的验后单位权中误差分别为σ0LS=0.60 mm、σ0CPLS=0.42 mm。这说明基于坐标参数化的平面基准转换精度高于普通的平面基准转换方法。结合LS、CPLS解算的参数估值与表 3中两种方法计算的坐标差可以看出，CPLS解算的坐标较差总体上更小，证明该方法解算的参数更加接近真值。为了比较两种方法的差异，采取了不同方法计算的Ⅱ套坐标差之差这一指标，如表 4所示，取绝对值来表示反算的坐标与样本点坐标的偏离程度，坐标较差的绝对值越大，偏离程度越大，反之就越小。用ΔXⅡCPLS-LS、ΔYⅡCPLS-LS可以比较出两种方法反算坐标偏离程度的大小，根据大量仿真实例计算的统计结果得出,CPLS反算坐标偏离程度小于LS反算坐标偏离程度的比例为70%左右，即CPLS有效地提高了平面基准转换的精度，表 2~表 4中的精度指标和坐标差等数据说明CPLS较LS可以求得更高精度的平面基准转换参数。
-
本文从理论推导和实例计算上验证了坐标参数化的平面基准转换方法要优于普通的平面基准转换方法，因为基于CPLS的方法同时考虑了公共点在两套坐标系下的坐标误差，在选择合适公共点的前提下，CPLS获得的基准转换参数更符合实际情况。
仿真实例的统计结果表明，CPLS验后单位权中误差小于LS验后单位权中误差的比例为100%，CPLS反算公共点在Ⅱ套坐标系下的坐标偏离程度小于LS反算坐标偏离程度的比例约为70%。CPLS验后单位权方差较LS有所改进，前者反算公共点在Ⅱ套坐标系下的坐标偏离程度较后者有所降低。
A Method of Plane Datum Transformation by Coordinate Parameterization
-
摘要: 工程测量中经常需要实现不同坐标系下成果的相互转换，而高精度的转换参数是完成这一工作的基础。获取基准转换参数的实质就是利用公共点在两套坐标系下的坐标，根据相似变换原理建立误差方程求解。传统的最小二乘（LS）相似变换法只考虑了公共点在一套坐标系下的误差，与实际情况不符。基于此，探讨了坐标参数化的平面基准转换方法，解决了考虑公共点在两套坐标系下坐标都含有误差时高斯-马尔科夫（Gauss-Markov，G-M）模型不成立的问题，以相似变换原理为基础，采取通用的最小二乘方法解算基准转换参数。Abstract: The mutual transformation of coordinates in different coordinate systems is a very common geoprocessing task. High-accuracy transformation parameters are the foundation for coordinate transformations, obtained by similarity transformations. Traditional least squares similarity transformation only considers common point coordinates in one coordinate system, and does not conform to many practical situations. Because common points have errors in two coordinate systems, the G-M model is not correct. The coordinates parameterization least squares, method overcomes this problem. Based on the principle of similar transformation, the traditional least squares method is used to calculate plane datum transformation parameters.
-
表 1 平面基准转换样本点/m
Table 1. Sample Points of Plane Datum Transformation/m
点号 无误差 加扰动误差 ${\hat{X}}$Ⅰ ${\hat{Y}}$Ⅰ ${\hat{X}}$Ⅱ ${\hat{Y}}$Ⅱ XⅠ YⅠ XⅡ YⅡ 1 1.000 00 8.000 00 11.771 14 12.376 26 1.000 88 8.000 27 11.772 68 12.376 43 2 1.000 00 9.000 00 12.537 19 13.019 04 1.001 77 9.001 94 12.538 21 13.020 22 3 1.000 00 10.000 00 13.303 23 13.661 83 1.001 06 10.001 81 13.304 65 13.663 61 4 3.000 00 8.000 00 13.056 72 10.844 17 3.001 42 8.000 09 13.058 37 10.845 02 5 3.000 00 9.000 00 13.822 76 11.486 96 3.000 51 9.000 20 13.823 73 11.487 82 6 3.000 00 10.000 00 14.588 81 12.129 74 3.001 99 10.000 83 14.590 30 12.130 81 
下载: 导出CSV
表 2 平面基准转换不同方法解算的精度和参数估值
Table 2. Estimation of Precision and Parameters from Different Plane Datum Transformation Methods
方法 σ0 /mm X0 /m Y0 /m θ/ (°) 尺度参数 LS 0.60 4.999 68 8.000 45 49.995 938 1.000 054 99 CPLS 0.42 5.002 47 7.996 23 49.996 494 1.000 000 00
下载: 导出CSV
表 3 平面基准转换不同方法计算的Ⅱ套坐标及坐标差
Table 3. Coordinates and Difference in Coordinate System Ⅱ from Different Plane Datum Transformation Methods
点号 XⅡLS/m YⅡLS/m XⅡCPLS/m YⅡCPLS/m ΔXⅡLS/mm ΔYⅡLS/mm ΔXⅡCPLS/mm ΔYⅡCPLS/mm 1 11.771 66 12.376 92 11.771 66 12.376 15 1.02 -0.49 1.02 0.28 2 12.539 55 13.020 20 12.539 16 13.019 80 -1.34 0.02 -0.95 0.42 3 13.305 04 13.663 53 13.304 30 13.663 51 -0.39 0.08 0.35 0.09 4 13.057 61 10.844 31 13.058 38 10.844 32 0.75 0.71 -0.01 0.70 5 13.823 16 11.487 96 13.823 53 11.488 35 0.58 -0.14 0.20 -0.53 6 14.590 63 12.130 11 14.590 65 12.130 87 -0.33 0.70 -0.35 -0.07
下载: 导出CSV
表 4 平面基准转换不同方法计算的Ⅱ套坐标差之差/mm
Table 4. Difference of Coordinate Difference from Different Plane Datum Transformation Methods/mm
点号 |ΔXⅡLS| |ΔYⅡLS| |ΔXⅡCPLS| |ΔYⅡCPLS| ΔXⅡCPLS-LS ΔYⅡCPLS-LS 1 1.02 0.49 1.02 0.28 0.00 -0.21 2 1.34 0.02 0.95 0.42 -0.39 0.39 3 0.39 0.08 0.35 0.09 -0.04 0.02 4 0.75 0.71 0.01 0.70 -0.74 -0.01 5 0.58 0.14 0.20 0.53 -0.37 0.39 6 0.33 0.70 0.35 0.07 0.02 -0.63
下载: 导出CSV
-
[1] 黄丁发.GPS卫星导航定位技术与方法[M].北京:科学出版社,2009 Huang Dingfa.GPS Satellite Navigation and Positioning Technologies and Methods[M].Beijing:Science Press,2009 [2] 陆珏,陈义,郑波.总体最小二乘方法在三维坐标转换中的应用[J].大地测量与地球动力学, 2008,28(5):77-81 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB200805016.htm Lu Jue,Chen Yi,Zheng Bo.Applying Total Least Squares to 3D Datum Transformation[J].Journal of Geodesy and Geodynamics, 2008, 28(5):77-81 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB200805016.htm [3] 陆珏,陈义,郑波.三维非线性基准转换的加权和加约束总体最小二乘解算模型[J].大地测量与地球动力学,2012,32(6):129-134 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201206029.htm Lu Jue,Chen Yi,Zheng Bo.Performing Nonlinear 3D Datum Transformation Using WTLS and CTLS Method[J].Journal of Geodesy and Geodynamics, 2012, 32(6):129-134 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201206029.htm [4] 李博峰,沈云中,李微晓.无缝三维基准转换模型[J].中国科学:地球科学,2012,42(7):1047-1054 http://www.cnki.com.cn/Article/CJFDTOTAL-JDXK201207012.htm Li Bofeng, Shen Yunzhong, Li Weixiao.The Seamless Model for 3D Datum Transformation[J]. Sci China Earth Sci,2012, 42(7):1047-1054 http://www.cnki.com.cn/Article/CJFDTOTAL-JDXK201207012.htm [5] 陈义,沈云中.非线性三维基准转换的稳健估计[J].大地测量与地球动力学,2004, 23(4):49-53 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB200304010.htm Chen Yi,Shen Yunzhong.Robust Estimation of 3D Nonlinear Datum Transformation[J].Journal of Geodesy and Geodynamics,2004, 23(4):49-53 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB200304010.htm [6] 孔建,姚宜斌,吴寒. 整体最小二乘的迭代解法[J].武汉大学学报·信息科学版, 2010, 35(6):711-714 http://ch.whu.edu.cn/CN/abstract/abstract970.shtml Kong Jian,Yao Yibin,Wu Han.Iterative Method for Total Least-Squares[J]. Geomatics and Information Science of Wuhan University,2010,35(6):711-714 http://ch.whu.edu.cn/CN/abstract/abstract970.shtml [7] 张鹏杰,邱卫宁,程进伟,等.加权整体最小二乘求解线性模型参数及精度估计[J].测绘信息与工程, 2012, 37(1):4-5 http://www.cnki.com.cn/Article/CJFDTOTAL-CHXG201201003.htm Zhang Pengjie,Qiu Weining,Cheng Jinwei,et al.Weighted Total Least Squares to Solve Linear Model Parameters and Accuracy Estimation[J].Journal of Geomatics,2012,37(1):4-5 http://www.cnki.com.cn/Article/CJFDTOTAL-CHXG201201003.htm [8] 杨仕平,范东明,龙玉春.加权整体最小二乘算法的改进[J].大地测量与地球动力学,2013,33(1):48-52 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201301012.htm Yang Siping,Fan Dongming,Long Yuchun.An Improved Weighted Total Least Squares Algorithm[J].Journal of Geodesy and Geodynamics,2013, 33(1):48-52 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201301012.htm [9] 葛旭明,伍吉仓.三维基准转换的约束加权混合整体最小二乘的迭代解法[J].武汉大学学报·信息科学版,2012, 37(2):178-182 http://ch.whu.edu.cn/CN/abstract/abstract124.shtml Ge Xuming,Wu Jicang.Iterative Method of Weight Constraint Total Least-Squares for 3D Datum Transformation[J].Geomatics and Information Science of Wuhan University,2012,37(2):178-182 http://ch.whu.edu.cn/CN/abstract/abstract124.shtml [10] 徐天河, 杨元喜.坐标转换模型尺度参数的假设检验[J].武汉大学学报·信息科学版,2001,26(1):70-74 http://ch.whu.edu.cn/CN/abstract/abstract5041.shtml Xu Tianhe,Yang Yuanxi. The Hypothesis Testing of Scale Parameter in Coordinate Transformation Model[J].Geomatics and Information Science of Wuhan University,2001, 26(1):70-74 http://ch.whu.edu.cn/CN/abstract/abstract5041.shtml [11] 曾文宪,陶本藻.三维坐标转换的非线性模型[J].武汉大学学报·信息科学版,2003, 28(5):566-568 http://ch.whu.edu.cn/CN/abstract/abstract4861.shtml Zeng Wenxian,Tao Benzao.Non-linear Adjustment Model of Three-Dimensional Coordinate Transformation[J].Geomatics and Information Science of Wuhan University,2003,28(5):566-568 http://ch.whu.edu.cn/CN/abstract/abstract4861.shtml [12] 曾怀恩,黄声享.三维坐标转换参数求解的一种直接搜索法[J].武汉大学学报·信息科学版,2008,33(11):1118-1121 http://ch.whu.edu.cn/CN/abstract/abstract1760.shtml Zeng Huai'en,Huang Shengxiang.A Kind of Direct Search Method Adapted to Solution of 3D Coordinate Transformation Parameters[J].Geomatics and Information Science of Wuhan University,2008,33(11):1118-1121 http://ch.whu.edu.cn/CN/abstract/abstract1760.shtml -
