航空重力測(cè)量中差分器的設(shè)計(jì)與比較

孫中苗， 孔　超， 溫　波

1.西安測(cè)繪研究所，陜西 西安，710054；2.地理信息工程國(guó)家重點(diǎn)實(shí)驗(yàn)室，陜西 西安，710054；3.信息工程大學(xué)地理空間信息學(xué)院，河南 鄭州，450052

?

航空重力測(cè)量中差分器的設(shè)計(jì)與比較

孫中苗1,2， 孔超3， 溫波3

1.西安測(cè)繪研究所，陜西 西安，710054；2.地理信息工程國(guó)家重點(diǎn)實(shí)驗(yàn)室，陜西 西安，710054；3.信息工程大學(xué)地理空間信息學(xué)院，河南 鄭州，450052

航空重力測(cè)量中垂直加速度的計(jì)算需要差分運(yùn)算。本文給出了5類常用差分器的系數(shù)及傳遞函數(shù)，以兩點(diǎn)中心差分、2階牛頓差分、3階牛頓差分、2階Lanczos、移動(dòng)平均差分和49階等波紋差分器為例，對(duì)差分器的性能做了分析比較。與理想差分器比較，結(jié)果表明，6種差分器的幅頻響應(yīng)在0～0.01Hz的低頻段均與理想差分器相一致，而在0.01～0.10Hz 頻段只有3階牛頓差分器和等波紋差分器吻合較好；與靜態(tài)GPS測(cè)量數(shù)據(jù)比較，結(jié)果表明，6種差分器獲得的垂直加速度經(jīng)20s低通濾波后標(biāo)準(zhǔn)差幾近相同，鑒于兩點(diǎn)中心差分器的邊界效應(yīng)最小，可在實(shí)踐中加以采用。

航空重力測(cè)量，垂直加速度，差分器，傳遞函數(shù)

1　引　言

航空重力測(cè)量在獲取地面重力測(cè)量難以到達(dá)地區(qū)的地球重力場(chǎng)信息方面已得到廣泛應(yīng)用[1-3]。目前實(shí)用的航空標(biāo)量重力測(cè)量系統(tǒng)包括平臺(tái)式和捷聯(lián)式兩類。平臺(tái)式將航空重力儀固定在阻尼平臺(tái)上進(jìn)行測(cè)量，如LCR航空重力儀[4]；捷聯(lián)式是捷聯(lián)式慣性導(dǎo)航和GPS組合的重力測(cè)量系統(tǒng)[5]。無論是平臺(tái)式系統(tǒng)還是捷聯(lián)式系統(tǒng)，獲取高質(zhì)量成果的前提是具備良好的飛行作業(yè)條件和完善的數(shù)據(jù)處理算法，其中，提高因飛機(jī)運(yùn)動(dòng)產(chǎn)生的垂直加速度的確定精度一直是重中之重。

垂直加速度目前一般利用GPS和數(shù)字濾波技術(shù)相結(jié)合進(jìn)行確定。慣用的做法是，先利用GPS觀測(cè)數(shù)據(jù)計(jì)算出運(yùn)動(dòng)載體的位置或速度，然后差分(二次或一次)獲取加速度，最后作低通濾波以消除高頻噪聲的影響。大量文獻(xiàn)研究表明[6-10]，在良好的環(huán)境條件下，對(duì)于90s～300s的濾波周期，利用上述方法可以1～2mGal的精度導(dǎo)出運(yùn)動(dòng)平臺(tái)的加速度。但由垂直加速度的頻譜特性分析表明[10]，垂直加速度不僅占據(jù)相當(dāng)寬的頻帶，而且頻譜特性隨測(cè)量環(huán)境(如不同的運(yùn)載平臺(tái)和大氣湍流條件)的變化也不盡相同。因此，若欲同時(shí)提高航空重力測(cè)量的精度和分辨率，在與航空重力相關(guān)的頻帶上，垂直加速度的精密確定仍是主要障礙。

大體而言，垂直加速度的確定精度取決于三個(gè)方面：一是GPS定位或測(cè)速的精度，主要體現(xiàn)在GPS各種誤差源的影響；二是差分方法；三是低通濾波器的性能，它反映了垂直加速度精度所對(duì)應(yīng)的分辨率。本文主要對(duì)其中差分方法的影響進(jìn)行分析，重點(diǎn)比較分析了六種常用差分器的性能，并以此給出了建議的差分方法。

2　垂直加速度計(jì)算方法

假設(shè)利用GPS觀測(cè)數(shù)據(jù)得到的等距離散位置序列為xn(n=1,…,NP)，則位置序列的一次差分即速度序列v(n)可表示為如下卷積形式：

(1)

式中，NP是位置序列的長(zhǎng)度；h(k)表示長(zhǎng)度為(2M+1)的有限沖擊響應(yīng)(FIR)差分器的系數(shù)。

同樣，由式(1)的速度序列可得加速度序列a(n)為：

(2)

顯然，由式(2)確定的加速度的精度除與位置或速度精度有關(guān)外，還取決于所用差分器的性能。

3　常用差分器的系數(shù)及傳遞函數(shù)

3.1理想差分器

一階理想差分器的傳遞函數(shù)為：

Hd(ejω)=jω

(3)

其幅頻響應(yīng)為過原點(diǎn)的一條斜線。

航空重力測(cè)量數(shù)據(jù)中含有許多高頻噪聲，故希望使用低通差分濾波器，其傳遞函數(shù)為一截?cái)嗪瘮?shù)：

(4)

式中απ是低通差分器的截止頻率。

3.2單純M次差分

這是最簡(jiǎn)單的差分算法，記：

Δk=x(n+k)-x(n-k)

(5)

則這種差分器的輸入輸出關(guān)系為：

(6)

其傳遞函數(shù)為：

(7)

若M=1，即為兩點(diǎn)中心差分：

(8)

(9)

3.3牛頓-柯斯特差分

這種差分器源自牛頓-柯斯特積分公式，其傳遞函數(shù)是：

(10)

當(dāng)M=2時(shí)，

(11)

(12)

當(dāng)M=3時(shí)，

(13)

(14)

3.4多項(xiàng)式擬合差分器(Lanczos差分器)

(15)

當(dāng)M=2時(shí)，

(16)

(17)

3.5平滑化差分

這是一種數(shù)據(jù)平滑和差分相結(jié)合的算法，目的在于在低頻段H(ejω)能更好地近似于jω，而在高頻段能有更大的衰減。

(18)

其頻率特性為：

(19)

可知H(ejω)是兩項(xiàng)的乘積，第一項(xiàng)表現(xiàn)了加權(quán)平滑特性，第二項(xiàng)體現(xiàn)了差分特性。

(20)

(21)

當(dāng)N=1，L=1時(shí)，

(22)

(23)

3.6等波紋差分器

在FIR濾波器的等波紋設(shè)計(jì)法中[11]，利用契比雪夫“交錯(cuò)點(diǎn)組定理”和數(shù)值分析中的Remez算法，可惟一地確定濾波器系數(shù)組{a(n)}，由此可得到濾波器的傳遞函數(shù)H(ejω)。設(shè)計(jì)低通濾波器時(shí)逼近的是幅頻響應(yīng)在通帶和阻帶的函數(shù)值，差分器的設(shè)計(jì)過程與低通濾波器完全一樣，只是其逼近的是頻率響應(yīng)在各帶的斜率。

以下討論中本文利用等波紋設(shè)計(jì)法設(shè)計(jì)了49階的等波紋差分器。

4　差分器的性能比較與分析

4.1與理想差分器的直觀比較

以常用的兩點(diǎn)中心差分、2階牛頓差分、3階牛頓差分、2階Lanczos、移動(dòng)平均差分和49階等波紋差分器為例，為直觀地評(píng)價(jià)上述各類差分器的性能，圖1(a)示出了各類差分器的幅頻響應(yīng)，其與理想差分器幅頻響應(yīng)的差值示于圖1(b)。顯而易見，在0～0.01Hz的低頻段，各個(gè)差分器的幅頻響應(yīng)均與理想差分器趨于一致，但在0.01～0.10Hz僅有3次牛頓-柯斯特差分器、等波紋差分器與其相吻合。因此，若目標(biāo)是差分獲取低于0.01Hz頻段內(nèi)的信號(hào)，可采用上述差分器中的任何一種；若需恢復(fù)0.01Hz以上頻譜，應(yīng)當(dāng)采用3次以上牛頓-柯斯特差分器或高階等波紋差分器。

圖1　各類差分器的幅頻響應(yīng)及與理想差分器的比較

4.2實(shí)測(cè)數(shù)據(jù)比較分析

按照航空重力測(cè)量作業(yè)要求，飛機(jī)起飛前需要進(jìn)行10～15m的靜態(tài)測(cè)量。靜態(tài)測(cè)量時(shí)，可以認(rèn)為飛機(jī)在作速度/加速度為零的運(yùn)動(dòng)，因此將所求速度/垂直加速度與真值'0'進(jìn)行比較，可以估算速度/垂直加速度的確定精度。實(shí)驗(yàn)數(shù)據(jù)源于2002年3月進(jìn)行的大同航空重力測(cè)量的靜態(tài)GPS觀測(cè)數(shù)據(jù)[10]。

為了對(duì)直觀比較結(jié)果作進(jìn)一步說明，圖2(a)、(b)分別示出了利用飛機(jī)起飛前靜態(tài)測(cè)量數(shù)據(jù)計(jì)算的高程和速度的幅度譜。表1統(tǒng)計(jì)了由高程數(shù)據(jù)導(dǎo)出的速度的誤差，顯然，差分器E、F的精度最高，差分器B、C精度次之，而差分器D、G精度最差。這是由于高程的能量主要集中在0.01Hz以內(nèi)，差分器E、F、B、C代之以在>0.01Hz頻段內(nèi)產(chǎn)生差分誤差，其所具有的低通特性部分地將高頻噪聲進(jìn)行了阻尼；而差分器D、G盡管在更寬的頻帶上接近理想差分器，卻在高頻段的差分運(yùn)算中放大了高頻噪聲。

表2統(tǒng)計(jì)了由垂直速度導(dǎo)出的垂直加速度的標(biāo)準(zhǔn)差，未濾波時(shí)可得出上述類似結(jié)論，20s低通濾波后，6種差分器的精度相當(dāng)。

航空重力測(cè)量中，經(jīng)差分運(yùn)算獲得加速度后，通常需要采用90～300s濾波尺度的低通濾波器進(jìn)一步消除高頻噪聲的影響。這一過程實(shí)質(zhì)上是差分器和低通濾波器的卷積運(yùn)算。從式(8)可以看出，兩點(diǎn)中心差分器的邊界效應(yīng)只有1，且其使用簡(jiǎn)便，建議在實(shí)踐中加以采用。

圖2　高程和速度的幅頻響應(yīng)(靜態(tài)測(cè)量)

表1靜態(tài)高程數(shù)據(jù)計(jì)算的速度的誤差(mm/s)

比較項(xiàng)目╲差分器BCDEFG差值均值-0.06-0.06-0.06-0.06-0.06-0.07標(biāo)準(zhǔn)偏差2.222.823.141.301.363.99

表2靜態(tài)垂直速度導(dǎo)出的加速度的誤差(mGal)

濾波尺度╲差分器BCDEFG0s175.6233.7261.868.677.7319.75s71.8783.1385.7446.7650.6486.8210s20.421.621.717.818.221.720s4.84.94.94.74.74.9

5　結(jié)束語

航空重力測(cè)量中垂直加速度的精確計(jì)算至關(guān)重要。垂直加速度通常由高程或速度通過差分運(yùn)算獲得。本文比較分析了六種常用差分器的性能，初步結(jié)論如下。

(1)與理想差分器幅頻響應(yīng)進(jìn)行直觀比較，可以看出，兩點(diǎn)中心差分、2階牛頓差分等6種差分器的幅頻響應(yīng)在0～0.01Hz的低頻段均與理想差分器趨于一致，但在0.01～0.10Hz僅有3階牛頓差分器和等波紋差分器吻合較好；如欲獲取0.01Hz以上高頻重力信號(hào)，應(yīng)當(dāng)采用3次以上牛頓差分器或高階等波紋差分器。

(2)從靜態(tài)GPS測(cè)量數(shù)據(jù)比較分析看，6種差分器獲得的垂直加速度經(jīng)20s低通濾波后，其標(biāo)準(zhǔn)差相差很小。

(3)航空重力測(cè)量中垂直加速度的計(jì)算實(shí)質(zhì)是差分器和低通濾波器的卷積運(yùn)算，鑒于兩點(diǎn)中心差分器的邊界效應(yīng)最小，可在實(shí)踐中加以采用。

[1]Brozena J, Childers V A. The NRL Airborne Geophysics Program [A]. In International Association of Geodesy Symposia, Geodesy Beyond 2000-The Challenges of the First Decade. Heidelburg: Springer Verlag,2000.

[2]Forsberg R, Olesen A V. Airborne Gravity and Geoid Surveys in the Arctic and Baltic Seas [A]. In Proceedings of IAG International Symposium on Gravity, Geoid and Geodynamics 2001. Banff, Canada,2001.

[3]Klingele E E, Cocard M, Halliday M E, Kahle H G. The Airborne Gravimetric Survey of Switzerland [M]. Contribution to the Geology of Switzerland Series, Swiss Geophysical Commission, Report No.31,1995.

[4]Valliant H D. The LaCoste and Romberg Air/Sea Gravity Meter: an Overview [M]. CRC Handbook of Geophysical Exploration at Sea, 2nd edition,1991.[5]Bruton A M. Improving the Accuracy and Resolution of SINS/GPS Airborne Gravimetry [D]. Calgary: Department of Geomatics Engineering at University of Calgary,2001.[6]Czompo J. Airborne Scalar Gravimetry Systems in the Spectral Domain [D]. Calgary: Department of Geomatics Engineering at the University of Calgary,2004.[7]Bruton A M. Improving the Accuracy and Resolution of SINS/GPS Airborne Gravimetry [D]. Calgary: Department of Geomatics Engineering at University of Calgary,2001.[8]李曉斌. 航空重力測(cè)量中GPS定位、速度及加速度解算技術(shù)與應(yīng)用[D]. 北京：中國(guó)地質(zhì)大學(xué),2010.

[9]歐陽永忠. ?？罩亓y(cè)量數(shù)據(jù)處理關(guān)健技術(shù)研究[D]. 武漢：武漢大學(xué),2013.

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

[11]胡廣書. 數(shù)字信號(hào)處理：理論,算法與實(shí)現(xiàn)(第二版)[M]. 北京:清華大學(xué)出版社,2003.

Design and Analysis of the Differentiator in Airborne Gravimetry

Sun Zhongmiao1,2，Kong Chao3, Wen Bo3

1 Xi’an Research Institute of Surveying and Mapping, Xi’an 710054, China 2 State Key Laboratory of Geo-information Engineering, Xi’an 710054, China 3.Institute of Geospatial Information, Information Engineering University, Zhengzhou 450052, China

Difference is needed to calculate the vertical acceleration calculation in airborne gravimetry. The coefficients and transfer functions of five types of common differentiators are provided in this paper. Six differentiators are set as examples to compare their performance, including differentiator of two centers, Newton differentiator of second order, Newton differentiator of third order, Lanczos differentiator of second order, moving average differentiator and equal-ripple differentiator of forty-ninth order. Compared with the ideal differentiator, it is indicated that the amplitude frequency responses of all the six differentiators are consistent with those of ideal differentiator in the frequency bands of 0～0.01Hz, while in the frequency bands of 0.01～0.10Hz only the Newton differentiator of third order and the equal-ripple differentiator will coincide well. The static GPS data indicates that the standard deviations of the vertical accelerations achieved by above six differentiators respectively are almost the same after using a lowpass filter with cutoff frequency of 20 seconds. Since the edge effect of the differentiator of two centers is the least, it can be applied to practice.

airborne gravimetry; vertical acceleration; differentiator; transfer function

2015-12-24。

孫中苗(1968—)，男，研究員，主要從事航空重力測(cè)量理論和方法研究。

P223

A