基于SINS/GNSS的航空矢量重力測量數據處理方法研究

寧津生

（1.武漢大學測繪學院，武漢 430079；2.武漢大學地球空間環境與大地測量教育部重點實驗室，武漢 430079）

基于SINS/GNSS的航空矢量重力測量數據處理方法研究

寧津生1，2

（1.武漢大學測繪學院，武漢 430079；2.武漢大學地球空間環境與大地測量教育部重點實驗室，武漢 430079）

航空矢量重力測量是未來高效測定中高頻地球重力場信息的先進技術，與地面重力測量、海洋重力測量和衛星重力測量技術相互補充。本文介紹了基于捷聯慣性導航系統/全球衛星導航系統（SINS/ GNSS）的航空矢量重力測量的原理，重點討論了航空矢量重力測量的數據預處理、數據歸算、帶限地形影響、向下延拓及大地水準面確定等問題，為發展我國航空矢量重力測量技術提供參考。

航空矢量重力測量；SINS/GNSS數據融合；地球重力場；大地水準面

1 前言

航空重力測量技術是以飛機為載體，快速測定近地空中重力加速度的重力測量方法。航空重力測量最早出現在20世紀50年代，經過數十年的沉淀與發展，已成為高效測定中高頻地球重力場信息的主要手段。在地面重力測量難以到達的地區進行航空重力測量，可以有效填補困難地區的重力空白，有利于改善重力場數據的精度和分辨率及精化局部重力場模型和區域大地水準面。航空重力測量結合其他技術手段獲取的高精度地球重力場信息，在國家經濟建設、軍事和國防建設、地球物理、海洋學、地球動力學、資源勘探等相關地球科學領域具有重要作用。

以美國、加拿大、丹麥等為代表的發達國家先后開展了航空重力測量研究，并進行了大量的飛行試驗[1～3]，我國也在2002年成功研制了自己的航空標量測量系統Chinese airborne gravimetry system（CHAGS）[4]。眾多研究結果表明，航空標量重力測量數據在波長5～10 km的尺度上能夠達到1～3 mGal（1 Gal=0.01 m/s2）的精度，這也意味著航空標量重力測量技術已經進入成熟階段。航空矢量重力測量技術不僅能獲取航空標量重力測量觀測值（即重力矢量的垂直分量），也能測得重力矢量的水平分量，是目前大地測量領域的研究熱點之一。Jekeli和Garcia[5]驗證了全球衛星導航系統（GNSS）相位加速度應用于航空矢量重力測量的可行性；Kwon和Jekeli[6，7]等率先啟動關于航空矢量重力測量的研究與試驗并提出了利用慣性導航系統（INS）/GNSS組合系統測定矢量重力的新方法，所得重力矢量垂直分量精度為3 mGal，水平分量精度為6～8 mGal，分辨率為10 km；Senobari[8]給出了基于捷聯慣性導航系統（SINS）/GNSS的航空矢量重力測量的試驗結果；美國政府也于2012年開始為期10年的GRAV-D計劃，并將航空矢量重力測量作為其中的一項關鍵技術，預期確定北美地區大地水準面的精度達到2 cm[9]。目前國內還沒有航空矢量重力測量系統，在航空矢量重力數據處理方面也尚處于起步階段。因此，充分吸收國內外現有研究成果，深入研究航空矢量重力測量的基礎理論和數據處理方法，對我國航空矢量重力測量技術的發展具有重要的科學意義和應用價值。

2 基本原理

依據牛頓第二定律，單位質點的加速度在數值上等于質點所受作用力的合力。在慣性系i系中質點的運動方程可寫為

ai=fi+gi（1）

式（1）中，ai為慣性加速度，由GNSS觀測值給出；fi為比力觀測值，由SINS給出；gi為重力加速度。將重力加速度矢量表示為正常重力矢量γn與擾動重力δgn之和，則在導航坐標系n系中，航空矢量重力測量的基本模型可表示為

式（2）中，右邊第三項為Coriolis加速度；Ωnie為地心地固坐標系e系相對于i系的旋轉角速度；Ωnen為n系相對于e系的旋轉角速度。對基于SINS/GNSS的航空矢量重力測量系統來說，加速度計和陀螺儀固定在載體上，比力是在載體坐標系b系中測得的，因此需要利用陀螺儀數據將比力轉換至導航坐標系n系，此時式可寫為

式（3）中，Rnb為b系到n系的旋轉矩陣。

3 數據處理方法

3.1 數據預處理

航空矢量重力測量是一個動態測量過程，各類觀測數據中均包含了大量的觀測噪聲，如何從低信噪比的加速度觀測信息中有效地提取重力信息是航空矢量重力測量數據處理的關鍵。這一過程又大致可分為數據去噪、誤差補償、數據融合、低通濾波等，具體數據處理流程如圖1所示。在數據預處理過程中對原始觀測數據進行去噪和誤差補償，可以提高數據的信噪比，進而提高數據處理精度。慣性測量單元（inertial measurement unit，IMU）是航空矢量重力測量系統的重要組成部分，主要由加速度計和陀螺儀構成。由于IMU自身的誤差特性，在外業測量前需對其進行外場標定，提高數據質量，減小系統誤差對后續數據處理帶來的影響[10]。

圖1 航空矢量重力測量的數據預處理Fig.1 Data preprocessing of airborne vector gravimetry

初始對準是SINS在進行各種應用前所必需的一項重要工作，計算出b系與n系或者i系之間的姿態轉移矩陣，其對準精度對后續的各種應用產生直接影響，而Kalman濾波技術是實現初始對準的一種有效方法。最優的Kalman濾波器是建立在函數模型和噪聲特性均精確已知的條件下，但噪聲統計特性一般難以獲得，自適應Kalman濾波是解決這一應用難題的主要方法之一。根據IMU的自身特性以及作業場地等因素，加速度計的噪聲特性通常是未知的，本文提出了基于簡化自協方差最小二乘（autocovariance least-squares，ALS）噪聲估計的SINS靜基座初始對準算法（SALS），計算流程見圖2。

圖2 基于簡化ALS噪聲估計的SINS靜基座初始對準流程圖Fig.2 Flow chart of the simplified ALS for SINS initial alignment on stationary base

圖3給出了數值仿真中東、北、天方向姿態角的初始對準誤差，從中可以看出，天方向姿態角比水平方向姿態角更易受觀測噪聲協方差矩陣的影響，當先驗觀測噪聲協方差矩陣存在較大偏差時，對天方向對準結果有較大影響；本文提出的基于簡化ALS噪聲估計的SINS SALS同時進行噪聲估計和姿態角修正，因此各方向姿態角的對準精度較常規Kalman濾波方法有著明顯得提高。靜基座初始對準的Kalman濾波模型并非完全可觀，三個失準角由于受陀螺常值漂移和加速度計常值漂移的影響，故均存在估計偏差。

圖3 初始對準誤差Fig.3 Errors of initial alignment

載體的運動加速度是最重要的改正項，其精度直接影響著航空重力測量的精度和空間分辨率[11，12]，但載體加速度并非直接觀測值，它是通過載體位置或速度差分計算得到的。理論上，已知載體位置的時間序列，采用數值差分的方法便可得到載體加速度的時間序列，但由于數值差分過程會放大高頻噪聲，在實際數據處理中，還需進行低通濾波，然而當噪聲和信號在同一頻帶上有交疊時，采用低通濾波方法必然會造成信號損失，且不同的低通濾波器及其參數的選擇也會對信號處理結果有一定的影響。Kalman濾波方法是一種基于模型的數字濾波方法，可以用來確定載體的加速度。在進行航空重力測量時，一般要求飛行過程比較平穩，因此可采用常加速度模型對載體運動狀態進行建模。由于Kalman濾波模型的噪聲協方差矩陣特別是狀態噪聲協方差矩陣難以獲取，不適當的噪聲參數的選取必然會對Kalman濾波的結果帶來影響。而顧及常加速度模型的狀態噪聲協方差矩陣滿足特定結構[13]，本文提出了一種基于ALS噪聲估計的載體加速度確定算法，其計算流程見圖4。

圖4 基于ALS噪聲估計的載體加速度確定算法流程圖Fig.4 Flow chart of the ALS method for acceleration estimation

表1中5種不同方案的狀態估計結果表明：a.采用9點差分法計算得到的速度和加速度精度最差；b.當先驗噪聲協方差存在較大偏差時，Kalman濾波精度較差，再采用RTS算法進行平滑處理，其狀態估計精度提高有限；c.先采用ALS進行噪聲估計，再進行RTS平滑，其狀態估計精度有著明顯提高。圖5說明差分法的差分誤差在高頻部分特別顯著，采用ALS及RTS平滑算法不僅在高頻段對誤差有較好的抑制作用，在低頻段對信號的擬合精度也明顯優于差分法。

表1 加速度計算結果比較Table 1 Results comparison of acceleration

圖5 加速度的頻譜圖Fig.5 Spectrum of acceleration

INS具有較高的短期穩定性，但誤差隨著時間累積，而GNSS具有較高的精度，但可能存在數據空白。基于兩者的互補性，通過差分GNSS得到載體的運動加速度，將載體加速度與比力及正常重力的差值作為觀測信息來建立Kalman濾波模型，并對SINS的各種誤差進行估計，此時擾動重力矢量存在于Kalman濾波的觀測殘差中，從而可提取到測線上的擾動重力。Kalman濾波是SINS/GNSS數據融合時最為常用的數據處理方法，但其有效性有待在實際數據處理中得到驗證。

數據融合后得到的測線擾動重力中仍含有大量的高頻噪聲，大量研究和試驗表明，低通濾波技術是消除高頻噪聲的有效方法，國內外學者針對航空標量重力測量數據的低通濾波做了大量的研究[13～18]。本文采用窗函數法（濾波器1）和最佳一致逼近法（濾波器2）設計了適用于航空矢量重力測量數據處理的FIR數字低通濾波器，結合GPS靜態測量數據進行了數值試驗。結果如圖6所示，說明兩類濾波器均能以±1～2 mGal的精度確定垂直加速度，以高于±1 mGal的精度確定水平加速度。

圖6 濾波后GPS測站的加速度Fig.6 Accelerations of GPS station after lowpass filtering

3.2 數據歸算

為滿足地球重力場有關研究及相關應用的需要，航空矢量重力觀測值往往被歸算成平均飛行高度面上規則格網化的擾動重力，這一過程包括平均高度面歸算、測線網平差和格網化數據生成。平均高度面歸算本質上可理解為實測數據相對于平均高度面的向上延拓或向下延拓改正，因此平均高度面歸算可采用重力延拓方法，如譜方法[19]、泊松積分法[20]和徑向改正方法[21]等，當實際航線相對于平均高度面的變化較小時，一般只考慮徑向改正方法。測線網平差實際上是一個系統誤差的補償問題，進行系統誤差補償的一般方法為自檢校平差法，但其設計矩陣是秩虧的，這時可選定某一條或幾條精度很高的測線作為固定測線，將這些測線的系統誤差作為虛擬觀測值，按秩虧自由網平差法進行求解[22]。但在很多情況下，因觀測區域的地面檢校數據不足，觀測數據的精度無法準確評定，此時可采用兩步平差法來進行系統誤差補償[23～25]。航空重力測量數據歸算的目的是獲得平均飛行高度面上具有一定空間分辨率、按經緯度劃分的規則格網重力數據，格網化方法的好壞對格網化輸出數據的質量、精度和可信度有著直接影響。重力數據處理中常用的格網化方法主要有加權平均法、Shepard曲面擬合法、Kriging方法和最小二乘配置法（LSC），這些格網化方法有著各自的優缺點，在實際應用中有所差異。

3.3 帶限地形影響

地形起伏反映了地球重力場的不規則變化，精確確定地形對地球重力場的影響具有十分重要的意義[26]。第二類Helmert凝集法是目前應用最為廣泛的地形歸算方法，分為解析形式和譜形式，前者提供全波段的地形質量改正，后者則是將牛頓積分核函數展開成Legendre多項式的收斂級數，可用來計算帶限地形影響。

航空重力測量的各類觀測數據中包含了大量的高頻觀測噪聲，消除這些噪聲的主要手段是使用低通濾波器，濾波后得到的重力信號為帶限重力信號[27]。因此，采用解析公式計算的全頻譜地形直接影響和間接影響對航空矢量重力觀測值進行改正是不合理的，解決這一問題的方法有兩類：a.將地形改正值也作相應的低通濾波處理，使得地形效應的頻譜范圍與航空重力觀測值一致[28]，這種方法需要輸入飛行高度上密集采樣點處的地形質量影響值，采用高分辨率的數字高程模型（digital elevation model，DEM）來進行計算時工作量相當龐大；b.基于第二類Helmert凝集法的譜形式，采用帶限公式來計算地形質量對引力位和重力擾動的直接影響公式和大地水準面的間接影響[29]。這種方法的優點是可依據重力觀測值的頻譜范圍來調整積分核函數的起始階數和最高階數，從而避免了對地形質量影響作低通濾波處理，原理簡單而且計算方便。Novák等[29]基于第二類Helmert凝集法的譜形式推導了地形質量對引力位和重力擾動的帶限直接影響改正公式以及對大地水準面的帶限間接影響改正公式，將基于級數核計算的帶限直接地形影響和間接地形影響與基于解析核計算并經低通濾波處理后得到的帶限直接地形影響和間接地形影響進行了比較，驗證了算法的有效性，該帶限地形影響改正公式可對重力擾動或重力異常進行地形歸算，針對的是航空標量重力測量觀測值。航空重力矢量測量獲取的是擾動重力的三個分量，而針對水平分量帶限地形影響改正的研究尚未有公開的文獻發表。本文在Novák研究的基礎上，推導了重力矢量水平分量的帶限直接地形影響改正公式，并在中國西部山區選取了兩條沿緯圈方向4 000 m（大地高）航線高度上的航線，采用3″×3″的數字高程模型基于解析核計算了水平分量的直接地形影響，經低通濾波后與基于帶限公式計算的地形影響進行了比較。數值計算結果如圖7～圖10所示，表明二者吻合較好。

3.4 向下延拓

航空矢量重力測量的觀測值為近地空中的中高頻重力信號，但地球物理和大地測量領域的許多應用都需要地形表面或大地水準面上的重力觀測值。因此，將航線上的重力測量觀測值進行向下延拓，是航空矢量重力測量數據處理的主要內容之一。

向下延拓是一個不適定的過程，屬于病態問題，會放大重力觀測數據中的噪聲，若不采用合適的延拓方法，將會導致延拓解的嚴重失真。在現有的向下延拓方法中，逆Possion積分法的應用最為廣泛，逆Possion積分法輔以Tikhonov正則化[20]、嶺估計[30]、廣義嶺估計[30]、迭代法[31]、濾波[32，33]等正則化方法才能得到穩定的延拓解，其關鍵在于選擇合理的正則化參數。Novák和Heck[34]認為航空重力為帶限信號，提出了基于帶限Possion核函數的直接延拓法，從而避免了逆Possion積分法面臨的大型法方程陣求逆問題。LSC兼具向上延拓和向下延拓的雙重功能，在延拓過程中可融入新的重力數據且能有效估計和改善系統偏差[32]，但LSC延拓結果的好壞很大程度上取決于協方差函數模型。解析延拓也是重力信號向下延拓的一種常用方法，同樣兼具向上延拓和向下延拓的能力，其基本思想是假設重力信號的垂直梯度已知，則延拓表達式可表示為一個Taylor展開級數；解析延拓在空域中的求解十分復雜，Wu[35]給出了解析延拓在頻域中的求解公式。此外，一些數學方法如樣條函數法、小波算法等的引入，使得重力信號向下延拓的手段更加豐富。

圖7 沿31.5°N緯線4 km飛行高度處北向水平分量的直接地形影響Fig.7 Direct topographical effect for northern component at 4 km flight height along the parallel of 31.5°N

圖8 沿32.5°N緯線4 km飛行高度處北向水平分量的直接地形影響Fig.8 Direct topographical effect for northern component at 4 km flight height along the parallel of 32.5°N

圖9 沿31.5°N緯線4 km飛行高度處東向水平分量的直接地形影響Fig.9 Direct topographical effect for eastern component at 4 km flight height along the parallel of 31.5°N

以上方法一般用于航空標量重力測量數據，即重力擾動（重力異常）的向下延拓，而航空矢量重力測量不僅能獲取重力擾動（重力矢量的垂直分量），還能測出重力矢量的水平分量。理論上水平分量的向下延拓原理與垂直分量類似，但現有公開發表的文獻中針對水平分量向下延拓的研究較少，且由于航空矢量重力測量技術自身的缺陷，空中獲取的水平分量的精度遠不如垂直分量，僅為6～8 mGal，如何實現低信噪比情況下水平分量的向下延拓是一個需要解決的關鍵問題。本文引入了頻域輸入輸出系統，在不同的噪聲水平下比較分析了單輸入單輸出系統（single-input single-output system，SISOS）和雙輸入單輸出系統（double-input singleoutput system，DISOS）對水平分量的向下延拓效果（見圖11和圖12）。研究結果表明：a.當水平分量的精度較高（如1.5 mGal）時，二者均能實現水平分量的穩定向下延拓；b.當水平分量的精度較低（如6 mGal）時，單輸入單輸出法的延拓結果較差，而雙輸入單輸出法能實現水平分量的穩定向下延拓。

圖10 沿32.5°N緯線4 km飛行高度處東向水平分量的直接地形影響Fig.10 Direct topographical effect for eastern component at 4 km flight height along the parallel of 32.5°N

圖11 東向分量的精度為1.5 mGal時延拓誤差的空間分布Fig.11 Distribution of continued errors of estern component with the accuracy of 1.5 mGal

圖12 東向分量的精度為6 mGal時延拓誤差的空間分布Fig.12 Distribution of continued errors of estern component with the accuracy of 6 mGal

4 航空矢量重力應用于大地水準面的確定

利用航空矢量重力測量數據確定大地水準面可歸結為求解航空矢量重力測量邊值問題。與傳統的重力學邊值問題不同，航空矢量重力測量邊值問題的邊界面為飛行高度面，可描述為[36]

式（4）中，Δ為拉普拉斯算子；T為擾動位；H為航線高度。由式（4）可知，擾動重力的水平分量和垂直分量均可用于精化區域大地水準面。垂直分量在全波長為5～10 km的分辨率尺度上能達到1～3 mGal的精度，已能滿足區域大地水準面確定的要求。

利用垂直分量確定大地水準面的方法主要有兩類。第一類為空域積分法，它可分為兩步法和一步法[37]。兩步法包括兩個計算步驟：a.求解位理論Dirichlet邊值問題，將飛行高度面上的垂直分量向下調和延拓至大地水準面；b.求解位理論Neumann邊值問題，按照Hotine積分公式由大地水準面上的垂直分量得到大地水準面上的擾動位，再由Bruns公式得到大地水準面高。一步法是一種比兩步法更穩定、精度更高的空域解法，它將向下延拓和求解大地水準面進行了綜合。由于航空重力觀測值為帶限信號，一步法采用頻譜范圍與觀測值一致的帶限積分核函數，由廣義Hotine積分公式結合Bruns公式即可得到大地水準面高。第二類為譜方法。通過某種形式的逼近，Hotine積分、Possion積分可化為空域卷積形式。研究表明，采用快速傅里葉變換（fast Fourier transform，FFT）將Hotine積分和Possion積分轉換至頻域進行計算可使積分法的計算速度得到極大的提升[38～40]。此外，還可采用球冠諧分析或矩諧分析將重力觀測數據表達為局域重力場的諧展開式。球冠諧分析[41]和矩諧分析[42]最早應用于地磁場表達，Li等[43]首次將球冠諧分析引入區域重力場研究中，邊少峰[44]和蔣濤[30]則將矩諧分析用于局部重力場表達。對于大范圍的區域重力場逼近，球冠諧分析是一種較為理想的方法，而矩諧分析則適用于區域重力場精細結構的逼近。

類似于天文水準原理，可采用逐測線剖面積分將擾動重力的水平分量轉化為航線上的擾動位，將擾動位向下延拓至大地水準面后，再采用Bruns公式計算大地水準面起伏[45]。值得注意的是，垂直分量確定的是絕對大地水準面，而水平分量確定的是相對大地水準面。

在聯合重力矢量的三個分量精化大地水準面方面，隨機解法、頻域輸入輸出法、多分辨率配置法、譜組合法等方法可將航線的重力矢量轉化為航線上的擾動位，將擾動位向下延拓至大地水準面，再采用Bruns公式計算大地水準面起伏。目前航空矢量重力測量水平分量的精度水平遠不如垂直分量，需要分析不同信噪比情況下的水平分量對三分量聯合求解大地水準面的影響，進一步研究聯合航空重力矢量的三個分量確定大地水準面的融合方法。

5 結語

基于SINS/GNSS的航空矢量重力測量是未來高效測定中高頻地球重力場信息的先進技術，聯合地面重力測量、海洋重力測量和衛星重力測量技術，可以獲取全球高分辨率高精度地球重力場信息，從而滿足國家經濟建設、軍事和國防建設、地球物理、海洋學、地球動力學、資源勘探等相關地球科學領域的需求。本文在國內外現有研究成果的基礎上，討論了航空矢量重力測量的數據預處理、數據歸算、帶限地形影響、向下延拓及大地水準面確定等問題。

航空矢量重力測量技術的實用化還需要解決一系列關鍵問題，如高精度SINS的自主研制、SINS/ GNSS的數據融合方法及滿足工程化應用的數據處理軟件、SINS和GNSS的誤差校正、地殼橫向密度擾動對向下延拓和精化大地水準面的影響、航空矢量重力及其他類型重力場信息的融合處理等。

[1]Wei M，Schwarz K P.Flight test results from a strapdown airborne gravity system[J].Journal of Geodesy，1998，72：323-332.

[2]Verdun J，Klingelé E，Bayer R，et al.The alpine Swiss-French airborne gravity survey[J].Geophysical Journal International，2003，152：8-19.

[3]Forsberg R，Olesen A，Vest A，et al.Gravity field improvements in the north atlantic region[C]//Proc.GOCE Workshop，ESA-ESRIN，2004.

[4]夏哲仁，石 磐，孫中苗，等.航空重力測量系統CHAGS[J].測繪學報，2004，33（3）：216-220.

[5]Jekeli C，Garcia R.GNSS phase accelerations for movingbase vector gravimetry[J].Journal of Geodesy，1997，71：630-639.

[6]Kwon J.Airborne vector gravimetry using GNSS/INS[D].USA：The Ohio State University，2000.

[7]Kwon J，Jekeli C.A new approach for airborne vector gravimetry using GNSS/INS[J].Journal of Geodesy，2001，74（10）：690-700.

[8]Senobari M S.New results in airborne vector gravimetry using strapdown INS/GNSS[J].J Geoid，2010，84：277-291.

[9]Smith D.The GRAV-D project：Gravity for the redefinition of the American Vertical Datum[OL].NOAA website：http：//www. ngs.noaa.gov/GRAV-D/pubs/GRAV-D_v2007_12_19.pdf，2007.

[10]Lee J G，Park C G，Park H W.Multiposition alignment of strapdown inertial navigation system[J].Aerospace and Electronic Systems，IEEE Transactions on，1993，29（4）：1323-1328.

[11]Kwon J H，Jekeli C.A new approach for airborne vector gravimetry using GNSS/INS[J].Journal of Geodesy，2001，74（10）：690-700.

[12]孫中苗.航空重力測量理論，方法及應用研究[D].鄭州：解放軍信息工程大學，2004.

[13]Singer R A.Estimating optimal tracking filter performance for manned maneuvering targets[J].IEEE Transactions on Aerospace and Electronic Systems，1970（4）：473-483.

[14]Hammada Y.A comparison of filtering techniques for airborne gravimetry[D].Calgary：Department of Geomatics Engineering at University of Calgary，1996.

[15]孫中苗.航空重力測量中FIR低通濾波器的設計與應用[J].解放軍測繪學院學報，1996，13（4）：247-250.

[16]柳林濤，許厚澤.航空重力測量數據的小波濾波處理[J].地球物理學報，2004，47（3）：490-494.

[17]孫中苗，翟振和.航空重力測量數據的小波閾值濾波[J].武漢大學學報（信息科學版），2009，34（10）：1222-1225.

[18]蔡劭琨，吳美平，張開東.航空重力測量中FIR低通濾波器的比較[J].物探與化探，2010，34（1）：74-78.

[19]蔣 濤，黨亞明，章傳銀，等.基于矩諧分析的航空重力向下延拓[J].測繪學報，2013，42（4）：475-480.

[20]王興濤，石 磐，朱非洲.航空重力測量數據向下延拓的正則化算法及其譜分解[J].測繪學報，2004，33（1）：33-38.

[21]鐘 波，劉華亮，羅志才，等.衛星重力梯度數據的歸算與格網化處理[J].大地測量與地球動力學，2011，31（3）：79-84.

[22]Hwang C，Hsiao Y S，Shih H C.Data reduction in scalar airborne gravimetry：Theory，software and case study in Taiwan [J].Computers&Geosciences，2006，32（10）：1573-1584.

[23]黃謨濤，翟國君，歐陽永忠，等.海洋重力測量誤差補償兩步處理法[J].武漢大學學報（信息科學版），2002，27（3）：251-255.

[24]周波陽，羅志才，林 旭，等.航空重力測量測線網平差中的粗差處理[J].大地測量與地球動力學，2012，32（2）：110-114.

[25]周波陽，羅志才，林 旭.航空重力測量測線網平差中的粗差處理方法[C]//中國地球物理學會第二十七屆年會論文集.合肥：中國科技大學出版社，2011.

[26]羅志才，寧津生，晁定波.衛星重力梯度向下延拓的譜方法[J].測繪學報，1997，26（2）：168-175.

[27]Novák P，Kern M，Schwarz K P.Numerical studies on the harmonic downward continuation of band-limited airborne gravity [J].Studia Geophysica et Geodaetica，2001，45（4）：327-345.

[28]Bayoud F A，Sideris M G.Two different methodologies for geoid determination from ground and airborne gravity data[J]. Geophysical Journal International，2003，155（3）：914-922.

[29]Novák P，Kern M，Schwarz K P，et al.Evaluation of band-limited topographical effects in airborne gravimetry[J].Journal of Geodesy，2003，76（11-12）：597-604.

[30]蔣 濤.利用航空重力測量數據確定區域大地水準面[D].武漢：武漢大學，2011.

[31]Serpas J.Local and regional geoid determination from airborne vector gravimetry[D].USA：The Ohio State University，2003.

[32]Hwang C，Hsiao Y S，Shih H C，et al.Geodetic and geophysical results from a Taiwan airborne gravity survey：Data reduction and accuracy assessment[J].Journal of Geophysical Research：Solid Earth（1978-2012），2007，112（B4）.

[33]周波陽，羅志才，許 闖，等.航空重力測量數據向下延拓最小二乘配置法和逆Possion積分法的比較[C]//中國地球物理.合肥：中國科技大學出版社，2012.

[34]Novák P，Heck B.Downward continuation and geoid determination based on band-limited airborne gravity data[J].Journal of Geodesy，2002，76（5）：269-278.

[35]Wu L.Spectral methods for post processing of airborne vector gravity data[D].Canada：University of Calgary，1996.

[36]Schwarz K P，Li Z.An introduction to airborne gravimetry and its boundary value problems[M]//Geodetic boundary value problems in view of the one centimeter geoid.Germany：SpringerBerlin Heidelberg，1997：312-358.

[37]Novák P.Optional model for geoid determination from airborne gravity[J].Studia geophysica et geodaetica，2003，47：1-36.

[38]Schwarz K P，Sideris M G，Forsberg R.The use of FFT techniques in physical geodesy[J].Geophysical Journal International，1990，100（3）：485-514.

[39]Sideris M G，She B B.A new，high-resolution geoid for Canada and part of the US by the 1D-FFT method[J].Bulletin Géodésique，1995，69（2）：92-108.

[40]周波陽，羅志才，林 旭，等.航空重力數據向下延拓FFT快速算法比較[J].大地測量與地球動力學，2013，33（1）：64-68.

[41]Haines G V.Spherical cap harmonic analysis[J].Journal of Geophysical Research：Solid Earth（1978-2012），1985，90（B3）：2583-2591.

[42]Alldredge L R.Rectangular harmonic analysis applied to the geomagnetic field[J].Journal of Geophysical Research，1981，86（B4）：3021-3026.

[43]Li J C，Chao D B，Ning J S.Spherical cap harmonic expansion for local gravity field representation[J].Manuscripta Geodaetica，1995，20：265-277.

[44]邊少峰.大地測量邊值問題數值解法與地球重力場逼近[D].武漢：武漢測繪科技大學，1992.

[45]Jekeli C，Kwon J H.Geoid profile determination by direct integration of GNSS inertial navigation system vector gravimetry[J].Journal of Geophysical Research：Solid Earth（1978-2012），2002，107（B10）：ETG 3-1-ETG 3-10.

Research on the data processing methods of airborne vector gravimetry using SINS/GNSS

Ning Jinsheng1，2
（1.School of Geodesy and Geomatics，Wuhan University，Wuhan 430079，China；2.Key Laboratory of Geospace Environment and Geodesy，Ministry of Education，Wuhan University，Wuhan 430079，China）

Airborne vector gravimetry is an advanced and efficient technology to determine high frequency information of earth’s gravity field in the future，and is a complement to ground gravimetry，marine gravimetry and satellite gravimetry.The principle of airborne vector gravimetry using SINS（strapdown inertial navigation system）/GNSS（global navigation satellite systems）is introduced in the paper，and then the data preprocessing，data reduction，band-limited topographical effect and downward continuation are discussed，as well as the geoid determination from airborne vector gravimetry.It provides the support for the development of airborne vector gravimetry in our country.

airborne vector gravimetry；SINS/GNSS data fusion；earth gravity field；geoid

P227

A

1009-1742（2014）03-0004-10

2013-12-03

國家自然科學基金項目（41174062）；地球空間環境與大地測量教育部重點實驗室開放基金項目（12-02-05，12-02-09）

寧津生，1932年出生，安徽桐城縣人，中國工程院院士，教授，主要從事物理大地測量的教學和科研工作；E-mail：jsning@sgg.whu.edu.cn