劉　敏，黃謨濤，歐陽(yáng)永忠，鄧凱亮，翟國(guó)君，吳太旗

1. 信息工程大學(xué)地理空間信息學(xué)院，河南 鄭州 450001； 2. 海軍海洋測(cè)繪研究所，天津 300061； 3. 海軍工程大學(xué)導(dǎo)航工程系，湖北 武漢 430033

?

顧及地形效應(yīng)的重力向下延拓模型分析與檢驗(yàn)

劉敏1，黃謨濤2,3，歐陽(yáng)永忠2，鄧凱亮2，翟國(guó)君2,3，吳太旗2

1. 信息工程大學(xué)地理空間信息學(xué)院，河南 鄭州 450001； 2. 海軍海洋測(cè)繪研究所，天津 300061； 3. 海軍工程大學(xué)導(dǎo)航工程系，湖北 武漢 430033

Foundationsupport：TheNationalBasicResearchProgramofChina(973Program) (No. 613219)；TheNationalNaturalScienceFoundationofChina(Nos. 41474012； 41174062； 41374018)；TheGreatScientificInstrumentDevelopmentProjectofChina(No.2011YQ12004503)

向下延拓是航空重力測(cè)量數(shù)據(jù)實(shí)際應(yīng)用中必不可少的技術(shù)環(huán)節(jié)。向下延拓屬于不適定反問(wèn)題，其解算過(guò)程具有較大的不確定性，故該問(wèn)題一直是大地測(cè)量領(lǐng)域國(guó)內(nèi)外學(xué)者的研究熱點(diǎn)。本文深入分析研究了當(dāng)前國(guó)內(nèi)外最具代表性的3種向下延拓計(jì)算模型的技術(shù)特點(diǎn)和適用條件，提出了應(yīng)用超高階位模型、局部地形改正和移去—恢復(fù)技術(shù)顧及地形效應(yīng)，以及位場(chǎng)延拓結(jié)果球面化曲面的工程化方法，重點(diǎn)探討了計(jì)算模型的穩(wěn)定性及數(shù)據(jù)觀測(cè)誤差對(duì)延拓計(jì)算結(jié)果的影響。通過(guò)理論分析、數(shù)值仿真和實(shí)測(cè)數(shù)據(jù)計(jì)算等手段，定量評(píng)估了不同向下延拓模型的解算精度及其可靠性。其主要結(jié)論是：傳統(tǒng)逆Poisson積分模型解嚴(yán)重受制于輸入數(shù)據(jù)觀測(cè)噪聲的干擾，在現(xiàn)有作業(yè)條件下，該模型至多只能用于1km以下高度的延拓解算；頻譜截?cái)喾e分和位模型加地改兩種延拓新模型具有良好的計(jì)算穩(wěn)定性，完全適用于2′分辨率和5km飛行高度條件下的航空重力測(cè)量數(shù)據(jù)向下延拓解算，其延拓計(jì)算精度可達(dá)2×10-5m/s2，可滿(mǎn)足各方面實(shí)際應(yīng)用需求。

航空重力測(cè)量；向下延拓；地形效應(yīng)；Poisson積分；有限帶寬頻譜；超高階位模型

精化大地水準(zhǔn)面始終是現(xiàn)代物理大地測(cè)量學(xué)永恒的研究主題[1-2]。精密確定大地水準(zhǔn)面需要聯(lián)合利用衛(wèi)星、航空、地面、海洋等多源重力和地形高數(shù)據(jù)，其計(jì)算過(guò)程涉及多源數(shù)據(jù)的融合處理問(wèn)題。向下延拓是多源重力數(shù)據(jù)融合處理必不可少的技術(shù)環(huán)節(jié)，在數(shù)據(jù)準(zhǔn)備階段，需要將衛(wèi)星和航空重力測(cè)量成果延拓到地面作聯(lián)合處理；在邊值問(wèn)題解算階段，需要將地面重力數(shù)據(jù)延拓到大地水準(zhǔn)面。由于物理場(chǎng)向下延拓在數(shù)學(xué)上屬于不適定反問(wèn)題，求解此類(lèi)問(wèn)題存在很大的不確定性[3]，一直以來(lái)國(guó)內(nèi)外諸多學(xué)者為解決這一棘手問(wèn)題付出了不懈的努力[4-28]。目前解決向下延拓問(wèn)題主要有3種途徑，第1種是直接求逆泊松(Poisson)積分方程，是應(yīng)用比較廣泛的主流途徑[29]，國(guó)內(nèi)外學(xué)者主要圍繞求逆過(guò)程引起的不穩(wěn)定性問(wèn)題開(kāi)展不同形式的正則化方法研究[10,13,16,18-23,25-26]；第2種是間接求逆途徑，包括最小二乘配置方法[30-32]、虛擬點(diǎn)質(zhì)量方法[1-2]和矩諧分析法[24]等，此類(lèi)方法雖然避開(kāi)了求逆Poisson方程，但仍涉及矩陣求逆過(guò)程，因此也不可避免存在不穩(wěn)定性問(wèn)題[33-34]；第3種可統(tǒng)稱(chēng)為非求逆途徑，文獻(xiàn)[17]提出的利用航空重力測(cè)量和DEM確定地面重力場(chǎng)的直接代表法，文獻(xiàn)[27—28]提出的聯(lián)合使用超高階位模型和地形信息確定向下延拓改正數(shù)的方法(以下簡(jiǎn)稱(chēng)差分延拓法)等都屬于非求逆途徑的范疇，此類(lèi)方法不受傳統(tǒng)求逆過(guò)程不穩(wěn)定性的影響。文獻(xiàn)[29]基于帶限(band-limited)航空重力測(cè)量數(shù)據(jù)特有的頻譜特性，提出了向下延拓的直接積分公式(以下簡(jiǎn)稱(chēng)頻譜截?cái)喾?，并將其推廣應(yīng)用于大地水準(zhǔn)面的直接解算[35]，其計(jì)算過(guò)程也避開(kāi)了方程求逆問(wèn)題。文獻(xiàn)[13]將頻譜截?cái)喾椒ê突谇竽孢^(guò)程的各類(lèi)正則化方法作了全面的數(shù)值比較和分析，得出的結(jié)論是，前者的計(jì)算精度和效率都明顯優(yōu)于其他方法。

如前所述，困擾向下延拓問(wèn)題的關(guān)鍵是其本身固有的不適定性。文獻(xiàn)[9]曾對(duì)向下延拓的穩(wěn)定性問(wèn)題做過(guò)深入分析，同時(shí)提出了改善計(jì)算穩(wěn)定性的地形效應(yīng)補(bǔ)償方法；文獻(xiàn)[11，13，14，29，35]專(zhuān)門(mén)針對(duì)航空重力向下延拓問(wèn)題，開(kāi)展了大量有效的數(shù)值計(jì)算和比對(duì)工作，得出了一些極具參考價(jià)值的結(jié)論。由已有的研究成果得知，雖然向下延拓問(wèn)題本身是不適定的，但其解算結(jié)果的不穩(wěn)定程度取決于延拓計(jì)算高度和數(shù)據(jù)分辨率(即網(wǎng)格間距大小)兩個(gè)方面。在一定條件下，向下延拓解算方程可以是良態(tài)的，而在超越一定界限以后，即使采用正則化處理技術(shù)，也無(wú)法取得有效的解算結(jié)果[11,13]。因此，對(duì)于實(shí)際應(yīng)用來(lái)說(shuō)，最要緊的是預(yù)先掌控前面所指的一定“條件”和“界限”，盡可能在規(guī)定的“條件”下開(kāi)展航空重力測(cè)量作業(yè)，避免超越已知的“界限”，只有這樣才能獲得預(yù)期的測(cè)量成果。另一方面，地形效應(yīng)對(duì)重力向下延拓解算結(jié)果具有重要影響，因?yàn)榍竽孢^(guò)程的穩(wěn)定性不僅取決于積分方程系數(shù)矩陣的結(jié)構(gòu)，還取決于觀測(cè)向量的頻譜特性[9]，故可通過(guò)地形效應(yīng)的“移去-恢復(fù)”運(yùn)算來(lái)改變重力觀測(cè)量的頻譜特征，從而改善向下延拓解算的穩(wěn)定性。本文的目的是，通過(guò)理論分析、模擬仿真和實(shí)際數(shù)值計(jì)算等手段，對(duì)當(dāng)前國(guó)內(nèi)外最具代表性的3種向下延拓模型進(jìn)行全面的分析比較和適用性檢驗(yàn)，力爭(zhēng)為工程應(yīng)用提供具有可操作性的延拓計(jì)算方案。

1　計(jì)算模型及穩(wěn)定性分析

1.1逆Poisson積分迭代解延拓模型

由文獻(xiàn)[6]知，Poisson積分向上延拓公式為

(1)

(2)

向下延拓計(jì)算是求式(1)的逆問(wèn)題，即已知Δgp，要求Δgd。對(duì)式(1)作離散化處理，可得一線性方程組，用矩陣形式表示為

y=Ax

(3)

式中，y為由航空重力測(cè)量獲得的已知重力異常向量；x為待求的球面重力異常向量；A為由式(1)積分核函數(shù)確定的系數(shù)矩陣。可采用如下的Jacobi迭代法解算方程式(3)[9]，將A改寫(xiě)為

A=E-B

(4)

代入式(3)有

x=y+Bx

(5)

取

x0=y

(6)

求

xi=Bxi-1(i>0)

(7)

(8)

式中，E為單位矩陣；ε為給定的某個(gè)大于零的限差值；k為迭代次數(shù)。

在實(shí)際應(yīng)用中，由于觀測(cè)數(shù)據(jù)有限，一般將球面積分區(qū)域劃分為近區(qū)和遠(yuǎn)區(qū)，近區(qū)是以計(jì)算點(diǎn)為中心、ψ0為半徑的球冠區(qū)域σ0，近區(qū)影響直接由觀測(cè)數(shù)據(jù)計(jì)算；遠(yuǎn)區(qū)是球面上的剩余部分(σ-σ0)，遠(yuǎn)區(qū)影響由位系數(shù)模型按下式計(jì)算[11,13]

(9)

(10)

(11)

式中，GM為地球引力常數(shù)；L為位系數(shù)模型的最高階數(shù)；l為參考場(chǎng)位模型的最高階數(shù)；a為參考橢球長(zhǎng)半軸；Δgn代表重力異常n階面球諧函數(shù)；Qn為Poisson核截?cái)嘞禂?shù)；Rn,m(ψ0)可由已知的遞推公式計(jì)算，具體參見(jiàn)文獻(xiàn)[36—37]。此時(shí)，對(duì)應(yīng)于式(3)的離散形式為

(12)

式中，N為待求點(diǎn)個(gè)數(shù)，一般取與已知點(diǎn)個(gè)數(shù)相同。為了減小中心數(shù)據(jù)塊(即與計(jì)算點(diǎn)重合的數(shù)據(jù)塊)離散化誤差的影響，可按下式計(jì)算系數(shù)矩陣A的對(duì)角線元素[13]

(13)

矩陣A的非對(duì)角線元素為

(14)

式中，N0為位于積分球冠區(qū)σ0內(nèi)的數(shù)據(jù)點(diǎn)個(gè)數(shù)；Δσj為數(shù)據(jù)塊面積。

如前所述，式(1)只是近似的球面解模型，式(1)只能反解得到等高度面(即球面或近似為大地水準(zhǔn)面)上的重力異常，而非地形面上的重力異常。如果是向下延拓到大地水準(zhǔn)面上，那么由于大地水準(zhǔn)面外地形質(zhì)量的存在，使得其延拓解與物理意義上的現(xiàn)實(shí)性不相對(duì)應(yīng)，即這里的延拓解只是一組虛擬的重力異常[6]。這類(lèi)重力異?？蓡为?dú)應(yīng)用于物理大地測(cè)量參數(shù)計(jì)算，但不宜與其他類(lèi)型重力異常數(shù)據(jù)聯(lián)合使用。為了求得地形面上的重力異常，提出如下位場(chǎng)延拓球面化曲面方法：

(1) 將計(jì)算區(qū)域的地形高度變化范圍(Hmax-Hmin)，按整百米間隔(如100 m或200 m)劃分為若干個(gè)等高度面；

(2) 利用前面的迭代法依次計(jì)算出各個(gè)等高度面上的重力異常；

(3) 通過(guò)內(nèi)插方法計(jì)算得到位于兩個(gè)等高面之間的地形面上的重力異常。

此外，為了改善向下延拓解算的穩(wěn)定性，可通過(guò)地形效應(yīng)的“移去-恢復(fù)”運(yùn)算改變重力觀測(cè)量的頻譜特征[9]，即首先在航空重力測(cè)量成果中移去地形質(zhì)量對(duì)空間點(diǎn)的作用，然后在延拓計(jì)算結(jié)果中恢復(fù)地形質(zhì)量對(duì)地面點(diǎn)的影響。地形效應(yīng)改正的計(jì)算模型可參見(jiàn)文獻(xiàn)[1，28]。

1.2頻譜截?cái)喾e分延拓模型

由于受到飛行環(huán)境動(dòng)態(tài)效應(yīng)的影響，航空重力測(cè)量原始觀測(cè)數(shù)據(jù)中一般包含上千甚至上萬(wàn)毫伽(mGal=10-5m/s2)的干擾加速度信息[38-39]。為了消除高頻干擾噪聲的影響，通常需要對(duì)航空重力測(cè)量數(shù)據(jù)作低通濾波處理，以獲取所需的重力異常信息。但經(jīng)濾波處理后的重力測(cè)量成果必定會(huì)損失掉一部分有用的高頻重力信息，損失量值大小取決于濾波截?cái)囝l率的選擇，實(shí)際應(yīng)用中需要平衡好測(cè)量精度和分辨率兩方面的需求。由此可見(jiàn)，航空重力測(cè)量成果是經(jīng)過(guò)高頻截?cái)嗪蟮闹亓Ξ惓?chǎng)信息，不包含某個(gè)頻率(L)以上的高頻分量。而在對(duì)航空重力測(cè)量成果作進(jìn)一步處理(如向下延拓)時(shí)，通常還需要引入全球重力場(chǎng)模型(GGM)進(jìn)行移去-恢復(fù)運(yùn)算，即事先需要從航空重力測(cè)量成果中移去GGM的影響，然后在解算結(jié)果中作反向的恢復(fù)運(yùn)算。這說(shuō)明，實(shí)際參與向下延拓解算的航空重力測(cè)量數(shù)據(jù)是一類(lèi)有限帶寬的重力場(chǎng)信息(band-limited airborne gravity data)[29]。設(shè)GGM的最高階次為(l-1)，則有L=l+b，b稱(chēng)為航空重力測(cè)量數(shù)據(jù)的帶寬。如何合理利用航空重力測(cè)量數(shù)據(jù)有限帶寬的頻譜特性，是有效提高向下延拓解算過(guò)程穩(wěn)定性的關(guān)鍵。文獻(xiàn)[9、11、13—14、29、35]全面研究了等高飛行條件下的航空重力測(cè)量數(shù)據(jù)延拓問(wèn)題，提出了如下相對(duì)穩(wěn)定的頻譜截?cái)喾e分向下延拓公式。

設(shè)在恒等高度H的飛行面上獲取了一組有限帶寬航空重力異常數(shù)據(jù)Δgb(R+H)，飛行面與大地水準(zhǔn)面之間已經(jīng)不存在地形質(zhì)量，要求確定大地水準(zhǔn)面(近似為半徑R的球面)上相對(duì)應(yīng)的有限帶寬擾動(dòng)位Tb(R)。理論上可將上述問(wèn)題表述為如下的偽邊值問(wèn)題[9,29]

從表2可以看出，參試品種的株高從72.5-88.0cm之間。對(duì)照(4個(gè)品種平均值)為83.2cm，青海13號(hào)最低，為72.5cm。4個(gè)品種均為直立生長(zhǎng)型。4個(gè)品種的單株有效分枝數(shù)在1.3個(gè)-2.2個(gè)之間，對(duì)照為1.7個(gè)，以青海13號(hào)最高，為2.2個(gè)。4個(gè)品種的單株莢數(shù)在5.2個(gè)-7.5個(gè)之間，對(duì)照為6.0個(gè)，青海13號(hào)最多，為7.5個(gè)。4個(gè)品種的莢粒數(shù)在1.4個(gè)-2.2個(gè)之間，對(duì)照為1.8粒，青海13號(hào)最多，為2.2個(gè)。最少為青海12號(hào)，為1.4個(gè)。4個(gè)品種百粒重在79.8g-170.0g之間，對(duì)照為140.8g，青蠶14號(hào)最高，為170.0g。青海13號(hào)最低，為79.8g。

(15)

(16)

(17)

對(duì)式(16)中的雙邊值面作近似處理并轉(zhuǎn)換為單邊值面后[9]，可求得上述問(wèn)題的擾動(dòng)位解為[29]

Δgb(R+H,Ω′)dΩ′

(19)

Δgb(R+H,Ω′)dΩ′

(20)

(21)

(22)

(23)

式中，L0為位模型最高階次，且L0≤L，其他符號(hào)意義同前。對(duì)于近區(qū)計(jì)算，同樣需要對(duì)中心數(shù)據(jù)塊作精細(xì)化處理，以減小離散化誤差的影響。類(lèi)似于式(13)和式(14)，此時(shí)對(duì)應(yīng)于積分式(20)的系數(shù)矩陣Ab的對(duì)角線元素為

(24)

(25)

[Pn+1(cosψ0)-Pn-1(cosψ0)]/(2n+1)

(26)

矩陣Ab的非對(duì)角線元素為

(27)

式中其他符號(hào)意義同前。由式(20)可直接計(jì)算得到地形面上的重力異常，故無(wú)須作位場(chǎng)延拓平面化曲面處理，此時(shí)數(shù)據(jù)面高度保持不變(r=R+H)，計(jì)算面高度由R改變?yōu)?R+hi)，hi代表地面計(jì)算點(diǎn)的大地高。類(lèi)似于前一小節(jié)的做法，這里同樣需要通過(guò)地形效應(yīng)的“移去-恢復(fù)”運(yùn)算，來(lái)消除計(jì)算面與數(shù)據(jù)面之間的地形質(zhì)量影響，以改善向下延拓解算的穩(wěn)定性。

1.3基于位模型和地形改正差分的延拓模型

為了規(guī)避傳統(tǒng)逆Poisson積分向下延拓解算過(guò)程的不適定性問(wèn)題，文獻(xiàn)[27—28]借鑒導(dǎo)航定位中的“差分”概念，曾先后提出了利用超高階位模型和地形高信息，直接實(shí)施海域和陸部航空重力測(cè)量數(shù)據(jù)向下延拓計(jì)算的新方案。其中，陸部顧及地形效應(yīng)延拓方案的核心思想是：以飛行高度面與地面對(duì)應(yīng)點(diǎn)的位模型重力異常差分信息表征向下延拓總改正數(shù)的中長(zhǎng)波分量，以相對(duì)應(yīng)的局部地形改正差分修正量表征總改正數(shù)的中高頻成分，從而實(shí)現(xiàn)航空重力數(shù)據(jù)向地面點(diǎn)對(duì)點(diǎn)的全頻段延拓。具體延拓計(jì)算模型為[28]

(28)

(29)

(30)

1.4計(jì)算模型穩(wěn)定性分析

如前言所述，重力場(chǎng)向下延拓計(jì)算在數(shù)學(xué)上屬于不適定反問(wèn)題，其解算過(guò)程存在不穩(wěn)定性是該問(wèn)題本身固有的一種屬性[3]。理論上，由逆Poisson積分公式(1)組成的向下延拓解算方程屬于第一類(lèi)弗雷德霍姆(Fredholm)積分方程，當(dāng)且僅當(dāng)滿(mǎn)足如下的Picard條件時(shí)，該方程存在唯一的收斂解[9]

(31)

式中，fi代表已知航空重力觀測(cè)量Δgp的傅里葉(Fourier)展開(kāi)系數(shù)；λi代表Poisson核函數(shù)的特征值。式(31)說(shuō)明，要想確保逆Poisson積分方程有解，從某個(gè)階次開(kāi)始，F(xiàn)ourier系數(shù)fi的衰減速度必須快于特征值λi的衰減速度。

對(duì)于離散化形式的線性方程組即式(3)，通常采用如下的系數(shù)矩陣條件數(shù)指標(biāo)來(lái)度量方程解的病態(tài)性[9,20]

(32)

式中，λmax和λmin分別代表系數(shù)矩陣A之特征值的最大值和最小值。一般認(rèn)為[20]，當(dāng)0<κ<100時(shí)，方程組是良態(tài)的；當(dāng)100≤κ≤1000時(shí)，定義為中等程度病態(tài)；當(dāng)κ>1000時(shí)，定義為嚴(yán)重病態(tài)。對(duì)于在恒等飛行高度H條件下的航空重力向下延拓問(wèn)題，逆Poisson積分方程系數(shù)矩陣的條件數(shù)可近似表示為[9,11]

(33)

式中，π/ΔΩ稱(chēng)為奈奎斯特(Nyquist)頻率，ΔΩ為數(shù)據(jù)網(wǎng)格間距。取ΔΩ=2′,5′，可分別計(jì)算得到對(duì)應(yīng)于不同延拓高度H的條件數(shù)κ，具體結(jié)果如表1所示。必須指出的是，式(33)給出的條件數(shù)只是一種理論上的病態(tài)性度量指標(biāo)，受各種擾動(dòng)因素影響后，實(shí)際問(wèn)題解的變化特性則要復(fù)雜得多，其穩(wěn)定性除了取決于系數(shù)矩陣的結(jié)構(gòu)外，還取決于已知觀測(cè)量的頻譜特性[9]。文獻(xiàn)[11]在完成大量的數(shù)值計(jì)算和分析比較后發(fā)現(xiàn)，航空重力數(shù)據(jù)向下延拓解算的穩(wěn)定性與數(shù)據(jù)網(wǎng)格間距和延拓高度的比值(α=ΔΩ/H)密切相關(guān)，在受到觀測(cè)噪聲干擾的情況下，只有當(dāng)比值α>1.1時(shí)，直接利用航空重力數(shù)據(jù)確定的局部大地水準(zhǔn)面解才是有效的。而航空重力異常直接向下延拓到地面結(jié)果的精度要比前者悲觀得多，即使采用一些正則化方法進(jìn)行優(yōu)化處理，也無(wú)法徹底解決觀測(cè)噪聲的放大問(wèn)題，要想獲得可接受的地面重力異常延拓結(jié)果，必須盡可能降低航空重力測(cè)量的飛行高度[11]。本文將在第2節(jié)繼續(xù)對(duì)此問(wèn)題進(jìn)行數(shù)值計(jì)算和檢驗(yàn)。表1同時(shí)列出了數(shù)據(jù)網(wǎng)格間距分別取ΔΩ=2′,5′時(shí)，對(duì)應(yīng)于不同延拓高度的比值α。

表1　不同延拓高度對(duì)應(yīng)的條件數(shù)和網(wǎng)格間距與高度比值

2　數(shù)值計(jì)算與分析

為了檢驗(yàn)上述各個(gè)計(jì)算模型的延拓效果，分別采用模擬仿真和實(shí)測(cè)數(shù)據(jù)兩種方式開(kāi)展數(shù)值計(jì)算及分析比較研究。

2.1模擬仿真計(jì)算與分析

2.1.1試驗(yàn)數(shù)據(jù)與計(jì)算方案

以美國(guó)本土一個(gè)3°×3°區(qū)塊作為試驗(yàn)區(qū)，選用EGM2008位模型作為標(biāo)準(zhǔn)場(chǎng)模擬產(chǎn)生不同高度的2′×2′網(wǎng)格航空和大地水準(zhǔn)面(即高度為零)重力異常。之所以選擇美國(guó)本土作為試驗(yàn)區(qū)，是考慮到EGM2008位模型在美國(guó)地區(qū)具有更高的逼近度[40-41]，并與3.2節(jié)選用的實(shí)測(cè)數(shù)據(jù)試驗(yàn)區(qū)取得一致。該區(qū)塊屬于地形變化比較劇烈的大山區(qū)，試驗(yàn)效果具有一定的說(shuō)服力。分別以1km、3km和5km3個(gè)高度的位模型重力異常作為觀測(cè)量，依次采用逆Poisson積分(簡(jiǎn)稱(chēng)模型1)和頻譜截?cái)喾e分(簡(jiǎn)稱(chēng)模型2)兩種計(jì)算模型，將3個(gè)不同高度面上的網(wǎng)格重力異常延拓到零高度面，并將其同由位模型計(jì)算得到的零高度面觀測(cè)量(作為基準(zhǔn)值)進(jìn)行比較，從而獲取相應(yīng)的計(jì)算模型精度評(píng)估參數(shù)。由于基于位模型和地形改正差分的延拓模型(簡(jiǎn)稱(chēng)模型3)本身就是建立在超高階位模型之上的，因此該模型暫不參加本階段的數(shù)值計(jì)算檢驗(yàn)。

為了進(jìn)一步考察觀測(cè)噪聲對(duì)延拓計(jì)算結(jié)果的影響，特別設(shè)計(jì)在模擬觀測(cè)量中分別加入±1mGal和±3mGal的隨機(jī)噪聲，生成對(duì)應(yīng)高度面上的帶噪聲的兩組觀測(cè)量，并重復(fù)前面的計(jì)算過(guò)程和對(duì)比分析。上述誤差量值與當(dāng)前國(guó)內(nèi)外航空重力測(cè)量的精度水平相當(dāng)[11,14,29,35,42]，因此其檢驗(yàn)結(jié)論具有實(shí)用意義。本試驗(yàn)統(tǒng)一采用EGM2008位模型的前360階次作為參考場(chǎng)(GGM)，對(duì)應(yīng)于零高度面和3個(gè)不同高度面的EGM2008位模型(361～2160階次)殘差重力異常統(tǒng)計(jì)結(jié)果見(jiàn)表2。

表2　試驗(yàn)區(qū)位模型殘差重力異常統(tǒng)計(jì)結(jié)果

2.1.2計(jì)算結(jié)果與分析

本試驗(yàn)在一個(gè)高度面上共有網(wǎng)格觀測(cè)數(shù)據(jù)30×3×30×3=8100個(gè)，對(duì)應(yīng)于模型1，需要求解8100階線性方程組，采用Jacobi迭代法解算時(shí)，迭代終止參數(shù)取ε=0.1mGal；對(duì)于模型2，模型參數(shù)分別取為l=361、L=L0=2160、b=1799；積分半徑統(tǒng)一取為ψ0=1°。按照前面設(shè)計(jì)的計(jì)算方案，對(duì)3個(gè)高度和兩個(gè)模型分別解算無(wú)誤差和有誤差干擾條件下的零高度面延拓重力異常，同時(shí)將其與零高度面基準(zhǔn)值作比較，具體計(jì)算結(jié)果如表3所示。為了減小積分邊緣效應(yīng)對(duì)統(tǒng)計(jì)結(jié)果的影響，計(jì)算區(qū)域邊緣1°范圍內(nèi)的數(shù)據(jù)不參加對(duì)比分析。

從表3結(jié)果可以看出，對(duì)于無(wú)誤差干擾條件下的輸入數(shù)據(jù)，模型1和模型2都能給出比較理想的輸出結(jié)果，基本不受計(jì)算高度大小的影響。但對(duì)于有誤差干擾條件下的輸入數(shù)據(jù)，模型1和模型2計(jì)算效果則表現(xiàn)出非常顯著的差異。在3km以下計(jì)算高度，模型2對(duì)數(shù)據(jù)誤差具有一定的抑制作用，只有當(dāng)延拓高度增大到5km時(shí)，模型2才對(duì)數(shù)據(jù)誤差產(chǎn)生一定的放大效應(yīng)，這一結(jié)果充分體現(xiàn)了模型2作為正解模型及其截?cái)嗪撕瘮?shù)所具有的超強(qiáng)的抗干擾能力。模型1的解算結(jié)果則明顯受到數(shù)據(jù)誤差的干擾，在所有計(jì)算高度上，模型1對(duì)數(shù)據(jù)誤差都有放大作用，計(jì)算結(jié)果嚴(yán)重失真，當(dāng)數(shù)據(jù)誤差大于1mGal，計(jì)算高度超過(guò)1km時(shí)，模型1的解算結(jié)果都是不可靠的。這一結(jié)果說(shuō)明，雖然在理想情況下，通過(guò)迭代計(jì)算求解逆Poisson積分方程，也能收斂到問(wèn)題的理論解[6,9]。但當(dāng)觀測(cè)數(shù)據(jù)存在噪聲干擾時(shí)，由于反問(wèn)題自身固有的不適定性，觀測(cè)噪聲將在迭代過(guò)程中得到累積和放大，使得較小的噪聲干擾也會(huì)引起問(wèn)題解出現(xiàn)較大的變化，最終造成解算結(jié)果嚴(yán)重偏離真解。

表3　不同模型計(jì)算結(jié)果與基準(zhǔn)值比對(duì)

2.2地面實(shí)測(cè)數(shù)據(jù)計(jì)算與分析

2.2.1試驗(yàn)數(shù)據(jù)與計(jì)算方案

檢驗(yàn)向下延拓計(jì)算模型適用性最有效的方法是，直接使用實(shí)測(cè)航空重力數(shù)據(jù)完成向下延拓計(jì)算，將計(jì)算結(jié)果與地面實(shí)測(cè)重力基準(zhǔn)值進(jìn)行比較，由此可得到不同計(jì)算模型的精度評(píng)價(jià)指標(biāo)。遺憾的是，陸部同時(shí)擁有高精度高分辨率航空和地面實(shí)際重力觀測(cè)數(shù)據(jù)的區(qū)域并不多見(jiàn)，目前筆者還缺乏這方面的可靠資料。為此，本文改用向上與向下延拓比對(duì)方法對(duì)計(jì)算模型精度進(jìn)行外部檢核，即首先利用地面網(wǎng)格重力和地形高數(shù)據(jù)，通過(guò)向上延拓方法計(jì)算得到一定高度面上的空中重力異常，將其作為航空重力測(cè)量的觀測(cè)量，進(jìn)而使用不同的計(jì)算模型將其向下延拓到地面，求取延拓計(jì)算值與地面已知網(wǎng)格重力值的差異，便可獲得相應(yīng)延拓計(jì)算模型的精度評(píng)價(jià)。

這里繼續(xù)選用與2.1節(jié)完全相同的美國(guó)本土3°×3°區(qū)塊作為試驗(yàn)區(qū)，開(kāi)展航空重力向下延拓的實(shí)際數(shù)值計(jì)算和對(duì)比分析。該區(qū)塊同時(shí)擁有2′×2′網(wǎng)格地面觀測(cè)重力異常和30″×30″網(wǎng)格地形高數(shù)據(jù)，兩組數(shù)據(jù)的變化特征如表4所示。

表4　試驗(yàn)區(qū)塊實(shí)測(cè)重力和地形數(shù)據(jù)變化特征

使用“先向下后向上”方法[6,43]，首先將地面重力向下延拓到大地水準(zhǔn)面，然后將其向上延拓到5km高度面。為了考察地形效應(yīng)對(duì)向下延拓計(jì)算結(jié)果的影響，這里采用兩種途徑完成地面重力向下延拓解算，得到兩組大地水準(zhǔn)面上的重力值，一種途徑是直接采用原始觀測(cè)重力異常Δg0(1)和式(3)進(jìn)行純粹數(shù)學(xué)意義上的向下延拓計(jì)算，得到一組虛擬的重力異常Δg*(1)，此方法不擾亂外部重力場(chǎng)[6]；另一種途徑是首先從地面重力中扣除掉大地水準(zhǔn)面以外地形質(zhì)量引力的影響，得到調(diào)整后的重力異常Δg0(2)，然后將其向下延拓到大地水準(zhǔn)面，得到一組消除地形效應(yīng)影響后的延拓重力異常Δg*(2)。最后利用Δg*(1)和Δg*(2)可由式(1)向上延拓得到對(duì)應(yīng)5km高度面的兩組航空重力觀測(cè)量Δgp(1)和Δgp(2)，并根據(jù)模型1和模型2完成后續(xù)的向下延拓計(jì)算。但此時(shí)需要確定的是地形面上的重力異常，故向下延拓計(jì)算高度面是對(duì)應(yīng)于地面網(wǎng)格點(diǎn)高度的等高面。對(duì)于模型1，需要利用前面提出的位場(chǎng)延拓球面化曲面方法，將球面上的計(jì)算結(jié)果內(nèi)插到地面上的各個(gè)測(cè)點(diǎn)；模型2可直接計(jì)算得到測(cè)點(diǎn)上的重力異常；模型3通過(guò)局部地形改正差分方式顧及了地形效應(yīng)的影響，故只需參加由Δg0(1)生成的觀測(cè)量檢驗(yàn)?？紤]到本文關(guān)注的重點(diǎn)是計(jì)算模型的適用性和地形效應(yīng)的影響，不必過(guò)分追求計(jì)算結(jié)果的絕對(duì)精度，故可對(duì)向上和向下延拓模型參數(shù)設(shè)置作統(tǒng)一的近似處理。這里將大地水準(zhǔn)面近似為球面，統(tǒng)一取ri=R+hi，rp=R+H，hi為地面點(diǎn)高程，其他符號(hào)意義同前。此時(shí)模型2的計(jì)算參數(shù)取為l=361，L=5400，b=5039,L0=2160，即與2′×2′數(shù)據(jù)分辨率相一致，其他模型參數(shù)保持不變。

2.2.2計(jì)算結(jié)果與分析

按照前面設(shè)計(jì)的計(jì)算方案，分別采用上述3個(gè)計(jì)算模型完成5km高度面各兩組觀測(cè)數(shù)據(jù)的向下延拓解算，將其分別與相對(duì)應(yīng)的地面基準(zhǔn)值Δg0(1)和Δg0(2)作比較，具體計(jì)算結(jié)果如表5所示。

表5　不同模型5 km高度向下延拓結(jié)果與基準(zhǔn)值比對(duì)

表5結(jié)果顯示，模型1的檢核效果反而明顯好于模型2和模型3，這是由于作為觀測(cè)量的空中重力異常是由地面重力異常通過(guò)Poisson積分向上延拓得到的，它與模型1使用的逆Poisson積分形成閉環(huán)關(guān)系，其計(jì)算模型誤差有相互抵消作用，因此該比對(duì)結(jié)果不能作為評(píng)價(jià)模型一計(jì)算性能的依據(jù)，但相應(yīng)結(jié)果至少說(shuō)明本文提出的位場(chǎng)延拓球面化曲面方法和顧及地形效應(yīng)策略是可行有效的。由于模型2和模型3的計(jì)算原理與Poisson積分向上延拓過(guò)程無(wú)關(guān)，故表5中對(duì)應(yīng)于模型2和模型3的檢核結(jié)果是獨(dú)立有效的。從表5看出，當(dāng)不顧及地形效應(yīng)影響時(shí)，模型2的檢核精度為±2.63mGal；而當(dāng)使用移去-恢復(fù)方法顧及地形影響時(shí)，模型2的檢核精度提高到±1.73mGal?？紤]到使用Poisson積分向上延拓解作為觀測(cè)量可能存在1～2mGal的數(shù)據(jù)誤差[44]，與新型航空重力測(cè)量系統(tǒng)精度水平基本相當(dāng)[42]，對(duì)照表3的仿真檢驗(yàn)結(jié)果，不難看出，表5給出的實(shí)測(cè)數(shù)據(jù)檢核結(jié)果完全符合預(yù)期，同時(shí)說(shuō)明顧及地形效應(yīng)對(duì)提高模型2的計(jì)算精度具有顯著成效。模型3的計(jì)算輸入量與航空重力測(cè)量觀測(cè)量無(wú)關(guān)，計(jì)算過(guò)程穩(wěn)定可靠，且模型自身就具備顧及地形效應(yīng)的功能，因此具有較高的計(jì)算精度，達(dá)到±1.83mGal，與顧及地形效應(yīng)后的模型2計(jì)算精度水平相當(dāng)。

3　結(jié)　論

根據(jù)前面的模型適用性分析和數(shù)值計(jì)算比對(duì)結(jié)果，可得出以下基本結(jié)論：

(1) 基于逆Poisson積分的傳統(tǒng)向下延拓模型明顯受到數(shù)據(jù)分辨率、計(jì)算高度和觀測(cè)噪聲等多種因素的制約，即使在比較理想的數(shù)據(jù)精度條件下，最大延拓高度也只能達(dá)到1km。使用正則化方法可在一定程度上提高傳統(tǒng)延拓模型的穩(wěn)定性，但其解算過(guò)程仍存在較多的不確定性因素，不利于該模型的推廣應(yīng)用。

(2) 基于頻譜截?cái)喾e分的直接向下延拓模型具有良好的計(jì)算穩(wěn)定性，對(duì)高頻觀測(cè)噪聲干擾具有很好的抑制作用。在當(dāng)前航空重力測(cè)量通常采用的2′分辨率和3km飛行高度條件下，該模型幾乎可以不損失觀測(cè)數(shù)據(jù)精度的能力恢復(fù)地球表面的重力異常，因此具有良好的應(yīng)用前景。在實(shí)際應(yīng)用中需要注意把握的問(wèn)題是積分核函數(shù)截?cái)嚯A數(shù)與航空重力測(cè)量數(shù)據(jù)濾波截止頻率的匹配關(guān)系。解算結(jié)果存在一定的邊緣效應(yīng)，同時(shí)觀測(cè)數(shù)據(jù)網(wǎng)格化處理給測(cè)線布設(shè)帶來(lái)的三維空間約束等因素是影響該模型推廣應(yīng)用的主要障礙。

(3) 基于位模型和地形改正差分的直接向下延拓模型具有良好的應(yīng)用前景，其計(jì)算過(guò)程完全獨(dú)立于航空重力觀測(cè)數(shù)據(jù)，因此不受觀測(cè)噪聲的干擾。其計(jì)算精度取決于超高階位模型和局部地形改正差分的相對(duì)精度，在現(xiàn)有技術(shù)條件下，該模型延拓計(jì)算精度可達(dá)2mGal，與當(dāng)前航空重力測(cè)量精度水平相當(dāng)。與模型2相比較，不需要對(duì)觀測(cè)數(shù)據(jù)作網(wǎng)格化預(yù)處理，也不存在積分邊緣效應(yīng)，可實(shí)現(xiàn)空中和地面點(diǎn)對(duì)點(diǎn)延拓計(jì)算是該模型的主要優(yōu)勢(shì)。但關(guān)于該模型物理意義的嚴(yán)密解釋還有待進(jìn)一步的研究，并通過(guò)更多的實(shí)際觀測(cè)數(shù)據(jù)計(jì)算來(lái)驗(yàn)證其應(yīng)用效果。

[1]李建成, 陳俊勇, 寧津生, 等. 地球重力場(chǎng)逼近理論與中國(guó)2000似大地水準(zhǔn)面的確定[M]. 武漢: 武漢大學(xué)出版社, 2003.

LIJiancheng,CHENJunyong,NINGJinsheng,etal.TheoryoftheEarth’sGravityFieldApproximationandDeterminationofChinaQuasi-geoid2000[M].Wuhan:WuhanUniversityPress, 2003.

[2]黃謨濤, 翟國(guó)君, 管錚, 等. 海洋重力場(chǎng)測(cè)定及其應(yīng)用[M]. 北京: 測(cè)繪出版社, 2005.

HUANGMotao,ZHAIGuojun,GUANZheng,etal.TheDeterminationandApplicationofMarineGravityField[M].Beijing:SurveyingandMappingPress, 2005.

[3]王彥飛. 反演問(wèn)題的計(jì)算方法及其應(yīng)用[M]. 北京: 高等教育出版社, 2007.WANGYanfei.ComputationalMethodsforInverseProblemsandTheirApplications[M].Beijing:HigherEducationPress, 2007.

[4]PELLINENLP.AccountingforTopographyintheCalculationofQuasigeoidalHeightsandPlumb-lineDeflectionsfromGravityAnomalies[J].BulletinGéodésique, 1962, 63(1): 57-65.

[5]BJERHAMMARA.GravityReductiontoaSphericalSurface[R].Stockholm:RoyalInstituteofTechnology, 1962.

[6]HEISKANENWA,MORITZH.PhysicalGeodesy[M].SanFrancisco:FreemanWHandCompany, 1967.

[7]JEKELIC.TheDownwardContinuationofAerialGravimetricDatawithoutDensityHypothesis[J].BulletinGéodésique, 1987, 61(4): 319-329.

[8]WANGYM.TheEffectofTopographyontheDeterminationoftheGeoidUsingAnalyticalDownwardContinuation[J].BulletinGéodésique, 1990, 64(3): 231-246.

[9]MARTINECZ.StabilityInvestigationsofaDiscreteDownwardContinuationProblemforGeoidDeterminationintheCanadianRockyMountains[J].JournalofGeodesy, 1996, 70(11): 805-828.

[10]XUPeiliang.TruncatedSVDMethodsforDiscreteLinearIll-posedProblems[J].GeophysicalJournalInternational, 1998, 135(2): 505-514.

[11]NOVKP,KERNM,SCHWARZKP.NumericalStudiesontheHarmonicDownwardContinuationofBand-limitedAirborneGravity[J].StudiaGeophysicaetGeodaetica, 2001, 45(4): 327-345.

[12]MUELLERF,MAYER-GUERRT.ComparisonofDownwardContinuationMethodsofAirborneGravimetryData[C]∥ProceedingsoftheInternationalAssociationofGeodesyIAGGeneralAssemblySapporo.Berlin:Springer, 2005, 128: 254-258.

[13]KERNM.AnAnalysisoftheCombinationandDownwardContinuationofSatellite,AirborneandTerrestrialGravityData[D].Calgary:UniversityofCalgary, 2003.

[14]ALBERTSB,KLEESR.AComparisonofMethodsfortheInversionofAirborneGravityData[J].JournalofGeodesy, 2004, 78(1-2): 55-65.

[15]TZIAVOSIN,ANDRITSANOSVD,FORSBERGR,etal.NumericalInvestigationofDownwardContinuationMethodsforAirborneGravityData[C]∥GGSM2004IAGInternationalSymposiumPorto.BerlinHeidelberg:Springer, 2005, 129: 119-124.

[16]ALBERTSB.RegionalGravityFieldModelingUsingAirborneGravimetryData[C]∥PublicationsonGeodesy70.Delft:NetherlandsGeodeticCommission, 2009.

[17]石磐, 王興濤. 利用航空重力測(cè)量和DEM確定地面重力場(chǎng)[J]. 測(cè)繪學(xué)報(bào), 1997, 26(2): 117-121.

SHIPan,WANGXingtao.DeterminationoftheTerrainSurfaceGravityFieldUsingAirborneGravimetryandDEM[J].ActaGeodaeticaetCartographicaSinica, 1997, 26(2): 117-121.

[18]沈云中, 許厚澤. 不適定方程正則化算法的譜分解式[J]. 大地測(cè)量與地球動(dòng)力學(xué), 2002, 22(3): 10-14.

SHENYunzhong,XUHouze.SpectralDecompositionFormulaofRegularizationSolutionforIll-posedEquation[J].JournalofGeodesyandGeodynamics, 2002, 22(3): 10-14.

[19]王興濤, 石磐, 朱非洲. 航空重力測(cè)量數(shù)據(jù)向下延拓的正則化算法及其譜分解[J]. 測(cè)繪學(xué)報(bào), 2004, 33(1): 33-38.

WANGXingtao,SHIPan,ZHUFeizhou.RegularizationMethodsandSpectralDecompositionfortheDownwardContinuationofAirborneGravityData[J].ActaGeodaeticaetCartographicaSinica, 2004, 33(1): 33-38.

[20]王振杰. 測(cè)量中不適定問(wèn)題的正則化解法[M]. 北京: 科學(xué)出版社, 2006.

WANGZhenjie.RegularizationMethodsforIll-posedProbleminSurvey[M].Beijing:SciencePress, 2006.

[21]顧勇為, 歸慶明. 航空重力測(cè)量數(shù)據(jù)向下延拓基于信噪比的正則化方法的研究[J]. 測(cè)繪學(xué)報(bào), 2010, 39(5): 458-464.

GUYongwei,GuiQingming.StudyofRegularizationBasedonSignal-to-noiseIndexinAirborneGravityDownwardtotheEarthSurface[J].ActaGeodaeticaetCartographicaSinica, 2010, 39(5): 458-464.

[22]蔣濤, 李建成, 王正濤, 等. 航空重力向下延拓病態(tài)問(wèn)題的求解[J]. 測(cè)繪學(xué)報(bào), 2011, 40(6): 684-689.JIANGTao,LIJiancheng,WANGZhengtao,etal.SolutionofIll-posedProbleminDownwardContinuationofAirborneGravity[J].ActaGeodaeticaetCartographicaSinica, 2011, 40(6): 684-689.

[23]鄧凱亮, 黃謨濤, 暴景陽(yáng), 等. 向下延拓航空重力數(shù)據(jù)的Tikhonov雙參數(shù)正則化法[J]. 測(cè)繪學(xué)報(bào), 2011, 40(6): 690-696.

DENGKailiang,HUANGMotao,BAOJingyang,etal.TikhonovTwo-parameterRegulationAlgorithminDownwardContinuationofAirborneGravityData[J].ActaGeodaeticaetCartographicaSinica, 2011, 40(6): 690-696.

[24]蔣濤, 黨亞民, 章傳銀, 等. 基于矩諧分析的航空重力向下延拓[J]. 測(cè)繪學(xué)報(bào), 2013, 42(4): 475-480.

JIANGTao,DANGYamin,ZHANGChuanyin,etal.DownwardContinuationofAirborneGravityDataBasedonRectangularHarmonicAnalysis[J].ActaGeodaeticaetCartographicaSinica, 2013, 42(4): 475-480.

[25]孫文, 吳曉平, 王慶賓, 等. 航空重力數(shù)據(jù)向下延拓的波數(shù)域迭代Tikhonov正則化方法[J]. 測(cè)繪學(xué)報(bào), 2014, 43(6): 566-574.DOI: 10.13485/j.cnki.11-2089.2014.0080.

SUNWen,WUXiaoping,WANGQingbin,etal.WaveNumberDomainIterativeTikhonovRegularizationMethodforDownwardContinuationofAirborneGravityData[J].ActaGeodaeticaetCartographicaSinica, 2014, 43(6): 566-574.DOI: 10.13485/j.cnki.11-2089.2014.0080.

[26]劉曉剛, 李迎春, 肖云, 等. 重力與磁力測(cè)量數(shù)據(jù)向下延拓中最優(yōu)正則化參數(shù)確定方法[J]. 測(cè)繪學(xué)報(bào), 2014, 43(9): 881-887.DOI: 10.13485/j.cnki.11-2089.2014.0160.

LIUXiaogang,LIYingchun,XIAOYun,etal.OptimalRegularizationParameterDeterminationMethodinDownwardContinuationofGravimetricandGeomagneticData[J].ActaGeodaeticaetCartographicaSinica, 2014, 43(9): 881-887.DOI: 10.13485/j.cnki.11-2089.2014.0160.

[27]黃謨濤, 歐陽(yáng)永忠, 劉敏, 等. 海域航空重力測(cè)量數(shù)據(jù)向下延拓的實(shí)用方法[J]. 武漢大學(xué)學(xué)報(bào)(信息科學(xué)版), 2014, 39(10): 1147-1152.HUANGMotao,OUYANGYongzhong,LIUMin,etal.PracticalMethodsforDownwardContinuationofAirborneGravityDataintheSeaArea[J].GeomaticsandInformationScienceofWuhanUniversity, 2014, 39(10): 1147-1152.

[28]黃謨濤, 寧津生, 歐陽(yáng)永忠, 等. 聯(lián)合使用位模型和地形信息的陸區(qū)航空重力向下延拓方法[J]. 測(cè)繪學(xué)報(bào), 2015, 44(4): 355-362.DOI: 10.11947/j.AGCS.2015.20130751.HUANGMotao,NINGJinsheng,OUYANGYongzhong,etal.DownwardContinuationofAirborneGravimetryonLandUsingGeopotentialModelandTerrainInformation[J].ActaGeodaeticaetCartographicaSinica, 2015, 44(4): 355-362.DOI: 10.11947/j.AGCS.2015.20130751.

[29]NOVKP,HECKB.DownwardContinuationandGeoidDeterminationBasedonBand-limitedAirborneGravityData[J].JournalofGeodesy, 2002, 76(5): 269-278.

[30]MORITZH.AdvancedPhysicalGeodesy[M].Karlsruhe:HerbertWichmannVerlag, 1980.

[31]HWANGC,HSIAOYS,SHIHHC,etal.GeodeticandGeophysicalResultsfromaTaiwanAirborneGravitySurvey:DataReductionandAccuracyAssessment[J].JournalGeophysicalResearch, 2007, 112(B4): 93-101.

[32]FORSBERGR,OLESENAV.AirborneGravityFieldDetermination[M]∥XUGuochang.SciencesofGeodesy-I.Berlin:Springer-Verlag, 2010: 83-103.

[33]黃謨濤, 歐陽(yáng)永忠, 翟國(guó)君, 等. 融合多源重力數(shù)據(jù)的Tikhonov正則化配置法[J]. 海洋測(cè)繪, 2013, 33(3): 6-12.

HUANGMotao,OUYANGYongzhong,ZHAIGuojun,etal.TikhonovRegularizationCollocationforMulti-sourceGravityDataFusionProcessing[J].HydrographicSurveyingandCharting, 2013, 33(3): 6-12.

[34]黃謨濤, 歐陽(yáng)永忠, 劉敏, 等. 融合海域多源重力數(shù)據(jù)的正則化點(diǎn)質(zhì)量方法[J]. 武漢大學(xué)學(xué)報(bào)(信息科學(xué)版), 2015, 40(2): 170-175.

HUANGMotao,OUYANGYongzhong,LIUMin,etal.RegularizationofPoint-massModelforMulti-sourceGravityDataFusionProcessing[J].GeomaticsandInformationScienceofWuhanUniversity, 2015, 40(2): 170-175.

[35]NOVKP,KERNM,SCHWARZKP,etal.OnGeoidDeterminationfromAirborneGravity[J].JournalofGeodesy, 2003, 76(9-10): 510-522.

[36]PAULMK.AMethodofEvaluatingtheTruncationErrorCoefficientsforGeoidalHeight[J].BulletinGéodésique, 1973, 110(1): 413-425.

[37]陸仲連. 球諧函數(shù)[R]. 鄭州: 測(cè)繪學(xué)院, 1984.LUZhonglian.SphericalHarmonicFunction[R].Zhengzhou:SurveyingandMappingInstitute, 1984.

[38]OLESENAV.ImprovedAirborneScalarGravimetryforRegionalGravityFieldMappingandGeoidDetermination[D].Copenhagen:UniversityofCopenhagen, 2002.

[39]孫中苗. 航空重力測(cè)量理論、方法及應(yīng)用研究[D]. 鄭州:信息工程大學(xué), 2004.

SUNZhongmiao.Theory,MethodsandApplicationsofAirborneGravimetry[D].Zhengzhou:InformationEngineeringUniversity, 2004.

[40]PAVLISNK,HOLMESSA,KENYONSC,etal.TheDevelopmentandEvaluationoftheEarthGravitationalModel2008 (EGM2008)[J].JournalGeophysicalResearch, 2012, 117(B4): 531-535.

[41]章傳銀, 郭春喜, 陳俊勇, 等.EGM2008地球重力場(chǎng)模型在中國(guó)大陸適用性分析[J]. 測(cè)繪學(xué)報(bào), 2009, 38(4): 283-289.

ZHANGChuanyin,GUOChunxi,CHENJunyong,etal.EGM2008andItsApplicationAnalysisinChineseMainland[J].ActaGeodaeticaetCartographicaSinica, 2009, 38(4): 283-289.

[42]歐陽(yáng)永忠, 鄧凱亮, 陸秀平, 等. 多型航空重力儀同機(jī)測(cè)試及其數(shù)據(jù)分析[J]. 海洋測(cè)繪, 2013, 33(4): 6-11.

OUYANGYongzhong,DENGKailiang,LUXiuping,etal.TestsofMulti-typeAirborneGravimetersandDataAnalysis[J].HydrographicSurveyingandCharting, 2013, 33(4): 6-11.

[43]MORITZH.TheComputationoftheExternalGravityFieldandtheGeodeticBoundary-ValueProblem[C]∥ProceedingsoftheSymposiumontheExtensionofGravityAnomaliestoUnsurveyedAreas.GeophysicalMonographSeriesNumber9. [S.l.]:AmericanGeophysicalUnionoftheNationalAcademyofSciences, 1966: 127-136.

[44]ARGESEANUVS.UpwardContinuationofSurfaceGravityAnomalyData[C]∥ProceedingsoftheIAGSymposiumonAirborneGravityFieldDetermination.Boulder: [s.n.], 1995: 95-102.

(責(zé)任編輯：陳品馨)

TestandAnalysisofDownwardContinuationModelsforAirborneGravityDatawithRegardtotheEffectofTopographicHeight

LIUMin1,HUANGMotao2,3，OUYANGYongzhong2，DENGKailiang2,ZHAIGuojun2,3，WUTaiqi2

1.InstituteofGeospacialInformation,InformationEngineeringUniversity,Zhengzhou450001,China； 2.NavalInstituteofHydrographicSurveyingandCharting，Tianjin300061，China； 3.DepartmentofNavigation,NavalUniversityofEngineering,Wuhan430033,China

Downwardcontinuationisanessentialtechnicalstepofdataprocessinginairbornegravimetryforfurtherapplications.Itisknownthatthesolutionofdownwardcontinuationisuncertainduetoitsill-posedness.Soithasbeenatopicofgeneralinterestformanyscholarsathomeandabroadingeodesy.Themainpurposeofthispaperistogive3representativemodelsfordownwardcontinuationincludingtraditionalinversePoissonintegrationandtwomodernmethods,andmakeacomprehensivecomparisonandanalysisontheirpropertyandapplicabilityamongthedifferentmodels.Ultra-high-degreegeopotentialmodel,localtopographiccorrectionandremove-restoretechniquearesuggestedtobeusedforregardtotheeffectoftopographicheight,andfortherealizationofdownwardcontinuationcombiningwithatransformationfromsphericaltoundulatingsurface.Wepayourattentiontotheinfluenceofsurveyeddataerrorsonthestabilityofdownwardcontinuationsolutions.Theoreticalanalysis,simulateddataandrealnumericalcomputationsarecarriedouttoevaluatethestabilityandaccuracyofdownwardcontinuationmodels.Andsomeusefulconclusionsareobtained.Underexistingworkingconditions,thetraditionalinversePoissonintegrationmethodcanonlybeusedtothecontinuationcomputationunder1kmduetotheseriousdisturbingofsurveyingnoise.Excellentcomputationstabilitycanbeachievedbyusingtheband-limitedspectrumandthegeopotentialmodelplustopographiccorrectionmethods.Thetwonewmodelscanbeusedtothedownwardcontinuationofairbornegravitydataon5kmheightand2′dataresolution.Andtheaccuracyofcorrespondingcontinuationsolutionscanbereach2×10-5m/s2.Itcanmeettherequirementsfromdifferentapplications.

airbornegravimetry;downwardcontinuation;terraineffect;Poissonintegration;band-limitedspectrum;ultra-high-degreegeopotentialmodel

LIUMin,HUANGMotao，OUYANGYongzhong，etal.TestandAnalysisofDownwardContinuationModelsforAirborneGravityDatawithRegardtotheEffectofTopographicHeight[J].ActaGeodaeticaetCartographicaSinica,2016,45(5):521-530.DOI:10.11947/j.AGCS.2016.20150453.

P223

A

1001-1595(2016)05-0521-10

國(guó)家973計(jì)劃(613219)；國(guó)家自然科學(xué)基金(41474012；41174062；41374018)；國(guó)家重大科學(xué)儀器設(shè)備開(kāi)發(fā)專(zhuān)項(xiàng)(2011YQ12004503)

引文格式：劉敏，黃謨濤，歐陽(yáng)永忠，等.顧及地形效應(yīng)的重力向下延拓模型分析與檢驗(yàn)[J].測(cè)繪學(xué)報(bào)，2016,45(5)：521-530.DOI:10.11947/j.AGCS.2016.20150453.