航空重力測量中牛頓中心差分器對確定載體加速度的影響

周波陽 羅志才 鐘 波 鄭 凱 魏艷平

1 廣東工業大學測繪工程系，廣州市大學城外環西路100號，510006

2 武漢大學測繪學院，武漢市珞喻路129號，430079

3 武漢大學地球空間環境與大地測量教育部重點實驗室，武漢市珞喻路129號，430079

4 中國石油集團工程設計有限責任公司華北分公司，任丘市建設路，062552

對機載應用來說，確定載體加速度的3種基本方法是位置求導法、相位差分法、多普勒頻移法［1－4］，這3種方法都需要通過差分運算才能得出載體的加速度。牛頓中心差分器是運用最為廣泛的一類差分器，其關鍵在于點數的選取。文獻［5］分析了利用一階牛頓中心差分器進行GPS定速的主要誤差來源，認為一階牛頓中心差分器的最佳點數應使截斷誤差與觀測誤差之和最小。這一結論是在全頻譜范圍內進行考慮的，適合于航空攝影測量、GPS／INS組合導航等領域。航空重力測量最終獲取的是帶限重力信號，需要將載體加速度從比力觀測值中分離出來。由于重力信號表現為長波特性，在數據處理中需對原始觀測值進行低通濾波處理，因此對航空重力測量而言，傳感器原始觀測值的高頻部分都可認為是噪聲，載體加速度只需在低頻端保持足夠的精度即可。考慮到這一情況，在航空重力測量的數據處理中，在全頻譜范圍內考慮差分器的精度可能并不合適。本文將從頻譜的角度分析牛頓中心差分器的點數對航空重力測量中載體加速度的影響，并采用實測的機載GPS數據進行驗證。

1 理想差分器與理想低通差分器

利用載體的位置信息通過差分方式可以計算出載體的速度和加速度。為簡化問題，本文僅以垂直方向為例。由微分的定義可知：

對式（2）取極限，得：

因此，理想差分器的頻率響應為：

式中，ω＝2πf為圓頻率。由式（4）可知，理想差分器的頻率特性在0～π內是線性增長的。和信號相比，噪聲一般在高頻端，而差分過程會放大噪聲，特別是高頻噪聲。這一問題的解決辦法是采用低通差分器，理想低通差分濾波器的頻率響應為：

式中，απ為理想低通差分器的截止頻率，頻率歸一化后則稱α為截止頻率。

2 牛頓中心差分器

理想低通差分器在物理上無法實現，人們通常采用頻率響應函數H（ejω）來逼近H′d（ejω）。此時，H（ejω）一般具有以下形式：

相應的轉移函數為：

若差分器的輸入為x（n），輸出為y（n），則有：

表1 一階牛頓中心差分器的系數Tab.1 Coefficients of one－order Newton central differentiator

依據式（4）、（10），可給出理想差分器與牛頓中心差分器的幅頻響應，如圖1所示。從圖中可以看出，牛頓中心差分器是對理想差分器的逼近，點數越多，對理想差分器的逼近效果越好，但同時差分器的低通效果在理論上也越差。文獻［7］指出，低通差分器對理想差分器的逼近也可采用截止頻率α來描述：απ為曲線｜H（ejω）｜與直線的交點的橫坐標，如圖1所示。表2給出了不同點數牛頓中心差分器的截止頻率αi（i＝3，5，7，9，11），可知α3＜α5＜α7＜α9＜α11。當信號的有效頻段位于低頻段，且信號的截止頻率α?α3時，從理論上可以得出兩條結論：1）為避免高頻噪聲的影響，應對差分結果進行低通濾波；2）由于不同點數的牛頓中心差分器的幅頻響應在低頻端完全吻合，低通濾波后不同點數的牛頓中心差分器的輸出結果之間的差異應比較小。

圖1 理想差分器及牛頓中心差分器的幅頻特性Fig.1 Magnitude－frequency characteristic of ideal differentiator and Newton central differentiators

表2 一階牛頓中心差分器的截止頻率Tab.2 Cutoff frequencies of one－order Newton central differentiators

3 數值實驗與分析

3.1 靜態實驗

選取美國某地區2個IGS跟蹤站的GPS觀測數據，測站編號分別為A和B，基線長度為30 km，數據采樣間隔為1s，觀測時間為2h。設測站A的坐標固定，將B的坐標觀測值與A的固定坐標相減，可得到時間序列→ZAB，對→ZAB作二次牛頓中心差分后可得到測站B相對于測站A的加速度。由于差分過程會放大噪聲，必須對差分結果作低通濾波處理。本文選用文獻［9］給出的迭代高斯低通濾波器對加速度進行濾波，若載體的飛行速度為155 m／s，所需重力數據的半波長分辨率為10km，則高斯濾波器的長度為65s，對應的截止頻率為α＝0.015Hz。在短時間內可認為測站B相對加速度的真值為零，因此低通濾波后所得加速度的標準差即可反映數值處理的絕對精度，表3給出了相關的統計信息。從表中可以看出，3點牛頓中心差分器輸出結果的精度最優，但隨著差分器點數的增加，輸出加速度的精度也隨之降低，同時也可認為采用牛頓中心差分器結合低通濾波器的處理方法是比較可靠的，其絕對精度優于1.8mGal。

表3 基準站B 的加速度絕對誤差統計信息Tab.3 Statistics of absolute acceleration error of station B

3.2 動態實驗

數據來源于美國某地區的某次航空重力測量結果。測區位于路易斯安那州，航線高度約為10 668m，飛機的地面速度約為155m／s，數據采樣率為1 Hz，觀測時間為6.6h。采用武漢大學研制的精密單點定位軟件TriP V2.0［10］解算飛機的位置信息，飛行軌跡如圖2所示。航空重力測量作業的理想情況是飛機能保持直線平穩水平飛行，因此我們選取圖2中的淺色部分作為測試數據，該段數據包括3 200個歷元。

圖2 航跡圖Fig.2 Fight trajectory

對飛機z方向的位置信息作二次牛頓中心差分可得到加速度，直接差分所得加速度的頻譜包含所有頻段內的信息。由于觀測噪聲不可避免，而差分過程會對噪聲進行放大，如果采用理想差分器進行數據處理，理論上在全頻段內噪聲的頻率越高，對噪聲放大得越嚴重。但牛頓中心差分器只是對理想差分器進行逼近，以3點牛頓中心差分為例，由表2 可知其截止頻率為α3＝0.603Hz。由圖3可知，在［0，α3］內，噪聲頻率越高放大效果越明顯，甚至遠遠超過了信號本身的大小；而在［α3，1］范圍內差分器對噪聲的放大程度有所減緩，說明差分器也具有一部分低通濾波的效果。上述結果也與圖1中牛頓中心差分器的幅頻特性一致。

已有航空重力測量觀測值的頻譜分析表明，重力信號表現為長波特性，在高頻端的能量幾乎可以忽略不計，但圖3表明，直接差分后得到的加速度的能量主要集中在高頻端，因此需對直接差分得到的加速度作低通濾波處理。低通濾波器為靜態實驗中提及的迭代高斯低通濾波器，其對應的截止頻率為α0＝0.015 Hz，易知α0?α3。圖3給出了對位置信息進行兩次3點牛頓中心差分后所得加速度在低通濾波前后的PSD。從圖中可以看出，低通濾波后高頻噪聲已得到有效消除。

圖3 低通濾波前后加速度的PSDFig.3 PSD of acceleration before and after low－pass filtering

由于沒有真值作為比較，這里僅分別將5、7、9、11點牛頓中心差分器與3點牛頓中心差分器的輸出結果進行比較。表4給出了低通濾波前后不同點數的牛頓中心差分器輸出結果的互差統計信息，其中（3，5）表示5點牛頓中心差分器與3點牛頓中心差分器輸出結果的較差，依此類推。由表4可知，低通濾波前5、7、9、11點牛頓中心差分器與3點牛頓中心差分器的輸出結果較差的標準偏差在102mGal級，而濾波后較差的標準偏差已經降低至2.1mGal以內。一方面，不同點數的牛頓中心差分器輸出結果之間差異較小，這是因為高斯濾波器的通帶［0，α0］位于所有牛頓中心差分的通帶［0，αi］（i＝3，5，7，9，11）內，且圖1表明，不同點數的牛頓中心差分器的幅頻響應在［0，α0］內完全一致，有差異的部分則位于高斯低通濾波器的阻帶范圍內；另一方面，因不同點數的牛頓中心差分器的截止頻率不一樣，在高頻段對噪聲的放大情況也不一致，同時低通濾波器不能濾除所有的高頻噪聲，這是不同點數的牛頓中心差分器的輸出結果之間存在較差的主要原因。此外，因3點牛頓中心差分器的截止頻率最小，對高頻噪聲的放大效果要弱于其他點數的牛頓中心差分器，同時牛頓3點中心差分器計算方便，使用簡單，邊界效應小。因此，在實際應用中建議使用牛頓3點中心差分器。

4 結 語

航空重力測量的有效頻段位于低頻端，采用牛頓中心差分器得到載體的加速度后還需對加速度作低通濾波處理。本文分析了不同點數的牛頓中心差分器的幅頻響應，從理論和數值實驗兩方面驗證了當低通濾波器的通帶位于所有牛頓中心差分器的通帶內時，不同點數的牛頓中心差分器所得結果經低通濾波后差異較小，建議航空重力測量數據處理時采用3點牛頓中心差分器。衛星重力測量中高低衛－衛跟蹤數據差分衛星加速度的過程也可采用類似的方法進行分析。

表4 低通濾波前后不同點數的牛頓中心差分器的輸出結果比較Tab.4 Result comparison of different point Newton central differentiator before and after low－pass filtering

［1］Ryan S，Lachapelle G，Gannon ME.DGPS Kinematic Carrier Phase Signal Simulation Analysis in the Velocity Domain［C］.ION，Kansas，1997

［2］Bruton A M，Gleeie G L，Schwarz K P.Differentiation for High－Precision GPS Velocity and Acceleration Determination［J］.GPS Solution，1999，2（4）：7－21

［3］肖云.利用GPS確定航空重力測量載體運動狀態的理論與方法［D］.武漢：武漢測繪科技大學，2000（Xiao Yun.Research on the Theories and Methods of Determination of Moving－Base Status for Airborne Gravimetry by Using GPS［D］.Wuhan：Wuhan Technical University of Surveying and Mapping，2000）

［4］肖云，夏哲仁.利用相位率和多普勒觀測值精確確定載體的加速度［J］.武漢大學學報：信息科學版，2003，28（5）：581－584（Xiao Yun，Xia Zheren.Comparison between Phase－Rate and Doppler to Determine Velocity［J］.Geomatics and Information Science of Wuhan University，2003，28（5）：581－584）

［5］王潛心，徐天河，龔佑興.利用中心差分法進行GPS定速時最佳點數的選取［J］.武漢大學學報：信息科學版，2012，37（3）：265－268（Wang Qianxin，Xu Tianhe，Gong Youxing.Optimal Points of Using Central Differences Method for GPS Velocity Determination［J］.Geomatics and Information Science of Wuhan University，2012，37（3）：265－268）

［6］歐陽永忠.海空重力測量數據處理關鍵技術研究［D］.武漢：武漢大學，2013（Ouyang Yongzhong.Study on Key Technologies of Data Processing for Air－Sea Gravity Survey［D］.Wuhan：Wuhan University，2013）

［7］胡廣書.數字信號處理理論、算法與實現［M］.北京：清華大學出版社，2012（Hu Guangshu.Theory，Algorithm and Implementation of Digital Signal Processing［M］.Beijing：Tsinghua University Press，2012）

［8］Cynar S J.Using Guassian Elimination for Computation of the Central Difference Coefficients［J］.Acm Signum Newsletter，1987，22（2）：12－19

［9］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：1 573－1 584

［10］張小紅，劉經南，Forsberg R.基于精密單點定位技術的航空測量應用實踐［J］.武漢大學學報：信息科學版，2006，31（1）：19－22（Zhang Xiaohong，Liu Jingnan，Rene Forsberg.Application of Precise Point Positioning in Airborne Survey［J］.Geomatics and Information Science of Wuhan University，2006，31（1）：19－22）