基于非迭代與迭代法聯(lián)合估計的七參數(shù)坐標轉(zhuǎn)換方法研究*

譚駿祥 李少達 楊容浩

(成都理工大學地球科學學院，成都 610059)

1 引言

空間直角坐標轉(zhuǎn)換常采用七參數(shù)模型［1］，針對七參數(shù)模型的參數(shù)求解，大致分為迭代解和非迭代解法兩類。根據(jù)旋轉(zhuǎn)矩陣求解方法不同，非迭代解法一般采用SVD 分解法［2］、正交矩陣法［3］、單位四元數(shù)法［4］、對偶四元數(shù)法［5］等;根據(jù)旋轉(zhuǎn)矩陣表達方式不同，迭代解法常使用歐拉角形式［6］、正交矩陣形式［7］和四元數(shù)形式［8］。

在已有七參數(shù)模型求解方法的基礎(chǔ)上，本文首先簡述了迭代解法與非迭代解法，然后通過實測數(shù)據(jù)對比分析了兩類解法的優(yōu)缺點，進而提出聯(lián)合兩類解法來估計七參數(shù)，以期實現(xiàn)任意旋轉(zhuǎn)角度的穩(wěn)健高精度空間直角坐標轉(zhuǎn)換。

2 模型的表示與算法評價

2.1 模型表示

據(jù)研究，一般化的三維坐標轉(zhuǎn)換非線性模型為:

式中，vS=［x，y，z］T和vT=［X，Y，Z］T分別表示待轉(zhuǎn)換坐標和目標坐標，vθ=［ΔX，ΔY，ΔZ］T為平移向量，x=［ΔX，ΔY，ΔZ，s，φ，θ，γ］T為七參數(shù)向量，s 為尺度參數(shù)，R 為旋轉(zhuǎn)矩陣。R 繞z、x、y 軸旋轉(zhuǎn)φ、θ、γ后組成的序列為:

假定e 為隨機誤差且e ～N(0，δ2)，則兩個坐標系中相對應的控制點的關(guān)系為:

式中i=1，2，…，M，M 為控制點的個數(shù)。

為了減小求解過程中坐標值的有效位數(shù)，提高計算過程的精度，對模型進行重心化處理，得

2.2 算法評價

使用均方根誤差(RMSE，Root Mean Square Error)度量轉(zhuǎn)換后坐標與目標坐標的接近程度，即

由于七參數(shù)的真值未知，我們根據(jù)轉(zhuǎn)換殘差來判斷求解的準確性。若RMSE 超過限差值e 時，即RMSE ＞e，則認為七參數(shù)的估計值錯誤，模型求解不準確。需要指出的是，實際坐標轉(zhuǎn)換中的尺度參數(shù)估值常接近1，否則，也會認為求解錯誤。

3 模型的最小二乘估計

3.1 非迭代解法

使用式(5)進行非迭代解法估計參數(shù)的分步求解，求解過程為:求尺度參數(shù)→旋轉(zhuǎn)矩陣→平移向量。其中尺度參數(shù)s 的求解式為:

求解旋轉(zhuǎn)矩陣R 常用的是單位四元數(shù)法，根據(jù)單位四元數(shù)法，先計算N

根據(jù)最大特征值對應的單位特征向量q=［q0，q1，q2，q3］，

3.2 迭代解法

式(4)按泰勒級數(shù)展開取至1 次項得:

式中v=Bx－l，其中v=［dvT，1，dvT，2，…，dvT，M］T，B=［B1，B2，…，BM］T，Bi=［I3，RvS，i，svS，i］，l=［l1，l2，…，lM］T，li=vT，i－v0－sRvS，i，I3表示3 階單位矩陣，x 為所求未知數(shù)改正數(shù)。由于R 有歐拉角、正交矩陣和四元數(shù)等3 種表達形式，而歐拉角形式有3個，正交矩陣形式有9 個，四元數(shù)形式有4 個。而表達旋轉(zhuǎn)只需要3 個必要參數(shù)，因此正交矩陣形式需要6 個約束條件，四元數(shù)形式只需1 個約束條件。

1)正交矩陣形式

2)四元數(shù)形式

設(shè)單位四元數(shù)q=［a，b，c，d］，構(gòu)造旋轉(zhuǎn)矩陣為:

約束條件為a2+b2+c2+d2=1。

式(8)表示的七參數(shù)線性化模型使用間接平差方法可求得x=(BTB)－1BTl。由于按泰勒級數(shù)展開進行線性化會產(chǎn)生誤差，因此采用迭代運算，迭代過程直至迭代前后兩次x 之差小于某一微小量為止。

4 算法對比

為了驗證七參數(shù)求解方法的精度和準確性，坐標轉(zhuǎn)換實驗采用兩組實測數(shù)據(jù)進行。第一組數(shù)據(jù)的控制點轉(zhuǎn)換前后坐標系分別為WGS84 和西安80 坐標系(表1)，覆蓋范圍比較大，最遠兩控制點間隔約80 千米，但兩個坐標的旋轉(zhuǎn)角度相對較小。第二組數(shù)據(jù)的控制點來自具有一定重疊度的兩站激光掃描數(shù)據(jù)(表2)，范圍很小，但旋轉(zhuǎn)角度比較大。

表1 第一組數(shù)據(jù)的控制點坐標(單位:m)Tab.1 Coordinates of control points of the first set data(unit:m)

表2 第二組數(shù)據(jù)的控制點坐標(單位:m)Tab.2 Coordinate of control points of the second set data(unit:m)

實驗中，非迭代解法不需要設(shè)置參數(shù);迭代求解中7 參數(shù)的近似值是未知的，分別設(shè)置v0、u 和旋轉(zhuǎn)角度的初始值為［0，0，0］T、1 和［0，0，0］T，設(shè)置最大迭代次數(shù)為300(表3、4)。

從表3 看，對于小角度的坐標轉(zhuǎn)換，非迭代解法與迭代解法均能求得轉(zhuǎn)換參數(shù)的正確估計值，其中非迭代解法的正交矩陣法的迭代次數(shù)最少，求解精度最高。從表4 看，對于角度較大的坐標轉(zhuǎn)換，非迭代解法仍然能求得轉(zhuǎn)換參數(shù)的正確估值，可迭代解法不能夠收斂至正確的值。由表3、4 可見:非迭代解法不需要迭代，均能收斂至正確的解，但精度可能不如迭代解法高;迭代解法需要迭代求解，可收斂至更高精度，但當旋轉(zhuǎn)角度過大時易收斂至錯誤解。

5 聯(lián)合估計七參數(shù)

非迭代解法采用分步求解七參數(shù)，不需要迭代，在大小旋轉(zhuǎn)角度下均能獲得正確估值，但沒有考慮同時包含七參數(shù)的誤差模型，精度可能較低;迭代解法采用同時求解七參數(shù)，并通過迭代的方式逐步逼近最佳待求解參數(shù)，不足之處在于迭代的收斂性易受初始值的影響，不適用于大旋轉(zhuǎn)角的坐標轉(zhuǎn)換。針對兩類解法的優(yōu)缺點，提出聯(lián)合兩類解法來估計七參數(shù)，具體步驟為:

表3 第一組數(shù)據(jù)的求解結(jié)果對比Tab.3 Comparison of solutions of the first set data

表4 第二組數(shù)據(jù)的求解結(jié)果對比Tab.4 Comparison of solutions of the second set data

1)使用非迭代解法求解轉(zhuǎn)換前坐標vS至轉(zhuǎn)換后坐標vT的初始轉(zhuǎn)換參數(shù)xini，包括旋轉(zhuǎn)矩陣Rini(或旋轉(zhuǎn)角度φini，θini，γini)、尺度參數(shù)sini和平移向量v0ini;

2)根據(jù)式(1)，利用xini將vS轉(zhuǎn)換至vTini，此時vTini與vT已經(jīng)相當接近;

3)設(shè)置v0、u 和旋轉(zhuǎn)角度的初始值為［0，0，0］T、1 和［0，0，0］T，利用迭代解法求解vTini至vT的轉(zhuǎn)換參數(shù)xend;

4)根據(jù)式(1)，利用xend將vTini轉(zhuǎn)換至vT，完成坐標轉(zhuǎn)換。

求解中，分別聯(lián)合常用的非迭代解法單位四元數(shù)法與三種迭代解法估計七參數(shù)，即將單位四元數(shù)法的解作為迭代解法的初始值進行迭代求解，于是產(chǎn)生三種聯(lián)合解法。利用表1、2 數(shù)據(jù)的聯(lián)合求解結(jié)果見表5、6。

表5 第一組數(shù)據(jù)的聯(lián)合求解結(jié)果Tab.5 The results of combined estimation for the first set of data

表6 第二組數(shù)據(jù)的聯(lián)合求解結(jié)果Tab.6 The results of combined estimation for the second set of data

由表5、6 可見，兩組數(shù)據(jù)均能獲得正確的結(jié)果，體現(xiàn)了此方法的正確性。其中，正交矩陣的迭代次數(shù)最少，并且能獲得最高的轉(zhuǎn)換精度，四元數(shù)形式次之，歐拉角的效果最差。與表3、4 相比較，聯(lián)合解法的精度幾乎高于非迭代解法，在大旋轉(zhuǎn)角度下也能獲得正確的估計值。

6 結(jié) 論

為了解決非迭代解法精度較低，收斂性易受初始值的影響、不適用于大旋轉(zhuǎn)角坐標轉(zhuǎn)換的問題，本文提出聯(lián)合兩類解法來估計七參數(shù)。通過算例驗證了聯(lián)合解法的準確性與有效性，其中MatrixUQ 的迭代次數(shù)最少，轉(zhuǎn)換后的坐標與目標坐標最接近。本文提出了一種不需迭代初始值而能應用于任意角度轉(zhuǎn)換的聯(lián)合估計七參數(shù)方法，但并沒有將其與非線性最小二乘［9］等求解方法作對比分析，還有待進一步完善。

1 王之卓.攝影測量原理［M］.北京:測繪出版社.1990.(Wang Zhizhuo.Principle of photogrammetry［M］.Beijing:Surveying and Mapping Press，1990)

2 Arun K S，Huang T S and Blostein S D.Least-square fitting of two 3D point sets［J］.IEEE Transactions on Pattern Analysis and Machine Intelligence，1987，9(5):698－700.

3 Horn B K P，Hilden H M nand Negahdaripour S.Closed-form solution of absolute orientation using orthonormal matrices［J］.Journal of the Optical Society of America A，1988，5(7):1 128－1 135.

4 Horn B K P.Closed-form solution of absolute orientation using unit quaternions［J］.Journal of the Optical Society of America A，1987，4(4):629－642.

5 Walker M W，Shao L and Volz R A.Estimating 3-D location parameters using dual number quaternions［J］.Image Understanding，1991，54(3):358－367.

6 曾文憲，陶本藻.三維坐標轉(zhuǎn)換的非線性模型［J］.武漢大學學報(信息科學版)，2003，28(5):566－568.(Zeng Wenxian and Tao Benzao.Non-linear adjustment model of three-dimensional coordinate transformation［J］.Geomatics and Information Science of W uhan University，2003，28(5):566－568)

7 陳義，沈云中，劉大杰.適用于大旋轉(zhuǎn)角的三維基準轉(zhuǎn)換的一種簡便模型［J］.武漢大學學報(信息科學版)，2004，29(12):1 101－1 105.(Chen Yi，Shen Yunzhong and Liu Dajie.A simplified model of three dimensional-datum transformation adapted to big rotation angle［J］.Geomatics and Information Science of Wuhan University，2004，29(12):1 101－1 105)

8 趙雙明，等.基于四元數(shù)的三維空間相似變換解算［J］.武漢大學學報(信息科學版)，2009，34(10):1 214－1 217.(Zhao Shuangming，et al.Quaternion-based 3D similarity transformation algorithm［J］.Geomatics and Information Science of Wuhan University，2009，34(10):1 214－1 217)

9 陳宇，白征東.基于非線性最小二乘算法的空間坐標轉(zhuǎn)換［J］.大地測量與地球動力學，2010，(2):129－132.(Chen Yu and Bai Zhengdong，An nonlinear least squares algorithm for spatial coordinate transformation［J］.Journal of Geodesy and Geodynamics，2010，(2):129－132)