韓 笑，蒲英霞

(1. 南京大學(xué) 地理與海洋科學(xué)學(xué)院，江蘇 南京 210023；2. 江蘇省地理信息技術(shù)重點實驗室，江蘇 南京 210023；3. 江蘇省地理信息資源開發(fā)與利用協(xié)同創(chuàng)新中心，江蘇 南京 210023)

0 引 言

空間基準和坐標系統(tǒng)是大地測量、地圖制圖和“3S”技術(shù)中非常重要的組成部分，為了獲取不同用途的地理空間數(shù)據(jù)，我們通常需要在不同空間基準（如北京1954、西安1980、WGS-84、CGCS2000等）和不同坐標系統(tǒng)（地理坐標、投影坐標、空間直角坐標等）之間相互轉(zhuǎn)換[1]。如當(dāng)采用GPS數(shù)據(jù)和地面控制點數(shù)據(jù)建立具有空間參考的衛(wèi)星影像時，需要將GPS數(shù)據(jù)和地面控制點數(shù)據(jù)統(tǒng)一到相同的坐標系中；當(dāng)合并不同數(shù)據(jù)來源的地圖數(shù)據(jù)時，也需要先統(tǒng)一他們的空間參考系然后進行合并處理。

坐標變換是將空間數(shù)據(jù)從一個坐標系統(tǒng)轉(zhuǎn)換到另一個坐標系統(tǒng)的過程，其實質(zhì)是根據(jù)兩個坐標系下的公共已知點求解坐標變換參數(shù)[2]。程新輝等人采用四參數(shù)法將北京1954轉(zhuǎn)換至西安1980坐標系，并進行了轉(zhuǎn)換誤差的分析[3]；王解先等人討論了將WGS-84地心坐標轉(zhuǎn)換到北京1954直角坐標的兩個模型——平面轉(zhuǎn)換模型和空間轉(zhuǎn)換模型，并且比較了空間轉(zhuǎn)換模型中三參數(shù)和七參數(shù)對模型轉(zhuǎn)換結(jié)果的影響[4]。黎舒等研究了從西安1980到CGCS2000坐標系的變換方法，并根據(jù)具體情況提出了二維四參數(shù)、七參數(shù)和三維七參數(shù)等模型[5]。在坐標變換的誤差處理方面，較為成熟且應(yīng)用最廣的是在求解參數(shù)的過程中采用最小二乘法[6]?？捉ɑ谡w最小二乘法推導(dǎo)了四參數(shù)仿射變換的參數(shù)估計公式，避免了常規(guī)矩陣分解方法計算的復(fù)雜性[7]；張鵬杰等提出利用加權(quán)整體最小二乘法進行參數(shù)求解，同時考慮原始坐標和目標坐標的誤差，并根據(jù)誤差的影響程度賦予不同的權(quán)值[8]?，F(xiàn)有的坐標變換方法研究中，大多是基于給定的兩個坐標系統(tǒng)來研究某一種坐標轉(zhuǎn)換方法，往往是獨立闡述，沒有試圖挖掘不同轉(zhuǎn)換方法之間的聯(lián)系。在利用最小二乘法推導(dǎo)參數(shù)估計表達式時，我們發(fā)現(xiàn)不同坐標變換方法的推導(dǎo)過程存在很大的相似性。本文在Allan研究的基礎(chǔ)上[9]，引入設(shè)計矩陣的概念，以此作為紐帶將不同的坐標變換方法聯(lián)系起來，推導(dǎo)出了基于最小二乘法的統(tǒng)一參數(shù)估計表達式。

1 不同坐標轉(zhuǎn)換方法

在進行坐標變換的過程中，若已知兩個不同坐標系統(tǒng)下足夠數(shù)量的坐標點，即觀測值個數(shù)大于必要觀測數(shù)時，則可以采用最小二乘法求解坐標轉(zhuǎn)換參數(shù)的最優(yōu)估計值。一般地，根據(jù)一組已知數(shù)據(jù)集，利用最小二乘法解算未知參數(shù)的過程如下：

首先，建立已知坐標數(shù)據(jù)與未知坐標轉(zhuǎn)換參數(shù)之間的一般函數(shù)關(guān)系式

在公式（1）中，當(dāng)達到最優(yōu)估計時，可得如下的微分關(guān)系式

將（3）式代入（4）式中，可以簡化為

A，C是雅克比矩陣，在這里稱之為設(shè)計矩陣[9]。

實際上，已知點的觀測值總會因為量綱或者觀測條件的不同而產(chǎn)生系統(tǒng)誤差，為減小誤差，我們在求解過程中計算殘差平方的加權(quán)和，因而引入權(quán)矩陣W，可以推導(dǎo)出以下表達式：

式中，W是一個對角矩陣，且對角線元素對應(yīng)著各自方差的倒數(shù)[9]。

基于公式（6）待求參數(shù)的一般表達式，問題的關(guān)鍵就變成求解不同坐標變換模型對應(yīng)的設(shè)計矩陣。下面針對幾個常用的坐標轉(zhuǎn)換模型，給出設(shè)計矩陣的推導(dǎo)過程。

1.1 相似變換

對于二維坐標系下的相似變換，我們通常采用四參數(shù)模型，即經(jīng)過坐標系的平移、旋轉(zhuǎn)和比例變化實現(xiàn)[10-11]。該方法適用于坐標軸正交且各個方向的尺度因子相同的平面坐標系之間的轉(zhuǎn)換，需要至少兩個控制點來估計變換參數(shù)，具體表達如下：

設(shè)（Xs,Ys）表示源坐標系下的坐標值，（XT,YT）表示目標坐標系下的坐標值，則有函數(shù)關(guān)系式：

式中，a0、b0分別為X、Y方向的平移參數(shù)，a=μcosα，b=μsinα，α為旋轉(zhuǎn)角度，μ為兩坐標軸的尺度參數(shù)。

設(shè)計矩陣A，C可以表示為：

進一步解算出：

以上的演算過程是針對單一坐標點的。通常情況下，會提供多個已知點的坐標值，這時我們可以用Ai，Ci來表示第i個點的設(shè)計矩陣，將所有已知點的設(shè)計矩陣組合在一起，可得表達式：

1.2 仿射變換

仿射變換是基于仿射坐標系而建立的一種坐標變換數(shù)學(xué)模型，是經(jīng)過坐標系的平移、比例、旋轉(zhuǎn)、對稱和錯切等復(fù)合變換得到的[12]。仿射變換適用于坐標軸不正交或者各個方向尺度因子不同的參考系，要求至少3個控制點來求解變換參數(shù)，其基本關(guān)系式為：

式中，a0、b0分別為X、Y方向的平移參數(shù)，a1=μxcosα，b1=μxsinα；a2=μysinα ，b2=μycosα，α為旋轉(zhuǎn)角度，μx、μy分別為X、Y軸的尺度參數(shù)。特別地，當(dāng) μx=μy時，有a1=b2，a2=b1，此時的六參數(shù)仿射變換關(guān)系式與四參數(shù)相似變換相同，即四參數(shù)相似變換是六參數(shù)仿射變換的特殊形式。

第i個點的設(shè)計矩陣Ai，Ci：

將所有已知點的設(shè)計矩陣組合在一起：

1.3 七參數(shù)相似變換法

三維地心參考系下的坐標變換，通常采用七參數(shù)變換模型[13]，即3個平移量、3個旋轉(zhuǎn)角和1個尺度因子，該變換需要至少3個以上的控制點，其模型表達式為：

式中，ΔX、ΔY、ΔZ分別為XS、YS、ZS方向的偏移量，R為旋轉(zhuǎn)矩陣，αx、αy、αz分別為XS、YS、ZS坐標軸的旋轉(zhuǎn)角。

一般地，比例因子接近于1，所以我們可以用（1+κ）表示比例因子，κ為極小量。第i個點的設(shè)計矩陣為：

這種轉(zhuǎn)換方法最普遍的應(yīng)用是一個大地基準向另一個大地基準的轉(zhuǎn)化，這些情況下旋轉(zhuǎn)角都非常小（幾乎只是幾秒），因此，XS＞＞ YS，則有：

同理，

將以上代入（16）、（18）式中，有：

閉合差向量：

將所有已知點的設(shè)計矩陣組合在一起，有：

2 算例驗證

本文選取了某地5個控制點分別在北京1954和西安1980坐標系統(tǒng)下的坐標值，采用二維坐標系統(tǒng)下的四參數(shù)相似變換和六參數(shù)仿射變換模型，推導(dǎo)設(shè)計矩陣，求解變換參數(shù)。表1給出了上述5個控制點分別在兩個不同坐標系統(tǒng)下的高斯平面直角坐標值。

表1 控制點在北京1954和西安1980坐標系下的坐標值Tab.1 Coordinates of control points in Beijing 1954 and Xi'an 1980

應(yīng)用公式(8)，可以計算出相似變換下的設(shè)計矩陣A、C，代入公式（6）求得相似變換參數(shù)：a0=-229.386 7，b0=-366.942 8，a=1.000 008 3，b=0.000003 84。

同理，應(yīng)用公式(12)，求出仿射變換下的設(shè)計矩陣A和C，代入公式（6）求得變換參數(shù)：a0=-200.398 66，a1=1.000 006 24，a2=0.000 003 306，b0=-303.194 67，b1=-0.000 008 182，b2=1.000 007 11。

選取CGCS2000和西安1980坐標系下的F、G、H、I 4個控制點坐標，需要將CGCS000的大地坐標和西安1980的高斯平面坐標轉(zhuǎn)化到空間直角坐標系中（見表2），然后采用三維坐標系下的七參數(shù)變換模型，推導(dǎo)出設(shè)計矩陣，求解變換參數(shù)。

表2 控制點在CGCS2000和西安1980下的坐標值Tab.2 Coordinates of control points in CGCS2000 and Xi'an 1980

根據(jù)公式(19)、（20）,列出七參數(shù)變換的設(shè)計矩陣A和C，由公式（6）計算得到相應(yīng)的變換參數(shù)：k=-2.219E-06，αx=-9.506 7E-07，αy=9.218E-06，αz=-1.295E-05，ΔX=190.248，ΔY=103.828，ΔZ=30.551 6。

接下來，將3種變換方法計算得到的變換參數(shù)應(yīng)用到該區(qū)域內(nèi)的11個控制點，分別計算這些點在北京1954轉(zhuǎn)西安1980的相似變換、北京1954轉(zhuǎn)西安1980的仿射變換、CGCS2000轉(zhuǎn)西安1980的七參數(shù)相似變換下的坐標值，比較計算值與真實值的差異(表3，表4，表5)。相似變化的X、Y均方根誤差分別為0.1353和0.1755，仿射變化的X、Y均方根誤差分別為0.1407和0.0127，七參數(shù)相似變化的X、Y、h上的均方根誤差分別為0.0033、0.0022和0.5178，均在可接受范圍內(nèi)。比較相似變換和仿射變換的均方根誤差可以看出，相似變換與仿射變換在X方向的均方根誤差比較接近，但在Y方向的均方根誤差相差較大。如我們所知，假設(shè)坐標軸正交，且各個方向的尺度因子相同，則仿射變換的表達式與相似變換相同，即a1=b2=a，a2=-b1=b。回顧上面參數(shù)估計的結(jié)果發(fā)現(xiàn)，a1、b2的估計值與a的估計值相近，a2的估計值與b的估計值也較為接近，但-b1的值與b的值相差較大，因此，相似變換與仿射變換在Y坐標的計算上差異較大。不同于前面兩者，七參數(shù)相似變換是適用于三維直角坐標系的，其考慮的參數(shù)更加充分，因此，在X、Y方向的誤差值相對較小，計算結(jié)果更接近真實值。

表3 相似變換的計算值和真值比較Tab.3 Comparison between calculated values and true values in similarity transformation

表4 仿射變換計算值與真值比較Tab.4 Comparison between calculated values and true values in aきne transformation

表5 七參數(shù)相似變換計算值與真值比較Tab.5 Comparison between calculated values and true values in sevenparameter similarity transformation

3 結(jié)束語

本文在常用坐標轉(zhuǎn)換模型研究的基礎(chǔ)上，引用設(shè)計矩陣的概念，根據(jù)最小二乘法原理推導(dǎo)了坐標變換統(tǒng)一的參數(shù)估計表達式。以四參數(shù)相似變換、六參數(shù)仿射變換和三維七參數(shù)相似變換為例，分別闡述了應(yīng)用設(shè)計矩陣求解變換參數(shù)的過程，并將求得的參數(shù)應(yīng)用于其他控制點，比較計算值與真實值的誤差，有效驗證了引用設(shè)計矩陣的合理性。與前人的工作相比，本文的特點在于以設(shè)計矩陣為紐帶，將不同的坐標變換方法聯(lián)系起來，對于給定的坐標變換方法，只要求出其設(shè)計矩陣，代入?yún)?shù)估計表達式進行求解，簡化了計算過程。需要注意的是，在進行坐標變換時應(yīng)該根據(jù)點的分布情況，合理地選擇控制點。在滿足控制點合理分布的條件下，控制點個數(shù)越多，計算結(jié)果往往更加精確。特別地，自2008年啟用CGCS2000坐標系以來，隨著國家經(jīng)濟建設(shè)的發(fā)展，地方坐標系測繪成果向CGCS2000轉(zhuǎn)換的需求越來越大，本文的研究有助于提高地理坐標系轉(zhuǎn)換的效率，推廣CGCS2000坐標系的使用。