姜英明,紀曉雨,侯曉偉,馮彥同
(1.千尋位置網絡有限公司,上海 200438;2.濰坊職業學院農林科技學院,山東濰坊 262737)
坐標參考框架是坐標系統的具體實現,通常通過一組核心站點的測站坐標和速度表示.坐標參考框架是測繪活動的參考基準[1],同時,高精度的坐標參考框架也是進行空間大地測量學研究所必須的前提[2].各個坐標參考框架之間往往通過Helmert轉換參數來聯系.隨著我國北斗衛星導航系統(BDS)的建成和BDS坐標系的啟用,使用BDS進行精密定位定軌等研究和生產,都需要關注坐標框架的概念[3].為了準確使用各種參考框架,需要分析其轉換的性質.
對于瞬時坐標系或精度要求不高的坐標轉換,一般采用Helmert七參數進行.而對于高精度的坐標參考框架如國際地球參考框架(ITRF),固定的Helmert七參數不能描述其框架的時變特性,因此Helmert七參數被視為是時變的,此時Helmert七參數被擴展成為Helmert十四參數[4].坐標參考框架在建立時,一般都會公布其到其他主流參考框架如ITRF之間的轉換參數,以方便用戶使用.
但在實際使用時,某些參考框架之間可能沒有直接的轉換參數,尤其對于某些區域性的坐標參考框架,此情形并不罕見.此時可通過某些“中間框架”如ITRF進行過渡,即首先將坐標或速度轉換到具有直接轉換參數的某個ITRF框架,再通過ITRF到目標框架的轉換參數,將坐標或速度轉換到目標框架.這造成了實際應用的復雜和不便.基于此類問題,本文推導了由間接的過渡參數計算直接轉換參數的公式,以簡化轉換過程,提高計算效率.
在不考慮框架轉換參數變率的情況下,一般采用Helmert七參數進行坐標轉換.假設從框架1到框架3之間沒有直接的轉換參數,但可以使用框架2進行過渡.記框架1到框架2的轉換參數為T1、D1、R1;框架2到框架3的轉換參數為T2、D2、R2.一般情況下,平移參數T的單位為mm,尺度縮放參數D的單位為10?9(ppb),旋轉矩陣R的單位為千分之一秒(mas).若坐標或角度分別以m或rad為單位,則T、D、R的量級分別為10?3、10?9、10?9.
依據布爾莎模型,將框架1下的某點P的坐標X1先轉換到框架2再轉換到框架3下的計算公式[5]為:


式中:X2為點P在框架2下的坐標;X3為點P在框架3下的坐標.將式(1)代入式(2),可得

整理式(3),可得式(4):

式(4)中,D2(T1+D1X1+R1X1)和R2(T1+D1X1+R1X1)兩項的量級皆為10?12m,在目前的測量精度下都可以忽略,因此在實際計算中,式(4)可以簡化為

對比布爾莎模型的轉換公式,可知框架1到框架3的直接轉換參數T3、D3、R3可用式(6)計算:

ITRF給出的各ITRF參考框架之間的轉換參數也是時變的,因此其轉換參數一般有14個,包括給定歷元t0下的七參數以及其各自的變化率分別為:T、D、R、T˙、D˙、R˙.變化率的單位分別為mm/a、ppb/a和mas/a,其數量級分別為:10?3、10?9、10?9、10?3、10?9、10?9.


整理式(11),可得

后兩項的量級皆為10?12m/a,在計算時可忽略不計,因此速度的轉換公式(12)可簡化寫為

同理,將式(7)、(8)代入式(9),可得坐標的轉換公式:

整理式(14)并忽略微小量,可得

作為對照,式(7)亦可寫成

歐洲參考框架網絡(EUREF)由多個持續觀測的永久性全球衛生導航系統(GNSS)站點構成,用于歐洲地區的參考框架維持.EUREF的數據由法國國家地理學院(IGN)或德國波茲坦地學中心(GFZ)統一處理,并基于其坐標和速度實現ETRS89坐標系.ETRS89是一種地心地固坐標系,其坐標系實現包括ETRF2000、ETRF2005和ETRF2014等[6].EUREF的站點分布由圖1所示。

圖1 EUREF站點分布圖
EUREF公布了其 714個站點分別位于ITRF2014、ETRF2014和ETRF2000框架下的坐標和速度,同時公布了ETRF到各ITRF框架的轉換參數.但是,其公布的成果中并未包含ETRF2014到ETRF2000的轉換參數.本文以將EUREF站點坐標從ETRF2014轉換到ETRF2000為例,驗證由式(17)計算直接轉換參數的可行性.此時坐標轉換有兩種方案:
1)首先將坐標和速度從ETRF2014轉到ITRF2014,然后再從ITRF2014轉換到ETRF2000,記為方法A;
2)使用式(17)計算ETRF2014到ETRF2000的直接轉換參數,然后對ETRF2014框架下的坐標和速度進行轉換,記為方法B.
EUREF公布的2010.0歷元下,ETRF2014轉換到ITRF2014的轉換參數如表1所示:
2010.0歷元下,ITRF2014到ETRF2000的轉換參數如表2所示:

表1 ETRF2014到ITRF2014的轉換參數

表2 ITRF2014到ETRF2000的轉換參數
由式(17)計算得,方法B所使用的2 010.0歷元下ETRF2014到ETRF2000的直接轉換參數如表3所示:

表3 ETRF2014到ETRF2000的轉換參數
作為驗證,首先使用方法A將714個站點在ETRF2014下的坐標和速度轉換到ETRF2000,并將轉換成果與EUREF官方公布的數據進行對比.統計方法A與EUREF官方公布的ETRF2000框架下坐標和速度的差異,統計結果如表4所示:

表4 方法A計算結果與官方數據對比統計
其坐標和速度的差異分布分別如圖2、圖3所示:

圖2 方法A與ETRF2000公布的坐標的差異分布

圖3 方法A與ETRF2000公布的速度的差異分布
從統計可以看出,方法A與EUREF官方公布的站點在ETRF2000框架下的成果,坐標差異均優于1mm,速度差異均優于0.1mm/a,且差異有正有負均勻地分布在0兩側,這說明方法A所進行的轉換是正確的,其轉換差異主要受到偶然誤差的影響.
然后通過方法B使用表3中的轉換參數進行坐標轉換,并將其轉換成果與方法A作為對比.差異分布統計如表5所示:
差異分布如圖4、圖5所示:

表5 方法A與方法B計算結果對比統計

圖4 方法A與方法B坐標差異分布圖

圖5 方法A與方法B速度差異分布圖
由以上統計可以看出,基于方法A和方法B得到的坐標差異均為10?8m量級,速度差異為10?9m/a量級.并且,其差異均大于0,說明二者之間存在系統差.結合方法B的推導過程,易知這是因為方法B的計算過程中忽略了一些微小量所造成的.但由于二者差異足夠小,因此可認為方法A和方法B在一般使用中是等價的.
為驗證式(18),使用參考框架ETRF2014到ITRF2014的轉換參數,將ETRF2014下的基站速度轉換到ITRF2014.作為對比:
1)方法C使用站點在ETRF2014下的坐標;
2)方法D假定站點在ETRF2014下的坐標未知,使用其在ITRF2014下的坐標近似.
兩種方法的差異統計如表6所示。

表6 使用精確和近似坐標進行速度轉換對比
其差異分布如圖6所示。

圖6 使用精確坐標與概略坐標轉換速度的差異分布圖
由以上統計可發現,由方法C和方法D所計算的基站速度差異在10?9m/a量級.二者雖然存在一定的系統差,但在實際使用時完全可以忽略.
利用Helmert轉換公式進行坐標轉換時,對于沒有直接轉換參數的坐標框架,使用多步“過渡”的方式進行轉換時公式復雜低效.本文對傳統的間接坐標轉換方法進行改進,推導了坐標框架間的直接轉換公式,并使用714個EUREF站點在ETRF2014和ETRF2000的坐標和速度進行驗證,得出以下結論:
1)對于Helmert七參數坐標轉換,兩個框架的直接轉換參數可由作為過渡的框架轉換參數求和計算,計算方法如式(6);
2)對于Helmert十四參數坐標轉換,兩個框架的直接轉換參數,其Helmert七參數的變率可由作為過渡的框架轉換參數的變率求和計算,Helmert七參數可將作為過渡的轉換參數歸算到某一相同歷元后求和計算,計算方法如式(17);
3)在對速度進行框架轉換時,若與原框架下速度對應的坐標未知,可使用其他框架下相差不大的坐標代替.
總體而言,坐標轉換的直接參數法和間接參數法差異較小,但本文的轉換方法可簡化計算步驟,各參考框架之間的轉換參數可提前計算后保存,方便程序實現,也可以節約計算資源.