航空和地面重力測量數據的一體化處理方法

翟振和，孫中苗，李迎春，郭飛霄

(西安測繪研究所，陜西 西安 710054)

The Integrative Method of Airborne Disturbing Gravity Data and

Land Gravity Anomaly Data

ZHAI Zhenhe，SUN Zhongmiao，LI Yingchun，GUO Feixiao

?

航空和地面重力測量數據的一體化處理方法

翟振和，孫中苗，李迎春，郭飛霄

(西安測繪研究所，陜西 西安 710054)

The Integrative Method of Airborne Disturbing Gravity Data and

Land Gravity Anomaly Data

ZHAI Zhenhe，SUN Zhongmiao，LI Yingchun，GUO Feixiao

摘要：在三維經典擾動位協方差模型基礎上，推導獲得局部重力異常協方差模型、局部擾動重力協方差模型，以及擾動重力和重力異常的互協方差模型，這3個模型在忽略小項后是一致的，即航空擾動重力數據與地面重力異常數據的聯合處理采用一個模型即可，大大簡化了處理流程。以最小二乘配置理論及本文推導的協方差模型為基礎，構建了航空擾動重力測線數據與地面重力異常數據一體化處理的具體方法，利用該方法可推估地面待求點的重力異常數據，且處理過程自動完成了航空擾動重力數據的向下延拓，同時也大大削弱了航空數據中存在的系統偏差影響。在某試驗區開展了一體化處理試驗，結果表明，采用一體化處理得到的重力異常數據比單獨使用5′分辨率航空重力數據推估獲得的數據精度提高2.4 mGal，比單獨使用2.5′分辨率地面重力數據推估獲得的數據精度提高1.2 mGal。

關鍵詞：航空擾動重力數據; 地面重力異常數據；協方差模型；一體化處理

一、引言

目前，利用航空重力測量技術獲取局部區域的高分辨率重力場信息成為一個非常有效的手段。在航空測量數據處理中，傳統的技術流程是先將觀測得到的航空擾動重力數據進行濾波、格網化，然后轉化為重力異常數據，最后向下延拓至對應的地面區域獲得最終的地面重力異常數據[1]。上述數據處理過程往往會使測量數據受到不當的處理(如采用不準確的模型)而影響最終數據的質量(如存在一定的系統誤差)。如果能找到一種直接利用測線擾動重力數據進行處理的方法則可以在減少處理流程的同時提高產品數據的精度。在航空重力測量實施的區域，地面上或多或少會有已知的重力異常測量點，這些測量點的重力數據較之航空測量精度更高，因此，如何將地面已知重力數據有效利用起來構建航空和地面一體化的處理方法成為國內外學者關注的重點。最小二乘配置法具有聯合多種數據類型的處理優勢[2-5]，然而在運用最小二乘配置法進行處理時，多數學者都是在固定的協方差模型及模型參數下進行解算，即所有觀測數據都采用同一個先驗信息，而后按照該先驗信息完成數據的推估和濾波[6-8]。針對上述問題，本文將在最小二乘配置法基礎上，構建一種能夠對航空和地面重力數據進行一體化處理的方法，該方法應顧及不同高度和不同數據類型，且兼具格網化和向下延拓等功能。

二、局部重力異常的協方差模型

構建能夠反映局部重力場空間變化的重力異常協方差模型是進行一體化處理的基礎工作。根據文獻[3]的研究成果，采用擾動位的三維協方差經驗模型如下

(1)

考慮到局部區域的應用環境，將式(1)用大地坐標替換得到如下形式的局部擾動位協方差模型

(2)

式中，f、α具體含義同上式；Hi=Hp+Hq+Di，Hp、Hq代表p、q兩點的大地高。

根據擾動場元間的線性關系，按照協方差傳播定律可得p、q兩點重力異常間協方差為

(3)

式中，CTT為擾動位協方差模型；ρp、ρq分別是p、q兩點的地心向徑。

(4)

對擾動位協方差模型的通用形式CTT求導得

(5)

二階導數為

(6)

將式(5)、式(6)代入式(3)得到局部重力異常協方差的一般形式為

(7)

其完整形式為

(8)

(9)

(10)

(11)

(12)

模型的完整形式為

(13)

式(9)—式(13)即為簡化后的局部重力異常的協方差實用模型。

三、局部擾動重力和重力異常的互協方差模型

擾動重力與擾動位間的關系為

(14)

考慮到橢球面上大地高方向H與n方向相同，因此

(15)

按照協方差傳播定律可得空間擾動重力觀測點p和地面重力異常待求點q的互協方差為

(16)

將局部擾動位協方差模型CTT代入上式并忽略小項得到互協方差模型為

(17)

式(17)即擾動重力與重力異常互協方差模型的完整形式。

同理，局部擾動重力的自協方差模型按照上述原理推導如下

(18)

將式(18)、式(17)與式(13)對比發現，擾動重力與重力異常互協方差模型與重力異常協方差模型、擾動重力協方差模型形式上完全一致，這樣在航空擾動重力數據和地面重力異常數據處理時只采用同一個協方差模型即可。

四、一體化處理方法

假設航空擾動重力數據觀測值為L1，地面重力異常數據觀測值為L2，對應的觀測信號為t1、t2，對應的觀測噪聲分別為Δ1、Δ2，航空重力數據對應的系統性參數為X，此時得到兩組觀測方程分別為

L1=A1X+t1+Δ1

(19)

L2=t2+Δ2

(20)

式中，A1為系數陣；觀測值的協方差陣為

(21)

式中，C11、C22表示兩組數據的自協方差；C21、C12表示兩組數據的互協方差。

未測信號S(地面未測點的重力異常)與已測信號t的互協方差為

(22)

按照經典理論可推導獲得未測信號的推估公式[9-11]如下

(23)

(24)

(25)

式中

(26)

(27)

(28)

(29)

系統參數的估值為

(30)

根據大量的試驗分析和理論推導，采用以下計算步驟。對于兩類數據的自協方差，分別采用各自的參數計算即Δg1的自協方差C11為

(31)

Δg2的自協方差為

(32)

對于兩類數據的互協方差的計算，首先加入兩個調整因子Pa、Pb對模型參數進行調整，為了保證加入調整因子后的參數符合實際情況，Pa、Pb應滿足如下關系

Pa+Pb=1

(33)

利用Pa、Pb對模型參數進行調整得到

Dcom=D1×Pa+D2×Pb

(34)

Ucom=U1×Pa+U2×Pb

(35)

(36)

(37)

(38)

由于未測信號的先驗信息未知，因此未測信號與已測信號的互協方差也按照調整后的模型參數計算

(39)

(40)

在實際應用中，Pa、Pb需要進行迭代計算，迭代循環的截止條件則以每次計算獲得的單位權中誤差作為評判依據，當前后兩次獲得的單位權中誤差之差小于某個限值時，循環結束，如圖1所示。

五、試驗分析

試驗區位于我國陸地某區域，該區域航空重力測量與陸地重力數據分布如圖2所示，整個區域有分辨率為5′×5′的航空擾動重力數據(飛行平均高度3600 m，測線數據采樣頻率為1 Hz)和2.5′分辨率的地面重力異常數據，地形數據采用SRTM 90 m分辨率的地形模型。為了檢驗一體化處理算法，本文假設參與地面重力異常數據ΔgT布設形式為均勻分布，如圖3所示，其中“★”代表參與一體化處理的地面重力數據，其余空格區域代表需要求解的重力數據，同時也用來對一體化處理結果進行檢核。

首先利用航空擾動重力數據進行格網化，然后利用直接代表法向下延拓[12]獲得地面待求點數據。其次利用地面重力數據按照經典的最小二乘配置法推估該區域待求點重力數據，兩種方法分別與已知的地面重力數據進行比較，見表1。

圖1　航空和地面重力數據一體化處理流程

圖2　航空重力測量與陸地重力數據分布

圖3　地面重力數據呈均勻分布

最后按照本文提出的一體化處理方法，利用航空擾動重力的測線數據求解地面待求點重力數據，統計結果見表2，最終處理后的地面重力異常數據如圖4所示。通過表1、表2對比可以看出，一體化處理方法比單獨使用航空和地面重力數據進行推估的效果都要好，這也說明了一體化處理方法對于聯合處理航空和地面重力數據是非常適用的。

表1　航空和地面重力數據分別推估待求點數據統計結果　mGal

表2　利用一體化處理方法得到的結果　mGal

圖4　航空和地面重力數據一體化處理后的三維重力圖

六、結論

本文深入研究了航空擾動重力數據和地面重力異常數據的聯合處理問題，構建了一種一體化處理方法，主要結論和創新點如下：

1) 在三維擾動位協方差模型基礎上，推導獲得局部重力異常協方差模型、局部擾動重力協方差模型，以及擾動重力和重力異常的互協方差模型。這3個模型在忽略小項后是一致的，即航空擾動重力數據與地面重力異常數據的聯合處理采用一個模型即可，大大簡化了處理流程，提高了計算效率。

2) 以最小二乘配置法為基礎，構建了航空擾動重力數據與地面重力異常數據一體化處理的具體方法和實施步驟，利用該方法可推估地面任意待求點的重力數據，并可有效估計航空重力數據中的系統偏差。

3) 在某試驗區開展了一體化處理試驗，結果表明，采用一體化處理獲得的重力異常數據比單獨使用5′分辨率航空重力數據推估獲得的數據精度提高2.4 mGal，比單獨使用2.5′分辨率地面重力數據推估獲得的數據精度提高1.2 mGal，且一體化處理過程自動完成了航空擾動重力數據的向下延拓，同時也大大削弱了航空數據中存在的系統偏差。

本文所構建的航空擾動重力數據和地面重力異常數據的一體化處理方法將為未來我國大規模的航空重力測量及磁力測量數據處理提供很好的參考和借鑒。

參考文獻：

[1]孫中苗.航空重力測量理論、方法及應用研究[D].鄭州: 中國人民解放軍信息工程大學,2004.

[2]OLESEN A V, ANDERSEN O B, TSCHERNING C C. Merging of Airborne Gravity and Gravity Derived from Satellite Altimetry: Test Cases Along the Coast of Greenland[J]. Studia Geophysica et Geodaetica, 2002, 46(3): 387-394.

[3]STRYKOWSKI G, FORSBERG R. Operational Merging of Satellite, Airborne and Surface Gravity Data by Draping Techniques[C]∥Geodesy on the Move: Gravity, Geoid, Geodynamics and Antarctica. Berlin: Springer Verlag, 1997: 243-248.

[4]BAYOUD F A. Some Investigations on Local Geoid Determination from Airborne Gravity Data[D]. Calgary, Canada:University of Calgary,2001:33-35.

[5]HWANG C, GUO J Y, DENG X L, et al. Coastal Gravity Anomalies from Retracked Geosat/GM Altimetry:Improvement, Limitation and the Role of Airborne Gravity Data[J]. Journal of Geodesy, 2006, 80(4): 204-216.

[6]郝燕玲，成怡，劉繁明, 等. 融合多類型海洋重力數據算法仿真研究[J]. 系統仿真學報, 2007, 19(21): 4897-4899.

[7]KERN M, SCHWARZ K K P P, SNEEUW N. A Study on the Combination of Satellite, Airborne, and Terrestrial Gravity Data[J]. Journal of Geodesy, 2003, 77(3-4): 217-225.

[8]楊元喜, 曾安敏. 大地測量數據融合模式及其分析[J]. 武漢大學學報: 信息科學版, 2008, 33(8): 771-774.

[9]楊元喜, 張菊清, 張亮. 基于方差分量估計的擬合推估及其在GIS誤差糾正的應用[J]. 測繪學報, 2008, 37(2): 152-155.

[10]陸仲連. 地球重力場理論與方法[M]. 北京: 解放軍出版社, 1996: 316-320.

[11]鄒賢才.最小二乘配置在物理大地測量中的應用[D]. 武漢: 武漢大學,2003.

[12]王興濤,夏哲仁,石磐, 等.航空重力測量數據向下延拓方法比較[J].地球物理學報,2004,47(6):1017-1022.

本刊開通微信公眾平臺服務賬號

[本刊訊]為了更全面、及時地傳播測繪地理信息學科相關研究的最新成果，方便廣大讀者、作者獲取本刊信息，更大地發揮刊物的作用，本刊自2015年1月起正式開通微信公眾平臺服務賬號，定期發布本刊每期論文目錄、重點文章和新聞報道，歡迎廣大讀者、作者關注。

《測繪通報》微信號：CN11-2246

(本刊編輯部)

引文格式： 翟振和，孫中苗，李迎春，等. 航空和地面重力測量數據的一體化處理方法[J].測繪通報，2015(1)：39-43.DOI:10.13474/j.cnki.11-2246.2015.0007

作者簡介：翟振和(1980—)，男，博士，助理研究員，主要從事物理大地測量方面的研究。E-mail： zhaizhenhe1980@163.com

基金項目：國家自然科學基金(41174017)

收稿日期：2014-07-01

中圖分類號：P223

文獻標識碼：B

文章編號：0494-0911(2015)01-0039-05