李　勇，楊德宏

(1.大理州國(guó)土資源規(guī)劃研究院，云南 大理　671000； 2.昆明理工大學(xué) 國(guó)土資源工程學(xué)院，云南 昆明　650093)

?

坐標(biāo)轉(zhuǎn)換求解算法分析與研究*

李勇1，楊德宏2

(1.大理州國(guó)土資源規(guī)劃研究院，云南 大理671000； 2.昆明理工大學(xué) 國(guó)土資源工程學(xué)院，云南 昆明650093)

在大地測(cè)量或工程測(cè)量實(shí)踐中，不可避免的需要解決測(cè)量坐標(biāo)系統(tǒng)建立及坐標(biāo)轉(zhuǎn)換的問(wèn)題。文章以布爾沙轉(zhuǎn)換模型為基礎(chǔ)，對(duì)坐標(biāo)轉(zhuǎn)換的5種求解算法進(jìn)行分析，并以某工程為例，采用MATLAB軟件編程進(jìn)行計(jì)算比對(duì)分析，結(jié)果表明：坐標(biāo)重心化算法求解轉(zhuǎn)換參數(shù)和轉(zhuǎn)換坐標(biāo)，其結(jié)果出現(xiàn)異常，計(jì)算結(jié)果和理論分析發(fā)現(xiàn)，該算法不適合布爾沙轉(zhuǎn)換模型，而其余4種算法求出的坐標(biāo)轉(zhuǎn)換參數(shù)和坐標(biāo)轉(zhuǎn)換結(jié)果相當(dāng)，能保證精度，具有可靠性。

布爾沙模型；坐標(biāo)轉(zhuǎn)換算法；對(duì)比分析

0　引言

實(shí)施坐標(biāo)轉(zhuǎn)換，如果已知轉(zhuǎn)換參數(shù)，那么應(yīng)用轉(zhuǎn)換公式直接進(jìn)行坐標(biāo)轉(zhuǎn)換計(jì)算即可，但實(shí)際面臨的情況往往是轉(zhuǎn)換參數(shù)未知，為此需要使用公共點(diǎn)求解轉(zhuǎn)換參數(shù)，且一般都有多余方程，為得到唯一解，在滿(mǎn)足最小二乘原則的基本前提下有最小二乘算法、坐標(biāo)重心化算法、基于尺度因子算法、最小二乘分解算法等擴(kuò)展算法。本文以布爾沙模型為基礎(chǔ)，以某工程為例，分析不同算法的適宜性和解算精度。

1　布爾沙模型

布爾沙模型計(jì)算公式如下：

2　衡量坐標(biāo)轉(zhuǎn)換精度指標(biāo)

重合點(diǎn)殘差：

V=重合點(diǎn)轉(zhuǎn)換坐標(biāo)值-重合點(diǎn)已知坐標(biāo)值

(2)

空間點(diǎn)位中誤差：

(3)

3　幾種坐標(biāo)轉(zhuǎn)換求解算法分析

3.1傳統(tǒng)最小二乘算法

以重合點(diǎn)的原坐標(biāo)系坐標(biāo)和目標(biāo)坐標(biāo)系坐標(biāo)為觀測(cè)值，建立誤差方程，按VTPV=min的原則組成法方程，進(jìn)而求逆解算法方程得到作為未知數(shù)的轉(zhuǎn)換參數(shù)。

3.2最小二乘分解法

用傳統(tǒng)最小二乘算法，需要對(duì)7階矩陣求逆，坐標(biāo)轉(zhuǎn)換中坐標(biāo)數(shù)據(jù)量較大，可能出現(xiàn)求逆的數(shù)值不穩(wěn)定。為此，可采用最小二乘分解算法，其基本思想是把布爾沙7參數(shù)法模型作線(xiàn)性表達(dá)[2]：

AS=B

(4)

S=[ΔxΔyΔzkεxεyεz]T；

B=[XTYTZT]T。

把系數(shù)矩陣A分解為以下形式：

A=QR

(5)

式中：Q為N×7列的正交矩陣；R為7×7階的上三角矩陣。式(4)表達(dá)為：

RS=QTB

(6)

運(yùn)用以上方法求解坐標(biāo)轉(zhuǎn)換參數(shù)過(guò)程中不需要對(duì)矩陣求逆也可得到坐標(biāo)轉(zhuǎn)換參數(shù)S的最小二乘解。

3.3坐標(biāo)重心化求解算法

工程實(shí)際應(yīng)用中，重合坐標(biāo)點(diǎn)在局部區(qū)域，可能僅分布在空間直角坐標(biāo)系中的某一象限，并不是均勻的分布在空間直角坐標(biāo)系中，其轉(zhuǎn)換過(guò)程中會(huì)影響布爾沙7參數(shù)的精度，增加轉(zhuǎn)換模型的誤差。為此，提出原點(diǎn)坐標(biāo)重心化算法。該方法是先將兩個(gè)坐標(biāo)系的原點(diǎn)平移至公共點(diǎn)的重心位置，重新建立空間直角坐標(biāo)，使坐標(biāo)點(diǎn)在新坐標(biāo)系中均勻分布，進(jìn)而得到7參數(shù)，得到過(guò)渡的布爾沙模型，解算出重心化坐標(biāo)系中的坐標(biāo)，然后再平移換算到相應(yīng)坐標(biāo)系中[3]。

3.4基于尺度因子求解算法

坐標(biāo)轉(zhuǎn)換模型本身是一個(gè)非線(xiàn)性函數(shù)，通常是將非線(xiàn)性轉(zhuǎn)化為線(xiàn)性后求解，這種線(xiàn)性處理增加旋轉(zhuǎn)矩陣中各項(xiàng)旋轉(zhuǎn)參數(shù)的誤差。針對(duì)尺度因子，從兩個(gè)空間直角坐標(biāo)系中同一線(xiàn)段的長(zhǎng)度比出發(fā)，通過(guò)多個(gè)公共點(diǎn)組成的多條公共線(xiàn)段，首先通過(guò)條件平差，得出具有高精度的尺度因子m，然后將尺度因子m作為已知值代入7參數(shù)轉(zhuǎn)換模型中，變成含有6參數(shù)的轉(zhuǎn)換模型，通過(guò)公共點(diǎn)以及最小二乘配置法得出剩余6參數(shù)[4]。

3.5總體最小二乘算法

在高斯-馬爾科夫模型中采用最小二乘配置法對(duì)其誤差方程進(jìn)行求解過(guò)程中，首先假設(shè)觀測(cè)向量L中含有偶然誤差，而其系數(shù)矩陣A中是不含有偶然誤差的。如果系數(shù)矩陣A中存在偶然誤差或者擾動(dòng)時(shí)，基于統(tǒng)計(jì)學(xué)觀點(diǎn)看來(lái)，基于最小二乘配置法求出的估值將不再是最優(yōu)解，是存有偏差的。最小二乘也就是ATPA=min，其中由于A中存在偏差，則其偏差的協(xié)方差因噪聲誤差的影響而增加。所以，在考慮系數(shù)矩陣A中存在偏差時(shí)，應(yīng)考慮采用其他的推廣最小二乘配置方法[5]。由于總體最小二乘算法中，解算方法多樣，本文采用傳統(tǒng)的迭代算法[6]。

4　不同坐標(biāo)轉(zhuǎn)換求解算法的計(jì)算結(jié)果及對(duì)比分析

4.1應(yīng)用實(shí)例

本文采用某測(cè)區(qū)11個(gè)工程實(shí)測(cè)重合點(diǎn)的WGS-84和1980西安坐標(biāo)系，利用MATLAB軟件編程，進(jìn)行計(jì)算對(duì)比分析，坐標(biāo)數(shù)據(jù)見(jiàn)表1，由于測(cè)繪數(shù)據(jù)的保密性，坐標(biāo)數(shù)據(jù)的前5位省略。

表1　工程實(shí)例數(shù)據(jù)

4.2不同坐標(biāo)轉(zhuǎn)換求解算法計(jì)算結(jié)果對(duì)比與分析

在此針對(duì)最小二乘算法，坐標(biāo)重心化算法，最小二乘分解算法，基于尺度因子求解算法，總體最小二乘算法這5種算法進(jìn)行分析。因測(cè)區(qū)呈面狀，1、2、3、4、5號(hào)點(diǎn)為四等點(diǎn)，保存較好，1、2、3號(hào)點(diǎn)相對(duì)分布均勻，6、7、8、9、10、11號(hào)點(diǎn)為一級(jí)點(diǎn)，為此設(shè)計(jì)兩種方案：計(jì)算方案一，以1、2、3號(hào)點(diǎn)作為公共點(diǎn)，4、5號(hào)點(diǎn)作為檢驗(yàn)比對(duì)點(diǎn)；計(jì)算方案二，以1、2、3、4、5號(hào)點(diǎn)作為公共點(diǎn)，10、11號(hào)點(diǎn)作為檢驗(yàn)比對(duì)點(diǎn)。分別以5種方法計(jì)算轉(zhuǎn)換參數(shù)并將WGS-84坐標(biāo)系坐標(biāo)轉(zhuǎn)換為1980西安坐標(biāo)系坐標(biāo)，根據(jù)計(jì)算結(jié)果進(jìn)行對(duì)比分析。結(jié)果顯示：兩種方案計(jì)算結(jié)果一致，因篇幅限制僅列出方案1的計(jì)算成果。

應(yīng)用MATAB軟件編程，得到幾種不同轉(zhuǎn)換算法求解的轉(zhuǎn)換參數(shù)，見(jiàn)表2。

表2　坐標(biāo)轉(zhuǎn)換參數(shù)

表2(續(xù))

轉(zhuǎn)換參數(shù)傳統(tǒng)最小二乘算法重心化求解算法最小二乘分解算法尺度因子求解算法總體最小二乘算法ΔZ4.35329687375.45866666674.3532968564.37606401774.3729663268m-0.0000011261-0.0000089755-0.0000011261-0.0000011343-0.0000011341εX0.00000291052-0.00000032390.000002910520.00000291050.0000029010εY-0.00000668760.0000105875-0.000006687639-0.0000066876-0.00000668751εZ0.000017087480.99999571460.00001708750.00001708750.00001708671

求得轉(zhuǎn)換參數(shù)后，對(duì)4、5號(hào)點(diǎn)進(jìn)行坐標(biāo)轉(zhuǎn)換，得到的坐標(biāo)殘差及點(diǎn)位偏差，見(jiàn)表3。

表3　坐標(biāo)點(diǎn)位殘差及偏差

通過(guò)表2及表3對(duì)比發(fā)現(xiàn)：坐標(biāo)重心化算法求解轉(zhuǎn)換參數(shù)偏差較大，其余4種算法坐標(biāo)轉(zhuǎn)換參數(shù)結(jié)果相當(dāng)；坐標(biāo)重心化求解方法得出的4、5號(hào)點(diǎn)位偏差分別為215.590 7 mm與389.243 4 mm，而其余4種算法得出的4、5號(hào)點(diǎn)位偏差平均值分別為6.550 9 mm與6.321 0 mm。

5　結(jié)論

通過(guò)前文的分析，筆者得出如下結(jié)論：

1)傳統(tǒng)最小二乘算法和最小二乘分解算法結(jié)果一致，尺度因子求解算法和總體最小二乘算法結(jié)果一致，這4種算法求得的坐標(biāo)轉(zhuǎn)換參數(shù)和坐標(biāo)轉(zhuǎn)換結(jié)果相當(dāng)，能保證精度，具有可靠性；

2)坐標(biāo)重心化算法實(shí)施坐標(biāo)轉(zhuǎn)換計(jì)算誤差較大，從理論上分析其根本原因是布爾沙模型旋轉(zhuǎn)參數(shù)和平移參數(shù)相關(guān)性較強(qiáng)[7]，坐標(biāo)重心化算法中因原點(diǎn)重心化改變了平移性質(zhì)進(jìn)而影響旋轉(zhuǎn)參數(shù)求解精度，導(dǎo)致轉(zhuǎn)換誤差較大，實(shí)際應(yīng)用中坐標(biāo)重心化算法不適宜于布爾沙模型。

[1]孔祥元，郭際明，劉宗泉.大地測(cè)量學(xué)基礎(chǔ)[M].武漢：武漢大學(xué)出版社，2002.

[2]魯鐵定，寧津生，周世健，等.最小二乘配置的QR分解解法[J].遼寧工程技術(shù)大學(xué)學(xué)報(bào)：自然科學(xué)版，2009，28(4)：550-553.

[3]田紅平.一種利用坐標(biāo)重心求解坐標(biāo)轉(zhuǎn)換關(guān)系式的方法[J].山西建筑，2011，37(15)：197-198.

[4]李卓昌.坐標(biāo)轉(zhuǎn)換求解算法及模型誤差分析與研究[D].昆明：昆明理工大學(xué)，2015.

[5]陸玨，陳義，鄭波.總體最小二乘方法在三維坐標(biāo)轉(zhuǎn)換中的應(yīng)用[J].大地測(cè)量與地球動(dòng)力學(xué)，2008，28(5)：77-81.

[6]孔建，姚宜斌，吳寒.整體最小二乘的迭代解法[J].武漢大學(xué)學(xué)報(bào)：信息科學(xué)版，2010(6)：711-714.

[7]王解先.七參數(shù)轉(zhuǎn)換中參數(shù)之間的相關(guān)性[J].大地測(cè)量與地球動(dòng)力學(xué)，2007，27(2)：43-46.

Analysis and Research on Algorithm of Coordinate Transformation

LI Yong1，YANG De-hong2

(1.Dali Institute of Land and Resources Planning and Research，Dali Yunnan 671000，China； 2.College of Land and Resources Engineering，Kunming University of Science and Technology，Kunming Yunnan 650093，China)

In the practice of geodetic surveying and engineering surveying，it is inevitable to solve the problem of establishment of measuring coordinate system and coordinate transformation.This paper analyses and studies 5 kinds of coordinate transformation algorithm based on the Bursa model.The related calculation and comparison analysis are carried out by using engineering instance data and MATLAB software programming.The results indicate that using coordinate centralization algorithm to solve conversion parameters and coordinate appear abnormal.Through the calculation results and theoretical analysis，the coordinate centralization algorithm is not suitable for Bursa model.The coordinate conversion parameters and results obtained by the other 4 algorithms are very close,and can ensure the accuracy.Therefore，the other 4 algorithms are reliable.

Bursa model；coordinate transformation algorithm；comparative analysis

2016-04-01

P 226+.3

B

1007-9394(2016)03-0029-03

李勇(1975～)，男，云南大理人，學(xué)士，高級(jí)工程師，現(xiàn)主要從事工程測(cè)量與土地規(guī)劃設(shè)計(jì)方面的工作。