自適應(yīng)卡爾曼濾波在航空重力異常解算的應(yīng)用研究

鄭崴， 張貴賓

1 中國地質(zhì)大學(xué)(北京)地球物理與信息技術(shù)學(xué)院， 北京　100083 2 河南科技學(xué)院信息工程學(xué)院， 河南新鄉(xiāng)　453003

自適應(yīng)卡爾曼濾波在航空重力異常解算的應(yīng)用研究

鄭崴1,2， 張貴賓1*

1 中國地質(zhì)大學(xué)(北京)地球物理與信息技術(shù)學(xué)院， 北京100083 2 河南科技學(xué)院信息工程學(xué)院， 河南新鄉(xiāng)453003

摘要依據(jù)航空重力測量基本原理，構(gòu)建了航空重力異常解算的卡爾曼濾波模型，將新息自適應(yīng)卡爾曼濾波器(IAE，Innovation based Adaptive Estimation)應(yīng)用于量測噪聲未知的航空重力異常解算.針對IAE濾波器滑動窗口寬度難以準(zhǔn)確確定的問題，通過對多個不同滑動窗口新息協(xié)方差估計的加權(quán)平均，獲得改進(jìn)的IAE濾波器，該IAE濾波器不僅具有量測噪聲自適應(yīng)估計能力，還能實現(xiàn)滑動采樣窗口的優(yōu)化選取.試驗結(jié)果表明，IAE濾波器可以降低因量測噪聲統(tǒng)計信息不明引起的解算誤差，改進(jìn)IAE解算的重力異常誤差約為1 mGal.

關(guān)鍵詞航空重力異常解算； 自適應(yīng)卡爾曼濾波； 新息； 量測噪聲

1引言

航空重力測量是以飛機為載體，測量近地空中重力異常的重力測量技術(shù)，它綜合運用了航空重力儀，全球定位系統(tǒng)GPS(Global Position System)，慣性導(dǎo)航系統(tǒng)INS(Inertial Navigation System)，無線電等技術(shù)設(shè)備(張昌達(dá)，2005；熊盛青，2009a；王靜波等，2009，2010).航空重力測量數(shù)據(jù)經(jīng)過后處理，便可得到所需的航空重力異常，這個數(shù)據(jù)后處理的過程稱為航空重力異常解算(簡稱為航重解算).近年來，差分GPS技術(shù)的進(jìn)步使得航空重力測量精度大幅度提高，解算所得重力異常精度大部分可達(dá)到1～2 mGal(1 mGal=10-5m·s-2)，最高的已達(dá)0.5 mGal(郭志宏等，2008；熊盛青，2009b; 熊盛青等，2010；Cai et al.，2013).

航空重力測量是在高速飛行運動中進(jìn)行的，測量數(shù)據(jù)必然會受到包括飛機運動產(chǎn)生的擾動加速度的干擾，干擾噪聲的幅度最高可達(dá)到106mGal以上(孫中苗，2004)，遠(yuǎn)遠(yuǎn)高于幾十到幾百毫伽的重力異常信號.航空重力測量數(shù)據(jù)的噪聲主要體現(xiàn)為高頻干擾，依據(jù)測量數(shù)據(jù)的頻域特點設(shè)計低通濾波器可消除噪聲影響提取所需的重力異常信號.目前，能夠?qū)崿F(xiàn)線性相位的FIR(Finite Impulse Response)濾波器是航重解算普遍采用的低通濾波器，利用FIR低通濾波器可以獲得精度約1 mGal左右的航空重力異常(郭志宏等，2007, 2011).FIR濾波器的顯著優(yōu)點是具有成熟的設(shè)計和實現(xiàn)方法，但由于重力信號與噪聲在頻譜上沒有明顯的分界線，為獲得高精度的解算輸出，需設(shè)計高階的低通濾波器以分離重力異常信號與噪聲(蔡劭琨等，2010；Cai et al.，2010；羅峰等，2012)，這加重了邊界效應(yīng)對解算結(jié)果的影響.

部分學(xué)者致力于利用卡爾曼濾波器代替數(shù)字低通濾波器實現(xiàn)航空重力異常的解算(Bolotin et al.，2003；Bolotin and Popelensky，2007；Bolotin and Golovan，2013).卡爾曼濾波器通常是采用狀態(tài)空間模型在時間域?qū)ο到y(tǒng)的運動狀態(tài)進(jìn)行描述，因此卡爾曼濾波器是基于系統(tǒng)運動模型的濾波方法.狀態(tài)空間模型的應(yīng)用將不可觀測的變量(狀態(tài)變量)并入可觀測的模型，并與其一起處理得到相應(yīng)的估計輸出，所以卡爾曼濾波可實現(xiàn)多維信號的并行處理(王靜波，2010；王靜波等，2011).在進(jìn)行航重解算時，還能夠?qū)娇罩亓x校準(zhǔn)參數(shù)進(jìn)行實時估計(Zhou and Cai，2013；Zheng et al.，2014).此外，航重解算是在線下進(jìn)行的，可對卡爾曼濾波的估計結(jié)果進(jìn)行平滑處理，獲得高精度的航空重力異常輸出(潘炎冰，2012；李瑞，2014).但是，卡爾曼濾波能夠獲得準(zhǔn)確解算輸出的重要前提，是能夠掌握系統(tǒng)噪聲和量測噪聲的統(tǒng)計信息，而這一前提往往是難以滿足的.解決這一問題的辦法就是用自適應(yīng)卡爾曼濾波器(AKF, Adaptive Kalman filter)代替?zhèn)鹘y(tǒng)的標(biāo)準(zhǔn)卡爾曼濾波器(SK, Standard Kalman filter)進(jìn)行航重解算.基于新息協(xié)方差估計的自適應(yīng)卡爾曼濾波器(IAE, Innovation-based Adaptive Estimation)已證明可以提高噪聲統(tǒng)計特性不明時的估計精度(Mohamed，1999；Mohamed and Schwarz，1999)，IAE的優(yōu)點在于實現(xiàn)方法簡單且能夠?qū)α繙y噪聲協(xié)方差進(jìn)行實時估計，但其新息協(xié)方差估計器的滑動采樣窗口寬度缺乏明確的確定方法.

航重解算中，卡爾曼濾波所需的量測量由GPS測定的位置或速度數(shù)據(jù)提供.由于測量過程中飛機運動狀態(tài)變化劇烈且具有不規(guī)律性，量測噪聲參數(shù)難以準(zhǔn)確預(yù)設(shè)，因此這里著重討論利用IAE濾波器實現(xiàn)量測噪聲自適應(yīng)的航重解算方法.本文將依據(jù)航空重力測量的基本原理，建立航空重力異常解算的卡爾曼濾波模型，設(shè)計適用于G型航空重力測量系統(tǒng)的IAE濾波器.另外，針對IAE濾波器的滑動窗口寬度難以準(zhǔn)確確定的問題，采用多個不同滑動窗口寬度的新息協(xié)方差估計器，通過對多個新息協(xié)方差估計的加權(quán)平均，獲得滑動窗口寬度優(yōu)化的改進(jìn)IAE濾波器.最后，將IAE濾波器和改進(jìn)的IAE濾波器用于實測數(shù)據(jù)的航重解算試驗.

2新息自適應(yīng)卡爾曼濾波器及改進(jìn)

2.1新息自適應(yīng)卡爾曼濾波器

設(shè)在時刻k，離散線性系統(tǒng)的狀態(tài)方程和量測方程可分別表示為

(1)

(2)

Xk為k時刻的狀態(tài)變量；Φk,k-1表示由k-1時刻到k時刻的一步轉(zhuǎn)移矩陣；Γk-1為系統(tǒng)噪聲驅(qū)動陣；Zk為時刻k的量測量；Hk為量測陣；Wk-1為系統(tǒng)激勵噪聲序列；Vk為量測噪聲序列.Qk是系統(tǒng)噪聲的協(xié)方差陣；Rk為量測噪聲協(xié)方差陣，則標(biāo)準(zhǔn)卡爾曼濾波的算法為(秦永元等，2012)

狀態(tài)一步預(yù)測：

(3)

一步預(yù)測誤差：

(4)

濾波增益：

(5)

狀態(tài)估計：

(6)

估計均方誤差：

(7)

新息(Innovation)定義為k時刻濾波器的量測估計值與實際量測值之差，它表示由量測量帶來的、濾波器無法預(yù)測的那部分信息，即

(8)

新息序列的理論協(xié)方差為Crk，則有

(9)

當(dāng)系統(tǒng)噪聲和量測噪聲服從不相關(guān)的高斯分布且新息序列滿足各態(tài)遍歷性時，可根據(jù)新息序列的滑動平均獲得新息協(xié)方差的極大似然估計，有(MohamedandSchwarz，1999；岳曉奎和袁建平，2005；卞鴻巍等，2006)

(10)

(11)

則量測噪聲協(xié)方差陣的自適應(yīng)估計為

(12)

2.2滑動采樣窗口寬度的優(yōu)化

實際應(yīng)用中，如建立的系統(tǒng)模型不準(zhǔn)確，不能真實反映物理過程時，易導(dǎo)致模型與獲得的量測值不匹配，將會出現(xiàn)濾波發(fā)散，這時需通過增大當(dāng)前時刻量測信息在狀態(tài)估計中的加權(quán)系數(shù)來抑制濾波發(fā)散(秦永元等，2012).

圖1　滑動采樣窗口優(yōu)化Fig.1　Structure of optimization for sliding sampling window

(13)

(14)

(15)

可見，相比于普通的IAE濾波器，改進(jìn)的IAE濾波器配備有多個不同滑動采樣窗口寬度的新息協(xié)方差估計器，通過對不同采樣窗口的新息協(xié)方差估計進(jìn)行加權(quán)平均，允許新息協(xié)方差估計器的滑動采樣窗口在一定范圍內(nèi)動態(tài)變化，一方面可以減弱不當(dāng)?shù)幕瑒硬蓸哟翱谠斐傻慕馑阏`差，另一方面也可以提高濾波器工作的穩(wěn)定性，其估計精度劣于取最優(yōu)滑動采樣窗口的情況，但優(yōu)于其他滑動采樣窗口寬度的估計結(jié)果.多個新息方差估計器的結(jié)構(gòu)如圖1所示.

2.3 RTS(Rauch-Tung-Striebel)平滑器

(16)

(17)

(18)

3航空重力異常解算的卡爾曼濾波模型

G型航空重力測量系統(tǒng)是三軸穩(wěn)定平臺型測量系統(tǒng)，它包括一個具有水平穩(wěn)定作用的慣性平臺(INS)，一個重力傳感器(航空重力儀)，以及機載GPS接收器(移動站)和地面GPS接收器(接收站).本部分將根據(jù)航空重力儀運動方程及GPS測量方程構(gòu)建航重解算的卡爾曼濾波模型.

3.1航空重力測量的基本原理

航空重力測量的基本原理是從航空重力儀測得的比力數(shù)據(jù)(飛機運動加速度與重力加速度之和)中扣除飛機運動加速度的影響，獲得所需要的重力加速度(熊盛青等，2010).根據(jù)牛頓第二定律，航空重力儀(垂直加速度計)傳感器探測質(zhì)量(SM，Sensitive Mass)在地理垂直方向上的動力學(xué)方程可表示為(Bolotin，2009)

(19)

厄特缶斯改正是為了補償因地球自轉(zhuǎn)產(chǎn)生的離心力和飛機運動產(chǎn)生的離心力附加在重力儀測量質(zhì)點造成的測量誤差，可計算為

(20)

uE為地球自轉(zhuǎn)速度；φ為飛機所處地理緯度；vE、vN分別是飛機的東向速度和北向速度；RE、RN分別是地理東向和北向的曲率半徑.

重力正常場和高度改正項G0(φ,h)可由Helmert公式近似計算(Bolotin，2009)：

(21)

式中，G=9.7803 m·s-2，β=0.005302，

β1=-0.000007，β0=-0.00014 m·s-2，

ω0=0.00123 s-1.

不難看出，為獲得式(19)中的重力異常Δg，需要對航空重力儀測量數(shù)據(jù)進(jìn)行一系列的改正.其中，一部分改正項有嚴(yán)密的解析公式，可以直接計算得到，如厄特缶斯改正、正常場和高度改正，它們統(tǒng)稱為解析類改正；而另一類改正無法直接用解析公式來計算，需借助濾波方法來實現(xiàn)，即垂直加速度改正.因此，利用卡爾曼濾波實現(xiàn)重力異常解算實際上是在消除噪聲干擾的同時實現(xiàn)垂直加速度改正的過程.3.2卡爾曼濾波模型的建立

將航空重力儀SM的實時高程，作為卡爾曼濾波的量測方程，寫為

(22)

h′為重力儀SM所處的實時高度；h為SM的實際高度；δh為測量誤差.

航空重力儀及兩個水平加速度計的測量方程可表示為(Bolotin，2009)

(23)

將式(23)代入式(19)聯(lián)立并整理得

(24)

其中

(25)

fΣ為已知的作用在SM上的合力產(chǎn)生的加速度總和；qΣ為系統(tǒng)噪聲項.為實現(xiàn)航重解算，將重力異?？醋餮販y線分布的隨機過程，并用濾波器表示，有

(26)

將式(24)和式(26)合并作為卡爾曼濾波的系統(tǒng)方程，式(22)為卡爾曼濾波的量測方程，并用狀態(tài)空間形式表示為

(27)

(28)

其中，Xf為狀態(tài)向量，重力異常Δg為其分量，即

(29)

Af為動力學(xué)矩陣；Bf為系統(tǒng)噪聲驅(qū)動矩陣；Cf為控制項驅(qū)動矩陣；qf和δh分別為系統(tǒng)噪聲和量測噪聲項；uf為系統(tǒng)的固定控制項.將式(27)和(28)進(jìn)行離散化，再利用設(shè)計的IAE濾波器和RTS平滑器進(jìn)行處理，即可實現(xiàn)重力異常的估計.

4航空重力異常解算試驗

4.1測量基本情況

我們在2007年利用G型航空重力測量系統(tǒng)在某地進(jìn)行了實地測量.實驗所用數(shù)據(jù)為第五架次飛行測量數(shù)據(jù)，飛行時間為4 h.此次測量是為了驗證系統(tǒng)工作的穩(wěn)定性及重復(fù)測量結(jié)果的一致性.圖2為測線位置示意圖，圖中橫軸為相對經(jīng)度，縱軸為相對緯度.由圖可知，此次測量共有6條東西方向的飛行測線，其中1號、2號、3號及6號測線為四條短測線，它們是沿同一位置往返飛行4次的重復(fù)測線；4號測線和5號測線為兩條長飛行測線，其中5號測線相較其他測線有約0.002°的北向偏移.測量選取在氣象條件較好的夜間進(jìn)行，以減小氣流對測量結(jié)果的影響.測量中，設(shè)計飛行高度為400 m，實際飛行測量高度在設(shè)計飛行高度的±15 m范圍內(nèi)波動，圖3為飛行測量高度變化示意圖.此外，飛行測量中的設(shè)計飛行速度為60 m·s-1.

圖2　測線位置示意圖Fig.2　Position diagram of survey lines

圖3　飛行高度示意圖Fig.3　Chart of flight height

剔除實測數(shù)據(jù)中的非測線數(shù)據(jù)，并進(jìn)行解析類改正后，所得結(jié)果(未濾波)如圖4所示.不難看出，航空重力測量信號受噪聲干擾非常嚴(yán)重，幾十到幾百毫伽的重力異常完全被數(shù)萬毫伽的噪聲所淹沒.

關(guān)于解算試驗有以下兩點需要說明，第一，由于暫無法獲得實際地面重力測量數(shù)據(jù)，試驗過程中將航空重力測量系統(tǒng)配套軟件解算結(jié)果作為解算參考標(biāo)準(zhǔn)，將本試驗的解算結(jié)果與參考標(biāo)準(zhǔn)進(jìn)行對比研究；第二，由于RTS平滑后的重力異常精度高于只使用卡爾曼濾波器估計的重力異常精度(潘炎冰，2012；李瑞，2014)，因此后文解算結(jié)果如無特別說明均是指平滑處理后的重力異常.

4.2不同預(yù)設(shè)量測噪聲的航空重力異常解算試驗

將SK解算的重力異常與參考標(biāo)準(zhǔn)進(jìn)行對比，通過反復(fù)調(diào)整量測噪聲和系統(tǒng)噪聲，獲得與解算參考標(biāo)準(zhǔn)誤差最小時的量測噪聲和系統(tǒng)噪聲分別為R0=(0.0153)2和Q0=diag[0,(3.5×10-5)2,0,0,0,0,0,(0.85×10-7)2]，diag[·]為取對角線元素符號.實驗中SK和IAE的系統(tǒng)噪聲設(shè)定為Q0，初始量測噪聲R分別預(yù)設(shè)為R/R0=0.01、0.25、25、100進(jìn)行解算.同時，為確保新息協(xié)方差估計的準(zhǔn)確性，IAE的滑動采樣窗口寬度暫設(shè)為80.

SK和IAE在不同預(yù)設(shè)量測噪聲時，解算的重力異常如圖5所示.為方便分析，將解算所得6條飛行測線的重力異常按時間順序連接組成一條曲線，紅色曲線為解算參考標(biāo)準(zhǔn)，藍(lán)色曲線為SK解算的重力異常曲線，黑色曲線為IAE解算的重力異常曲線.由圖可見，SK解算的重力異常曲線隨著預(yù)設(shè)R的變化出現(xiàn)較大的波動，而且預(yù)設(shè)的R與R0的偏差越大，SK解算重力異常的誤差也越高；而在同樣的量測噪聲條件下，IAE解算的重力異常曲線基本保持穩(wěn)定，即IAE可降低因預(yù)設(shè)量測噪聲不當(dāng)引起的解算誤差.

統(tǒng)計SK和IAE解算輸出與參考標(biāo)準(zhǔn)的誤差，如表1所示.可以看到，在預(yù)設(shè)量測噪聲R=0.01R0時，SK解算結(jié)果的均方誤差高達(dá)3.423367mGal；當(dāng)預(yù)設(shè)量測噪聲R=0.25R0時，SK解算重力異常與參考標(biāo)準(zhǔn)的均方誤差為1.198765mGal.在IAE解算誤差一欄中，不同預(yù)設(shè)量測噪聲條件下，IAE解算重力異常的最大均方誤差為0.923542mGal，最小均方誤差為0.887827mGal，兩者的差距僅為0.035715mGal，而且IAE解算結(jié)果與解算參考標(biāo)準(zhǔn)的誤差都保持在1 mGal以下.所以，IAE濾波器在解算過程中對量測噪聲協(xié)方差進(jìn)行實時的估計和修正，減弱了解算精度對先驗量測噪聲統(tǒng)計信息的依賴.因此，相比基于SK的航重解算方法，利用IAE進(jìn)行航重解算更具有普遍性的意義.

圖4　未濾波航空重力異常Fig.4　Un-filtered airborne gravity anomaly

解算濾波方法比較點數(shù)誤差統(tǒng)計項預(yù)設(shè)量測噪聲(R/R0)0.010.2525100SK19518平均誤差(mGal)3.0827970.9786831.9484592.57978619518均方誤差(mGal)3.4233671.1987652.0226482.652140IAE(N=80)19518平均誤差(mGal)0.6822390.6834810.7147730.71157719518均方誤差(mGal)0.8878270.8891970.9235420.920205

圖5　SK、IAE與解算參考標(biāo)準(zhǔn)解對比圖(a) R/R0=0.01； (b) R/R0=0.25； (c) R/R0=25； (d) R/R0=100.Fig.5　Airborne gravity anomaly of SK、IAE VS reference

4.3不同滑動窗口的航空重力異常解算試驗

由IAE算法步驟可知，滑動采樣窗口寬度N直接關(guān)系到新息協(xié)方差估計的結(jié)果，不當(dāng)?shù)幕瑒硬蓸哟翱谝矔?dǎo)致IAE解算精度的下降.本部分將利用不同滑動窗口的IAE濾波器和改進(jìn)IAE濾波器進(jìn)行解算試驗，驗證改進(jìn)IAE濾波器的優(yōu)化效果，并討論N對IAE解算重力異常精度的影響.試驗中，IAE濾波器的滑動采樣窗口取N=40、80、120、160、200分別進(jìn)行解算；改進(jìn)IAE濾波器的新息協(xié)方差估計器設(shè)定為5個，各新息協(xié)方差估計器的采樣窗口分別為N1=40，N2=80，N3=120，N4=160，N5=200.

圖6為不同滑動窗口的IAE濾波器估計的量測噪聲協(xié)方差.由圖可見，在測線上進(jìn)行飛行測量時，飛機的運動狀態(tài)平穩(wěn)，此時量測噪聲協(xié)方差估計曲線變化平緩；當(dāng)飛機由一條測線轉(zhuǎn)飛向另一條測線時(A、B、C、D、E指示的部分)，飛機的運動狀態(tài)變化劇烈，量測噪聲協(xié)方差估計曲線出現(xiàn)較大的跳動，而且曲線的跳動頻率和幅度隨著滑動窗口寬度的增加而減小，即滑動窗口越寬，IAE對量測噪聲的估計越平穩(wěn)，反之，IAE估計的量測噪聲曲線跳變越劇烈.

統(tǒng)計不同滑動窗口的IAE與改進(jìn)IAE的解算誤差，如表2所示.不難看出，在N=80時，IAE的解算輸出與解算參考標(biāo)準(zhǔn)的誤差最小，此時的滑動采樣窗口寬度兼顧了系統(tǒng)的運動狀態(tài)變化和新息協(xié)方差估計的準(zhǔn)確性；而在其他采樣窗口寬度條件下，IAE的解算結(jié)果誤差均高于N=80的情況.這也印證了前述的滑動窗口寬度選取的不宜太大也不宜太小.

注意到利用改進(jìn)IAE解算的重力異常與參考標(biāo)準(zhǔn)的均方誤差在1.021047 mGal，解算精度劣于N=80時IAE的解算結(jié)果，但優(yōu)于其他的滑動窗口寬度的情況.因此，改進(jìn)的IAE可以減小由于滑動采樣窗口寬度選取不當(dāng)引起的解算誤差.

圖6　IAE估計R變化曲線(N=40、80、120、160、200)Fig.6　Estimated R by IAE (N=40,80,120,160,200)

濾波方法比較點數(shù)最大誤差(mGal)最小誤差(mGal)平均誤差(mGal)均方誤差(mGal)改進(jìn)IAE195184.660951-3.9856210.7810011.021047IAE(N=40)195186.116302-3.9087150.8110201.069920IAE(N=80)195183.872649-3.0882630.6840910.889861IAE(N=120)195184.070119-3.6601320.8671471.063740IAE(N=160)195184.310235-3.9103800.9236401.104894IAE(N=200)195184.645537-3.6042080.9348591.102127

5結(jié)語

本文利用新息自適應(yīng)卡爾曼濾波器及改進(jìn)的新息自適應(yīng)卡爾曼濾波器實現(xiàn)航空重力異常解算，并得出如下結(jié)論：

(1) 新息自適應(yīng)卡爾曼濾波器的應(yīng)用，可以降低量測噪聲協(xié)方差選取不準(zhǔn)確造成的解算誤差，減弱了解算精度對量測噪聲先驗統(tǒng)計信息的依賴性；

(2) 對新息自適應(yīng)卡爾曼濾波器進(jìn)行改進(jìn)，即在新息自適應(yīng)卡爾曼濾波器中設(shè)計多個不同采樣窗口的新息協(xié)方差估計器，實現(xiàn)新息協(xié)方差估計器的滑動采樣窗口優(yōu)化.解算結(jié)果表明，改進(jìn)的新息自適應(yīng)卡爾曼濾波的應(yīng)用能夠減小由滑動采樣窗口選取不當(dāng)造成的解算誤差.

為獲得精度更高、更有普遍意義的航空重力異常解算方法，建議進(jìn)一步加強對自適應(yīng)卡爾曼濾波器設(shè)計方法的研究.同時，還應(yīng)繼續(xù)加強對航空重力測量系統(tǒng)的研究，進(jìn)一步完善航空重力異常解算的卡爾曼濾波模型，提高解算精度.

致謝感謝中國國土資源航空物探遙感中心的熊盛青、周堅鑫、郭志宏和周錫華對本研究的大力支持和幫助.

References

Bian H W, Jin Z H, Wang J P, et al. 2006. The innovation-based estimation adaptive Kalman filter algorithm for INS/GPS integrated navigation system.JournalofShanghaiJiaotongUniversity(in Chinese), 40(6): 1000-1003, 1009.Bolotin Y V, Golovan A A, Parusniov N A. 2003. Mathematical models of airborne gravimetry systems with aided inertial platforms.MoscowUniversityMathematicsBulletin, 58(4):13-23.Bolotin Y V, Popelensky M Y. 2007. Accuracy analysis of airborne gravity when gravimeter parameters are identified in flight.JournalofMathematicalSciences, 146(3):5911-5519.

Bolotin Y V. 2009. Mathematics behind GTGRAV. Moscow: Moscow Lomonosov State University.

Bolotin Y V, Golovan A A. 2013. Methods of inertial gravimetry.MoscowUniversityMechanicsBulletin, 68(5): 117-125.

Cai S K, Wu M P, Zhang K D. 2010. A comparison of digital lowpass FIR-filters in airborne gravimetry.GeophysicalandGeochemicalExploration(in Chinese), 34(1):74-78.

Cai S K, Wu M P, Zhang K D, et al. 2013. The first airborne scalar gravimetry system based on SINS/DGPS in China.ScienceChinaEarthSciences, 56(12): 2198-2208.

Guo Z H, Luo F, An Z F. 2007. Experimental researches on FIR lowpass digital filters based on window functions of airborne gravity data.GeophysicalandGeochemicalExploration(in Chinese), 31(6):568-571, 576.

Guo Z H, Xiong S Q, Zhou J X, et al. 2008. The research on quality evaluation method of test repeat lines in airborne gravity survey.ChineseJournalofGeophysics(in Chinese), 51(5):1538-1543.

Guo Z H, Luo F, Wang M, et al. 2011. The design and experiment of IIR lowpass digital filters for airborne gravity data.ChineseJournalofGeophysics(in Chinese), 54(8):2148-2153.

Li R. 2014. Research of gravity anomaly solution in airborne gravimetry based on Kalman filtering [Master′s thesis] (in Chinese). Beijing: China University of Geosciences Beijing.

Luo F, Guo Z H, Luo Y, et al. 2012. Experimental researches on FIR lowpass filter based on equiripple.GeophysicalandGeochemicalExploration(in Chinese), 36(5): 856-860.

Mohamed A H. 1999. Optimizing the estimation procedure in INS/GPS integration for kinematic applications [Ph. D. thesis]. Calgary: Canada Calgary University.

Mohamed A H, Schwarz K P. 1999. Adaptive Kalman Filtering for INS/GPS.JournalofGeodesy, 73(4): 193-203.

Pan Y B. 2012. The Kalman filter method of free air anomaly solution in airborne gravimetry [Master′s thesis] (in Chinese). Beijing: China University of Geosciences Beijing.

Qin Y Y, Zhang H Y, Wang S H. 2012. Kalman Filter and Principle of Integrated Navigation 2nded (in Chinese). Xi′an: Northwestern Polytechnical University Press.Sun Z M. 2004. Theory, methods and applications of airborne gravimetry [Ph. D. thesis] (in Chinese). Zhengzhou: PLA Information Engineering University.

Tan F J, Xu J N, Li A, et al. 2009. Innovation-based adaptive Kalman filter for accelerometer signal de-noising.JournalofDataAcquisitionandProcessing(in Chinese), 24(2): 227-231.Wang J B, Xiong S Q, Zhou X H, et al. 2009. The advances in the study of the airborne gravimetry system.GeophysicalandGeochemicalExploration(in Chinese), 33(4): 368-373.

Wang J B, Xiong S Q, Guo Z H, et al. 2010. Estimation of the vertical acceleration for the airborne gravimetry using Kalman smoothing.ProgressinGeophysics(in Chinese), 25(3): 968-974, doi: 10.3969/j.issn.1004-2903.2010.03.035.

Wang J B. 2010. Methodologies and technology of data processing for airborne gravimetry [Ph. D. thesis] (in Chinese). Beijing: China University of Geosciences Beijing.

Wang J B, Xiong S Q, Guo Z H, et al. 2011. Research on methods for estimating airborne gravity anomalies.GeophysicalandGeochemicalExploration(in Chinese), 35(4): 493-498.

Xiong S Q. 2009a. The present situation and development of airborne gravity and magnetic survey techniques in China.ProgressinGeophysics(in Chinese), 24(1): 113-117.

Xiong S Q. 2009b. The strategic consideration of the development of China′s airborne geophysical technology.GeologyinChina(in Chinese), 36(6):1366-1374.

Xiong S Q, Zhou X H, Guo Z H, et al. 2010. Theory, Method and Application of Airborne Gravity Prospecting (in Chinese).Beijing: Geology Publishing.

Yue X K, Yuan J P. 2005. An adaptive Kalman filtering algorithm based on maximum-likelihood criterion.JournalofNorthwesternPolytechnicalUniversity(in Chinese), 23(4): 469-474.Zhang C D. 2005. On the subject of airborne gravimetry and airborne gravity gradiometry.ChineseJournalofEngineeringGeophysics(in Chinese), 2(4):282-291.Zheng W, Zhang G B, Li R. 2014. Application of Kalman filtering and smoothing for airborne gravimetry.AppliedMechanicsandMaterials, 548-549:1192-1195.

Zhou W, Cai T J. 2013. Study on filtering methods of airborne gravity.AppliedMechanicsandMaterials, 333-335:516-521.

附中文參考文獻(xiàn)

卞鴻巍, 金志華, 王俊璞等. 2006. 組合導(dǎo)航系統(tǒng)新息自適應(yīng)卡爾曼濾波算法. 上海交通大學(xué)學(xué)報, 40(6): 1000-1003, 1009.

蔡劭琨, 吳美平, 張開東. 2010. 航空重力測量中FIR低通濾波器的比較. 物探與化探, 34(1): 74-78.

郭志宏, 羅鋒, 安戰(zhàn)鋒. 2007. 航空重力數(shù)據(jù)窗函數(shù)法FIR低通數(shù)字濾波試驗. 物探與化探, 31(6): 568-571.

郭志宏, 熊盛青, 周堅鑫等. 2008. 航空重力重復(fù)線測試數(shù)據(jù)質(zhì)量評價方法研究.地球物理學(xué)報, 51(5):1538-1543.

郭志宏, 羅鋒, 王明等. 2011. 航空重力數(shù)據(jù)無限脈沖響應(yīng)低通數(shù)字濾波器設(shè)計與試驗研究. 地球物理學(xué)報, 54(8): 2148-2153.

李瑞. 2014. 基于卡爾曼濾波的航空重力異常解算的研究[碩士論文]. 北京: 中國地質(zhì)大學(xué)(北京).

羅峰, 郭志宏, 駱遙等. 2012. 航空重力數(shù)據(jù)的等波紋FIR低通濾波試驗. 物探與化探, 36(5): 856-860.

潘炎冰. 2012. 航重測量中解算自由空氣異常的卡爾曼濾波方法[碩士論文]. 北京: 中國地質(zhì)大學(xué)(北京).

覃方君, 許江寧, 李安等. 2009. 基于新息自適應(yīng)卡爾曼濾波的加速度計信號降噪. 數(shù)據(jù)采集與處理, 24(2): 227-231.

秦永元, 張洪鉞, 汪叔華. 2012. 卡爾曼濾波與組合導(dǎo)航原理(第二版). 西安: 西北工業(yè)大學(xué)出版社.

孫中苗. 2004. 航空重力測量理論、方法及應(yīng)用研究[博士論文]. 鄭州: 解放軍信息工程大學(xué).

王靜波, 熊盛青, 周錫華等. 2009. 航空重力測量系統(tǒng)研究進(jìn)展. 物探與化探, 33(4): 368-373.

王靜波, 熊盛青, 郭志宏等. 2010. 利用Kalman平滑技術(shù)估算航空重力測量中的載體垂直加速度. 地球物理學(xué)進(jìn)展, 25(3): 968-974, doi: 10.3969/j.issn.1004-2903.2010.03.035.

王靜波. 2010. 航空重力測量數(shù)據(jù)處理方法技術(shù)研究[博士論文]. 北京: 中國地質(zhì)大學(xué)(北京).

王靜波, 熊盛青, 郭志宏等. 2011. 航空重力異常估計方法研究. 物探與化探, 35(4): 493-498.

熊盛青. 2009a. 我國航空重磁勘探技術(shù)現(xiàn)狀與發(fā)展趨勢. 地球物理學(xué)進(jìn)展, 24(1): 113-117.

熊盛青. 2009b. 發(fā)展中國航空物探技術(shù)有關(guān)問題的思考. 中國地質(zhì), 36(6): 1366-1374.

熊盛青, 周錫華, 郭志宏等. 2010. 航空重力勘探理論方法及應(yīng)用. 北京: 地質(zhì)出版社.

岳曉奎, 袁建平. 2005. 一種基于極大似然準(zhǔn)則的自適應(yīng)卡爾曼濾波算法. 西北工業(yè)大學(xué)學(xué)報, 23(4): 469-474.

張昌達(dá). 2005. 航空重力測量和航空重力梯度測量問題. 工程地球物理學(xué)報, 2(4): 282-291.

(本文編輯胡素芳)

Application research on adaptive Kalman filtering for airborne gravity anomaly determination

ZHENG Wei1, 2, ZHANG Gui-Bin1*

1DepartmentofGeophysicsandInformationTechnology,ChinaUniversityofGeosciences,Beijing100083,China2DepartmentofInformationEngineering,HenanInstituteofScienceandTechnology,XinxiangHenan453003,China

AbstractBased on the basic principle of airborne gravimetry, Kalman filtering model for airborne gravity anomaly determination is established and innovation-based adaptive estimation Kalman filter (IAE) is applied to estimate the airborne gravity anomaly when measurement covariance is unknown. At the same time, improved IAE filter installing several innovation covariance estimators is designed to reduce the difficulty of sampling window width selection. It is can be seen from the experimental result that IAE is able to reduce gravity anomaly determination errors caused by lack of statistical information of measurement noise. Moreover, the improved IAE filter can decrease the rate of determination errors introduced by the improper selection of sampling window width, and the error of gravity anomaly determined by improved IAE is nearly 1 mGal.

KeywordsAirborne gravity anomaly determination； Adaptive Kalman filtering； Innovation； Covariance of measurement noise

基金項目國家高技術(shù)研究發(fā)展計劃(863計劃)(2011AA060501)資助.

作者簡介鄭崴，男，1984年生，博士，主要從事航空重力測量技術(shù)的研究. E-mail:13667060@qq.com *通訊作者張貴賓，男，1958年生，教授. 主要從事地球物理反演、重磁電勘探等研究. E-mail:gbzhang@cugb.edu.cn

doi:10.6038/cjg20160410 中圖分類號P223,P631

收稿日期2015-09-14，2015-12-24收修定稿

鄭崴， 張貴賓. 2016. 自適應(yīng)卡爾曼濾波在航空重力異常解算的應(yīng)用研究. 地球物理學(xué)報，59(4):1275-1283,doi:10.6038/cjg20160410.

Zheng W, Zhang G B. 2016. Application research on adaptive Kalman filtering for airborne gravity anomaly determination.ChineseJ.Geophys. (in Chinese),59(4):1275-1283,doi:10.6038/cjg20160410.