李 軍,高詠梅
(中國科學院,測量與地球物理研究所動力大地測量學重點實驗室,湖北武漢430077)
準確、可靠的地殼運動速度場的獲得是地殼運動及動力學研究的基礎,而GPS資料的解算主要依賴三大主力軟件,即GAMIT、Bernese和GIPSY。很多研究表明:采用不同的軟件、不同的解算方案,不同的作者處理同樣的GPS資料往往得到不同的解算結果。中國地殼運動觀測網絡1998~2000年GPS觀測資料計算的幾個結果顯示:不同研究機構采取不同的處理軟件和方案得到的速度場的差異絕對值平均為2~3 mm/a[1]。弄清這些差異產生的原因,并采用相應的解決處理方案對我們獲得一致的、可靠的地殼運動速度場至關重要。
目前,中國科學院測量與地球物理研究所與日本東京大學地震研究所和韓國天文及空間科學研究院正在開展東亞地區地殼運動的合作研究。根據合作研究計劃,中日韓三方分別采用GAMIT、Bernese和GISPY軟件獨立處理相同的GPS觀測資料,這為我們探討三種軟件處理結果差異提供了一個良好的機會[2]。
采用東亞地區的GPS觀測資料進行三種軟件數據處理的對比研究,日本東京大學采用GIPSY軟件,韓國天文研究院采用BERNESE軟件,而我們則采用GAMIT/GLOBK以及BERNESE軟件。三方對相同時段的觀測資料進行處理,得到統一參考框架下的解算結果。通過對處理方案以及處理結果的分析,研究差異性的原因,并探討得到一致的處理結果的合理方案。
此次對比研究,我們選取了東亞地區13個GPS站點,包括中國的3個站點,BJFS(北京)、WUHN(武漢)、SHAO(上海);韓國的2個站點,SUWN和DAEJ;日本的8個站點:0094、0202、0252、0603、0729、0767、USUD、TSKB。資料的時間段為2002年中每個月第一個星期的數據,共84天。在處理過程中,為了能客觀地比較三種軟件的處理結果,我們盡量采用相同的處理參數,如相同的參考框架(ITRF2000),相同的慣性系(J2000.0),相同的衛星截止角(15°)以及相同的精密衛星星歷。
GAMIT/GLOBK是由MIT和SIO共同研制的GPS綜合分析軟件[3]。我們分三個步驟采用GAMIT/GLOBK進行資料處理。首先,利用偽距和載波相位觀測資料的雙差組合求得臺站坐標和衛星軌道的單日松弛解。此步計算中,GPS衛星的軌道和測站坐標都不在一個穩定的參考框架里,整個GPS網作為一個剛體在平移和旋轉。其次,利用GLOBK將每天的單日松弛解和圣地亞哥海洋研究所軌道中心(SOPAC)所計算出的全球IGS跟蹤站的多個單日松弛解合并,得到一個包含全球GPS測站的合并松弛解。最后,對一系列的單日解進行卡爾曼濾波,估算出測站的位置和速度。在估算位置和速度時,通過對全球IGS核心站的約束來進行參考框架的定義和7參數轉換,從而得到ITRF2000參考架下的測站的坐標和速度[4]。
GIPSY是由美國噴氣推進實驗室(JPL)研制的GPS數據處理軟件,它和GAMIT以及BERNESE的最大區別在于直接處理載波非差觀測量[5]。GIPSY處理的方案是首先解算出測站的三維坐標,并給出坐標的方差-協方差,形成單日時段解。然后,將所有的單日時段解轉換到統一的坐標框架下。同樣,我們利用全球的IGS核心站進行7參數轉換,將所有單日時段解作整網平差,平差過程中對所有站都不給約束,最后得到ITRF2000參考框架下測站的坐標和速度。
BERNESE是由瑞士BERNE大學開發研制的GPS及GLONASS資料處理軟件[6]。我們利用BERNESE處理東亞GPS資料的步驟和方案是,利用雙差數據進行最小二乘參數估計,解算出各時段的基線和模糊度,然后,對同步測區各時段的觀測進行整體平差,最后,利用法方程堆積方法對非同步觀測結果進行綜合序貫平差。我們在運用BERNESE做網平差的時候采取兩種方案,一種是固定網內精度比較高的IGS站來進行方差估計及位移速度解算,得到ITRF2000參考架下測站的坐標和速度;另外一種方法就是不固定任何站而對整個觀測網進行方差估計以及位移速度解算,得到IT RF2000參考框架下測站的坐標和速度。
依照上面的方案步驟,三個單位分別用三種處理軟件包對東亞13個GPS站點2002年84天的觀測資料進行了解算,并對解算的結果進行了分析。首先,分析BERNESE按第一種方案(也就是固定某幾個參考站來定義參考框架)處理的結果。這里,我們引進IGS處理的結果來作為參考值,由于IGS中心分析了各個分處理中心的結果,并進行加權平均,坐標的精度和可靠性都非常高,所以,我們采用他們提供的坐標值來作為參考值。表1列出了這三種軟件的結果統計,我們將每個臺站每天的IT RF2000參考架下的坐標值分別與IGS公布的相同歷元相同參考架下的坐標值之差計算出來,分別對三種軟件得到的13個臺站84天的殘差求平均,按X,Y,Z三個方向列出來。按照上面介紹的方案,BERNESE軟件處理的結果與IGS的結果差在三個方向上分別為0.071 m,0.019 m和0.025 m,都差到厘米級。GIPSY的結果在 X,Y方向差為7 mm和2 mm,Z方向也達到了1.2個cm。而GAMIT的結果顯示差都在4個mm以下。

表1 三種軟件計算結果統計
在圖1中,我們列出了其中一個代表臺站TSKB的坐標序列圖,依次為坐標的三個方向X、Y、Z,每個序列圖的橫坐標為時間,縱坐標則是三種軟件解算出來的測站坐標與IGS公布的差,其中符合“+”表示是BERNESE處理的結果,符號“▲”表示的是GAMIT處理的結果,而符號“○”則代表GIPSY處理的結果。
從圖中可以看出,BERNESE方案存在著一個比較明顯的系統誤差,這也可以從上面的統計看出,而GIPSY方案和GAMIT方案則吻合得比較好,其中,GAMIT方案更與IGS吻合。我們認為,之所以BERNESE方案出現這個系統誤差,主要原因是在實現ITRF2000參考框架統一的時候引起的,由于本研究中BERNESE方案采取的是利用網內的IGS站的約束來進行參考架的7參數轉換,而GIPSY和GAMIT方案均是利用全球均勻分布的IGS核心站做為約束來進行參考架的7參數轉換[7],IT RF參考框架是全球范圍的參考框架,僅僅依靠少量的幾個IGS站的約束得到的參考框架肯定是不穩定和不可信的,這些站的資料質量不好引起的誤差都很容易影響到定義的參考架。相對BERNESE方案,GIPSY方案和GAMIT方案的參考架的統一更合理些,但它們之間也存在著毫米級的差別,這個是可以接受的,畢竟首先兩種軟件本身采取的是不同的算法結構,另外,在參考框架定義的時候,雖然都選取全球分布的IGS站做約束,但可能選的具體站稍有差別,所以存在4~5個mm的差別是可以接受的。

我們在統一了參考框架定義的情況下,發現BERNESE的結果有了很大的改善,見圖 2,以TSKB站為例,在三種軟件統一了參考框架定義的情況下,全部采用不固定任何參考站,和全球的IGS站一起進行方差平差,并定義參考框架ITRF2000,從圖中可以看出,三個方向上三種軟件的結果吻合得非常好,誤差都在4~5 mm以內。其中BERNESE和GAMIT的結果更接近,這也說明了在參考框架以及處理方案完全一致的情況下,不同軟件得到的結果的差異和軟件處理的算法有關系,因為BERNESE和GAMIT均是采用雙差模糊度解算,而GIPSY采用的是直接處理非差載波相位觀測量。
雖然處理的資料的時間尺度才一年,我們還是求出了一個大致的速度,我們列出了BERNESE方案和GAMIT方案得到的速度場,圖3列出了東亞地區13個測站在ITRF2000框架下的速度場,其中紅色的代表BERNESE方案的結果,藍色的則是 GAMIT的結果。從圖中可以看出,BERNESE與GAMIT的結果大致上是吻合的,由于用來求速度的資料時間跨度才一年,很多季節性的影響都無法消除[8],兩者速度上的差異是可以接受的,這說明了BERNESE方案雖然得到的絕對坐標有很大的系統差,但在求速度場的時候影響不是很大。

圖3 GAMIT和BERNESE速度比較圖
GPS數據處理在一定程度上是一項復雜且個性化的工作,即使是同一軟件,在數據的取舍,參考框架的約束,處理方案的制定,不同的操作人員之間都有細微的差別,從而會導致結果不同。然而,同一軟件之間的相互比較只能說明處理方案一致的情況下,軟件的計算過程的穩定性,不同軟件之前的比較,則能比較客觀的評價不同處理結果。通過這次研究,客觀地比較了利用三種軟件制定的三種處理方案的處理結果,確定了一個合理的固定參考框架的方案,采用此方案,GAMIT和GIPSY的結果是吻合的,驗證了方案的合理性。導致結果存在差異的最主要的因素還是參考架框架的歸算引起的,由于在GPS數據處理中,還沒有一個嚴格統一的參考框架的定義辦法,所以,不同的研究單位在處理GPS數據時,所采取的參考框架的定義也會有細微的差別,從而導致結果的差別。下一步,我們將考慮如何對三種軟件的合理結果進行融合,以便得到一個更加可靠穩定的計算結果。
致 謝:感謝韓國國立天文與空間科學研究院Pil-Ho Park博士和日本東京大學地震研究所Teruyuki Kato教授提供韓國和日本的GPS觀測數據。
[1] 賴錫安,黃立人.中國大陸現今地殼運動[M].北京:地震出版社,2004.
[2] 許厚澤,熊 熊.東北亞大陸地殼運動的GPS研究[J].大地測量與地球動力學,2004,24(4):1-6.
[3] KING R W,BOCK Y.Documentation forthe GAMIT GPS Analysis Software[M].America:M ass Inst of Technol.,Cambridge Mass,2002.
[4] ERIC C,MAT HILDE V.GPS measurements of crustal deformation in the Baikal-M ongolia area(1994-2002)[J].Journal of Geophysical Research,2003(8):14-1;14-13.
[5] WEBB F H,ZUM BERG J F.An introduction to GIPSY/OASIS2II precision software for analysis of data from Global Positioning System[M].America,Jet Propulsion Laboratory,California Institute of Technology,1993.
[6] HUGEN TOBLER U,SCHAER S,FRODEZ P,et al.Documentation for Bernese GPS Software Version 4.2[M].German,2001.
[7] BLEWIT T G,BOUCHER C,DAVIES P B H,et al.ITRF Densification and Continuous Realization by the IGS[R].Rio de Janeiro,IAG97,1997.
[8] 陳俊平,王解先.中國境內IGS站數據處理及地殼運動研究[J].同濟大學學報,2005,33(10):1414-1417.