-
由于GPS可以在同一时间获得点位的三维坐标, 测量工作的效率得到了很大的提高, 但是GPS获得的是基于WGS84参考椭球的大地高, 在实际工程应用中, 我国高程系统普遍采用正常高, 如果大地高能够准确地转换为正常高, 现在的测量工作所需要的平面控制测量和高程控制测量就能够同时进行。
对于大多数工程应用, 由于测区范围不大, 高程异常变化平缓, 多项式拟合便可满足精度要求[1], 可通过建立GPS高程异常值与平面坐标之间的多项式关系求解转换参数。在实际参数估计问题中, 对于大部分线性模型来说, 一般由于数据采样大小、模型化及测量等原因引起了系数矩阵的误差, 故其系数阵并非常数阵而是由其他方法计算得到的具有一定误差的变量矩阵[2], 通常可以采用总体最小二乘来进行解算。目前已有相关学者采用总体最小二乘方法进行GPS高程转换参数的估计:文献[3]考虑不同曲面拟合模型中的水平和竖直坐标中含有误差的情形, 分别基于二次多项式曲面、三次多项式曲面、移动二次多项式曲面, 利用最小二乘估计、总体最小二乘和加权总体最小二乘算法进行高程转换; 文献[4]引入基于加权总体最小二乘的二次多项式GPS高程转换, 并采用了一种简化的定权方法; 文献[5]考虑到二次多项式GPS高程转换模型系数矩阵的特点, 采用混合总体最小二乘的方法进行GPS高程转换, 在合理考虑系数矩阵和观测向量误差的情况下, 利用其提出的基于总体最小二乘平差的粗差探测方法进行粗差剔除后, 对GPS高程转换参数进行求解; 文献[6]考虑系数矩阵和观测数据的随机误差, 讨论基于总体最小二乘法建立的高斯-马尔科夫模型, 同时考虑由于粗差的存在, 观测数据精度不相等的情况, 提出稳健总体最小二乘的算法来计算转换参数。
但是, 在现有的关于总体最小二乘方法在GPS高程转换中应用的研究中, 通常都是将转换参数的计算与待求点高程异常计算独立进行, 并且在对高程异常待求点进行高程转换时忽略了待求点本身的坐标误差。因此, 本文在文献[7]的基础上提出GPS高程转换的总体最小二乘拟合推估模型, 该模型考虑到所有点的点位误差, 且将计算转换参数和计算待求点高程异常联合处理, 最后采用拟合推估法进行求解, 并通过模拟算例对该方法的可行性进行了验证。
-
基于二次多项式曲面的GPS高程转换模型为:
(1) 式中, ξi为坐标xi, yi的高程异常; λ=[a0a1a2a3a4a5]T为待估参数。式(1)的矩阵形式可写为:
(2) 其中, Ai=[1xiyixi2xiyiyi2]。
由于控制点的坐标值及高程异常不可避免地受观测误差的影响, 对观测向量及系数矩阵引入误差向量可得:
(3) 其中, EAi=[0ExiEyiExi2ExiyiEyi2]。
联合n个高程异常已知点求解高程异常转换参数, 则:
(4) 式中,
假设有m个高程异常待求点, 计算它们的高程异常:
(5) 式中, 系数矩阵B=
(6) 假设控制点的点位坐标误差相关, 且服从零均值正态分布, 则联合模型的各项误差的随机模型为:
(7) 式中, vec(·)为矩阵向量化算子; Q是协因数阵。
-
拟合推估模型的标准形式为[7]:
(8) 式中, y和ey分别为n×1观测向量和误差向量; β是k×1的参数向量; C是它的设计矩阵且列满秩; y0是n0×1未知随机信号向量; ey0和C0分别为它的误差向量和设计矩阵。通常要求n>k, 但并不要求n0 > k。模型(8)的求解包含平差和预报两个部分, 平差是指采用观测向量y估计参数β, 预报则是指根据信号y0与参数β的函数关系及其与观测向量y的随机关系对其做出合理的推估。若ey和ey0服从正态分布:
(9) 根据最小二乘准则:
(10) 求解的最优线性无偏参数估计值和预报值分别为:
(11) (12) -
GPS高程转换联合模型是非线性模型, 本文将根据高斯-牛顿法[7]和转换参数的数值特点, 分别将式(6)转化为线性模型, 并采用拟合推估法进行求解。
-
采用高斯-牛顿法将联合模型线性化后进行迭代计算, 设第j次迭代后的估值为λ(j), 则参数λ表示为:
(13) 式中, δλ(j)为第j+1次迭代的参数改正数。再将式(13)代回式(6)可得:
(14) 尽管EAδλ(j)和EBδλ(j)是二阶的小量, 但是完全忽略二阶以上的小量, 可能导致迭代发散或收敛错误。采用EA和EB的第j次迭代估值EA(j)和EB(j)代替EA和EB计算系数矩阵可避免该缺陷。记A(j)=A-EA(j), B(j)=B-EB(j), 则(14)可变为:
(15) 将上式改写成拟合推估模型的标准形式:
(16) 式中, In和Im分别表示n维和m维单位阵;
(17) 顾及式(17) ey和ey0及式(7), 根据协方差传播律推出Qyy, Qy0y分别为:
(18) (19) 根据上述中的无偏参数和预报值的公式可以推出:
(20) (21) (22) (23) 其中, 式(20)、(21)式用于估计高程异常转换参数, 式(22)、(23)用于计算待求点的高程异常值。由于预报和估计是相对独立的过程, 因此只需要按照式(20)、(21)计算转换参数, 迭代终止后, 再按式(22)、(23)计算最后的推估值。
-
求解两种模型需要已知QAA、QBA, 定义变量
由系数矩阵的表达式, 结合误差传播定律可以得到:
(24) 并可推出
则n个点系数矩阵误差可以写成:
(25) 分别令:
根据误差传播定律, 可以得到:
(26) 式中,
-
选取14个点, 将其中的A、B两点作为起始坐标点, 模拟距离观测值, 建立一个测边网, 计算其他12个点的坐标平差值及其协方差阵, 选取2、6、8、9号点作为检验点, 各点的坐标真值如表 1所示。图 1为这14个点的空间分布, 可以看出所有数据点的空间分布较为均匀。令(b0, b1, b2, b3, b4, b5)=(-11, 0.05, 0.1, 0.01, 0.009, -0.02)为二次多项式曲面高程异常转换参数的真值[6], 并根据坐标真值及参数真值计算高程异常的真值。
图 1 测边网点的空间分布
Figure 1. Spatial Distribution of Trilateration Net's Point
表 1 各点坐标、高程异常/m
Table 1. Coordinates and Height Anomaly of Points/m
点号 X Y ξ A 100 100 - B 120 120 - 1 120 140 -87.8 3 160 120 149.8 4 170 200 -187.5 5 150 190 -225 7 130 220 -524.1 10 110 180 -336.3 11 170 160 35.3 12 200 170 144 2 140 160 -102.4 6 190 140 220.9 8 180 100 294 9 200 220 -151 反算两点间的边长作为距离的真值S0, 并对S0添加服从正态分布N(0, σ2)的随机误差作为距离的观测值, 其中, σ=a+b×S, 本文中计算取常数a=5 mm, b=0.05 mm/km。联合所有距离误差方程, 权取距离方差的倒数, 计算得到这12个点的坐标平差值及其协方差阵。再给高程异常的真值添加服从正态分布N0, σξ2随机误差, 作为高程异常的观测值, 其中σξ2为测边网平差过程中计算得到的单位权方差, 得到的平差坐标值及含有误差的高程异常值如表 2所示。
表 2 各点平差坐标、加入误差的高程异常/m
Table 2. Adjusted Coordinates and Add the Abnormal Height of Error/m
点号 X Y ξ A 100 100 - B 120 120 - 1 119.997 9 140.000 6 -87.798 0 3 159.989 1 120.001 7 149.795 9 4 170.011 9 199.980 9 -187.495 4 5 150.013 1 189.987 9 -224.997 6 7 130.016 7 219.995 4 -524.100 7 10 110.008 0 180.000 5 -336.299 9 11 170.006 2 159.987 1 35.300 2 12 200.013 5 169.973 9 144.003 0 2 140.003 4 159.996 6 -102.396 2 6 190.003 3 139.983 0 220.899 8 8 179.962 3 100.003 3 293.999 2 9 200.005 6 219.983 8 -151.001 1 本文中采用的方案如表 3所示, 其中方案5所采用的定权方法, 是根据协因数传播定律求得协因数之间的相互关系, 并考虑到协因数会随点的不同而变化, 为了方便起见, 对定权的方法做了一定的简化[4]; 方案4则引用了文献[8]中采用的加权总体最小二乘迭代计算的方法, 该方法直接采用系数矩阵的协方差阵进行计算。
表 3 各方案所用方法列表
Table 3. List of Methods Used in Different Schemes
各个方案的计算结果如表 4所示, 从计算结果的比较中可以看出, 采用传统最小二乘方法和总体最小二乘方法的计算结果差别并不大; 两种加权总体最小二乘法解算精度较好, 方案5采用的是简化的定权方法, 而方案4根据协因数阵进行计算, 方案4的计算方法相对方案5而言更加合理, 且方案4的转换效果也优于方案5。但是, 不论采用何种加权总体最小二乘方法, 仍然存在着目前关于总体最小二乘方法在GPS高程转换的研究中普遍存在的不足, 即将转换参数的计算与待求点高程异常计算分两步进行处理, 并且只考虑由高程异常已知点的平面坐标组成的系数矩阵的误差, 忽略了高程异常待求点的坐标误差。方案1将计算转换参数和计算检验点高程异常联合处理, 且考虑到所有点的点位误差, 从理论上来说该方案要更为严谨与合理, 从计算的结果可以看出, 该方案相对其他4种方案, 转换精度有所提高。
表 4 检验点转换的结果及与真值的残差范数/m
Table 4. Checkpoint's Result of the Transformation and Norm of Residual Between the Estimated Value and True Value/m
方案 点号 
‖Δ 真值 -102.4 220.9 294 -151 - 方案1 -102.394 2 220.901 8 293.896 6 -151.149 4 0.181 7 方案2 -102.413 5 220.914 5 293.929 7 -151.186 5 0.200 3 方案3 -102.413 5 220.914 5 293.929 7 -151.186 5 0.200 3 方案4 -102.405 1 220.893 1 293.900 7 -151.154 6 0.183 9 方案5 -102.413 5 220.914 5 293.929 7 -151.186 5 0.199 3 -
本文考虑到由高程异常待求点组成的系数矩阵中的误差, 并且将计算高程异常转换参数与计算待求点的高程异常联合处理, 提出了GPS高程转换拟合推估模型, 结合模拟数据进行转换参数的计算以及解算结果精度的评定。模拟算例的计算结果表明, 本文方法能够得到更高精度的高程转换结果。
-
摘要: 在现有的关于GPS高程转换的总体最小二乘方法研究中，通常是将高程异常转换参数的计算与待求点高程异常的计算分两步进行处理，并且只考虑由已知高程异常点的平面坐标组成的系数矩阵的误差，忽略了高程异常待求点的坐标误差。针对以上问题，本文提出了GPS高程转换的总体最小二乘拟合推估模型，将计算高程异常转换参数和待求点高程异常联合处理，且考虑到所有点的点位误差，最后采用拟合推估法进行求解。实验结果表明，本文方法能够有效地提高高程转换的精度。Abstract: In the current research on the total least squares method in the conversion of GPS height, the calculation of the conversion parameter and elevation abnormities of the check points are generally performed in two steps, and only consider the error in the coefficient matrix used to calculate the parameters; errors in the coordinate of the check point are ignored. In view of this gap, we put forward a total least squares fitting estimation model of GPS height transformation, that combines the calculation of fitting parameters with the calculation of elevation abnormities at inspection points, and considers the position error of all points. Collocation calculation experiemental results verify the feasibility of this method. These test results show that the method can effectively improve the accuracy of elevation conversion.
-
Key words:
- total least squares /
- collocation /
- GPS elevation /
- elevation transformation /
- polynomial fitting
-
表 1 各点坐标、高程异常/m
Table 1. Coordinates and Height Anomaly of Points/m
点号 X Y ξ A 100 100 - B 120 120 - 1 120 140 -87.8 3 160 120 149.8 4 170 200 -187.5 5 150 190 -225 7 130 220 -524.1 10 110 180 -336.3 11 170 160 35.3 12 200 170 144 2 140 160 -102.4 6 190 140 220.9 8 180 100 294 9 200 220 -151 
下载: 导出CSV
表 2 各点平差坐标、加入误差的高程异常/m
Table 2. Adjusted Coordinates and Add the Abnormal Height of Error/m
点号 X Y ξ A 100 100 - B 120 120 - 1 119.997 9 140.000 6 -87.798 0 3 159.989 1 120.001 7 149.795 9 4 170.011 9 199.980 9 -187.495 4 5 150.013 1 189.987 9 -224.997 6 7 130.016 7 219.995 4 -524.100 7 10 110.008 0 180.000 5 -336.299 9 11 170.006 2 159.987 1 35.300 2 12 200.013 5 169.973 9 144.003 0 2 140.003 4 159.996 6 -102.396 2 6 190.003 3 139.983 0 220.899 8 8 179.962 3 100.003 3 293.999 2 9 200.005 6 219.983 8 -151.001 1
下载: 导出CSV
表 3 各方案所用方法列表
Table 3. List of Methods Used in Different Schemes
下载: 导出CSV
表 4 检验点转换的结果及与真值的残差范数/m
Table 4. Checkpoint's Result of the Transformation and Norm of Residual Between the Estimated Value and True Value/m
方案 点号 ‖Δ 真值 -102.4 220.9 294 -151 - 方案1 -102.394 2 220.901 8 293.896 6 -151.149 4 0.181 7 方案2 -102.413 5 220.914 5 293.929 7 -151.186 5 0.200 3 方案3 -102.413 5 220.914 5 293.929 7 -151.186 5 0.200 3 方案4 -102.405 1 220.893 1 293.900 7 -151.154 6 0.183 9 方案5 -102.413 5 220.914 5 293.929 7 -151.186 5 0.199 3
下载: 导出CSV
-
[1] 徐绍铨, 张华海, 杨志强, 等.GPS测量原理及应用[M].武汉:武汉大学出版社, 2008 Xu Shaoquan, Zhang Huahai, Yang Zhiqiang, et al. GPS Measurement Principle and Application[M]. Wuhan: Wuhan University Press, 2008 [2] 王乐洋, 许才军.总体最小二乘研究进展[J].武汉大学学报·信息科学版, 2013, 38(7): 850-856 http://www.cnki.com.cn/Article/CJFDTOTAL-WHCH201307022.htm Wang Leyang, Xu Caijun. Progress in Total Least Square[J]. Geomatics and Information Science of Wuhan University, 2013, 38(7): 850-856 http://www.cnki.com.cn/Article/CJFDTOTAL-WHCH201307022.htm [3] 丁海勇, 孙景领.GPS高程转换的总体最小二乘方法研究[J].大地测量与地球动力学, 2013, 33(3):52-55 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201303011.htm Ding Haiyong, Sun Jingling. Research on Total Least-Squares Methods for Transformation of GPS Elevation[J]. Journal of Geodesy and Geodynamics, 2013, 33(3): 52-55 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201303011.htm [4] 赵辉, 张书毕, 张秋昭.基于加权总体最小二乘法的GPS高程拟合[J].大地测量与地球动力学, 2011, 31(5): 88-90 http://www.cnki.com.cn/article/cjfdtotal-dkxb201105018.htm Zhao Hui, Zhang Shubi, Zhang Qiuzhao. GPS Height Fitting of Weighted Total Least-Squares Adjustment[J]. Journal of Geodesy and Geodynamics, 2011, 31(5):88-90 http://www.cnki.com.cn/article/cjfdtotal-dkxb201105018.htm [5] 龚循强, 陈磬, 周秀芳.总体最小二乘平差方法在GPS高程拟合中的应用研究[J].测绘通报, 2014(3): 6-8 http://www.cnki.com.cn/Article/CJFDTOTAL-CHTB201403002.htm Gong Xunqiang, Chen Qing, Zhou Xiufang. The Application Research of GPS Height Fitting Based on TLS Adjustment Method[J]. Bulletin of Surveying and Mapping, 2014(3): 6-8 http://www.cnki.com.cn/Article/CJFDTOTAL-CHTB201403002.htm [6] Tao Y Q, Gao J X, Yao Y F. TLS Algorithm for GPS Height Fitting Based on Robust Estimation[J]. Survey Review, 2014, 46(336):184-188 doi: 10.1179/1752270613Y.0000000083 [7] 李博峰, 沈云中, 李薇晓.无缝三维基准转换模型[J].中国科学(地球科学), 2012, 42(7): 1 047-1 054 http://www.cnki.com.cn/Article/CJFDTOTAL-JDXK201207012.htm Li Bofeng, Shen Yunzhong, Li Weixiao. The Seamless Model for Three-Dimensional Datum Transformation[J]. Sci China Earth Sci, 2012, 42(7): 1 047-1 054 http://www.cnki.com.cn/Article/CJFDTOTAL-JDXK201207012.htm [8] 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 [9] 鲁铁定, 宁津生.总体最小二乘平差理论及其应用[M].北京:中国科学技术出版社, 2011 Lu Tieding, Ning Jinsheng. Total Least Squares Adjustment Theory and Its Applications[M]. Beijing: China Science and Technology Press, 2011 [10] 王乐洋.测边网坐标的总体最小二乘平差方法[J].大地测量与地球动力学, 2012, 32(6): 81-85 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201206018.htm Wang Leyang. Trilateration Net's Coordinate Adjustment Based on Total Least Squares[J]. Journal of Geodesy and Geodynamics, 2012, 32(6): 81-85 http://www.cnki.com.cn/Article/CJFDTOTAL-DKXB201206018.htm [11] 孔建, 姚宜斌, 吴寒.整体最小二乘的迭代解法[J].武汉大学学报·信息科学版, 2010, 35(6): 711-714 http://www.cnki.com.cn/Article/CJFDTOTAL-WHCH201006024.htm 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://www.cnki.com.cn/Article/CJFDTOTAL-WHCH201006024.htm [12] 陶本藻, 邱卫宁, 姚宜斌.误差理论与测量平差基础[M].武汉:武汉大学出版社, 2009 Tao Benzao, Qiu Weining, Yao Yibin. Error Theory and Fundation of Surveying Adjustment[M].Wuhan: Wuhan University Press, 2009 [13] 李丽华, 高井祥, 刘晋斌.协方差函数拟合高程异常方法探析[J].测绘工程, 2005, 14(4):33-35 http://www.cnki.com.cn/Article/CJFDTOTAL-CHGC200504009.htm Li Lihua, Gao Jingxiang, Liu Jinbin. Discussion of the Methods of the Covariance Function Prediction for the Fitting of Abnormal Height[J]. Engineering of Surveying and Mapping, 2005, 14(4):33-35 http://www.cnki.com.cn/Article/CJFDTOTAL-CHGC200504009.htm [14] 王增利, 黄腾, 邓标.基于二次曲面的拟合推估法在GPS高程测量中的应用[J].测绘工程, 2009, 18(1):50-52 http://www.cnki.com.cn/Article/CJFDTOTAL-CHGC200901016.htm Wang Zengli, Huang Teng, Deng Biao. Application of Collocation Model Based on Quadric Surface in GPS Leveling Surveying[J]. Engineering of Surveying and Mapping, 2009, 18(1):50-52 http://www.cnki.com.cn/Article/CJFDTOTAL-CHGC200901016.htm [15] Teunissen P J G, Simons D, Tiberius C. Probability and Observation Theory[M]. The Netherlands: Delft University of Technology, 2008 -
点击查看大图
图(1) / 表(4)
计量
- 文章访问数: 1175
- HTML全文浏览量: 43
- PDF下载量: 418
- 被引次数: 0
