楊 曄，李 達(dá)，高 巍

（天津航海儀器研究所，天津 300131）

旋轉(zhuǎn)加速度計(jì)式重力梯度儀輸出解調(diào)與濾波

楊 曄，李 達(dá)，高 巍

（天津航海儀器研究所，天津 300131）

重力梯度儀是對(duì)地球表面微小重力梯度變化進(jìn)行連續(xù)測(cè)量的儀器。由于核心敏感元件加速度計(jì)工藝與性能水平限制，以及多環(huán)節(jié)安裝誤差等因素導(dǎo)致系統(tǒng)實(shí)際輸出信號(hào)中包含了大量噪聲，且信噪比極低，為了能夠在強(qiáng)噪聲中有效提取真實(shí)的重力梯度信號(hào)，需在信號(hào)解調(diào)過程中降低諧波干擾引起的測(cè)量偏差。結(jié)合誤差產(chǎn)生機(jī)理，分析比較了不同的解調(diào)方法對(duì)重力梯度信號(hào)解調(diào)的影響，明確了合理的信號(hào)解調(diào)手段。同時(shí)，在數(shù)據(jù)處理過程中增加帶通濾波環(huán)節(jié)，進(jìn)一步降低動(dòng)態(tài)噪聲對(duì)系統(tǒng)的影響。結(jié)合動(dòng)態(tài)搖擺實(shí)驗(yàn)數(shù)據(jù)的仿真驗(yàn)證結(jié)果表明：上述措施對(duì)重力梯度儀原始觀測(cè)數(shù)據(jù)的噪聲抑制具有明顯效果，測(cè)量偏差在一定程度上得到解決，解調(diào)后信號(hào)的噪聲幅度下降至20%，提高了系統(tǒng)的空間分辨率，具有工程化應(yīng)用意義。

重力梯度儀；信號(hào)提??；帶通濾波；同步解調(diào)

重力梯度儀是用于測(cè)量重力場(chǎng)分布梯度的精密設(shè)備，基于Bell Aerospace公司提出的旋轉(zhuǎn)加速度計(jì)測(cè)量原理的重力梯度儀（后文中簡稱梯度儀）是迄今唯一實(shí)用的近地表動(dòng)態(tài)測(cè)量重力梯度儀[1-3]。礙于重力梯度信號(hào)的微弱和加速度計(jì)性能的限制，梯度儀輸出信號(hào)的信噪比極低，需后期數(shù)據(jù)處理進(jìn)一步提高測(cè)量精度[4-6]。由于世界范圍內(nèi)能夠從事該技術(shù)研究的國家和單位很少，相關(guān)技術(shù)文獻(xiàn)鮮有報(bào)道，本文結(jié)合科研工作和實(shí)測(cè)數(shù)據(jù)來驗(yàn)證方法的有效性，優(yōu)化重力梯度信號(hào)數(shù)據(jù)處理過程。

1 重力梯度測(cè)量原理與存在問題

1.1 測(cè)量原理

如圖1所示，旋轉(zhuǎn)加速度計(jì)式重力梯度儀基于加速度計(jì)位置差分測(cè)量原理，通過機(jī)械旋轉(zhuǎn)的方式將重力梯度張量分量調(diào)制到系統(tǒng)旋轉(zhuǎn)頻率的二倍頻處。加速度計(jì)輸出與重力梯度張量分量之間的關(guān)系可表示為

式中：a1、a2、a3、a4是4個(gè)加速度計(jì)的輸出；R是加速度計(jì)檢測(cè)質(zhì)心到旋轉(zhuǎn)中心的距離，Γxx、Γyy、Γxy是重力梯度張量分量；ω是旋轉(zhuǎn)裝置的旋轉(zhuǎn)角速度。

圖1 旋轉(zhuǎn)加速度計(jì)式重力梯度儀測(cè)量原理Fig.1 Measurement principle of rotating accelerometer gravity gradiometer

以往的信號(hào)解調(diào)方案如圖2所示，其中方波解調(diào)公式為

式中：Aout(i)是經(jīng)位置同步的加速度計(jì)輸出采樣值；N是重力梯度儀圓盤旋轉(zhuǎn)半周時(shí)采樣數(shù)據(jù)個(gè)數(shù)；R是旋轉(zhuǎn)半徑；π2是將解調(diào)值轉(zhuǎn)化為幅值的比例系數(shù)。

圖2 原始重力梯度數(shù)據(jù)處理框圖Fig.2 Flow chart for processing the original gravity gradiometer data

1.2 存在問題

梯度儀在進(jìn)行水平搖擺試驗(yàn)（搖擺幅度為2°，搖擺周期為20 s）時(shí)的典型解調(diào)出的信號(hào)如圖3所示，可以看出在4 892 s由靜態(tài)開始搖擺時(shí)，信號(hào)均值出現(xiàn)600E(1E=10-9s-2)的偏差“臺(tái)階”，并且原始信號(hào)信噪比極低，噪聲峰峰值高達(dá)10000E，相關(guān)數(shù)據(jù)見表1。

圖3 搖擺試驗(yàn)時(shí)梯度儀Γxx-Γyy解調(diào)信號(hào)（方波解調(diào)）Fig.3 Demodulation signal of gravity gradiometer Γxx-Γyyin swing test (square-wave demodulation)

表1 搖擺試驗(yàn)梯度信號(hào)Γxx-Γyy的均值和標(biāo)準(zhǔn)差(方波解調(diào))Tab.1 Mean values and standard deviations of gravity gradient signals in swing test (square-wave demodulation)

2 梯度儀數(shù)據(jù)處理方法

在此分析梯度儀在搖擺試驗(yàn)中出現(xiàn)的測(cè)量偏差的機(jī)理，系統(tǒng)可采用的解調(diào)方式主要有方波解調(diào)和正弦解調(diào)兩種，其中正弦解調(diào)的公式為

式中：ts是梯度儀的采樣時(shí)間間隔；分子上的2是將解調(diào)值轉(zhuǎn)化為幅值的比例系數(shù)。兩種解調(diào)方式在解調(diào)2ω頻率上的信號(hào)幅頻響應(yīng)對(duì)比如圖4所示。

圖4 方波解調(diào)和正弦解調(diào)方式幅頻響應(yīng)對(duì)比Fig.4 Comparison on amplitude-frequency responses between square-wave demodulation and sinusoidal demodulation

從圖4可以看出，在解調(diào)頻點(diǎn)附近，兩種解調(diào)方式效果相當(dāng)，但在解調(diào)頻率的奇數(shù)倍頻處，方波解調(diào)方式會(huì)引入高頻諧波項(xiàng)，而正弦解調(diào)方式的結(jié)果不被高頻諧波項(xiàng)影響。引入一個(gè)模擬梯度信號(hào)y= sin(2πf1t)+sin(2πf2t )進(jìn)行仿真驗(yàn)證（其中，f1=2ω，f2=6ω），該信號(hào)時(shí)域圖見圖5（T=12ω），解調(diào)結(jié)果對(duì)比見表2。

圖5 模擬梯度信號(hào)時(shí)域圖Fig.5 Time domain of simulated gravity gradiometer signal

表2 兩種解調(diào)方式模擬梯度信號(hào)解調(diào)結(jié)果對(duì)比Tab.2 Comparison of two demodulation methods for simulated gradient signals

由此可知，方波解調(diào)值會(huì)以一定比例引入高頻諧波分量，由于搖擺試驗(yàn)中在機(jī)械旋轉(zhuǎn)調(diào)制頻率的高頻諧波頻點(diǎn)處存在較大噪聲耦合，方波解調(diào)的方式會(huì)將這些諧波頻點(diǎn)噪聲轉(zhuǎn)化為解調(diào)值，匯入梯度信號(hào)，導(dǎo)致梯度信號(hào)出現(xiàn)測(cè)量偏差，而正弦解調(diào)值不會(huì)受到高頻諧波噪聲的影響。從圖4可以看出，解調(diào)頻點(diǎn)附近的噪聲也會(huì)對(duì)解調(diào)結(jié)果產(chǎn)生影響，如果原始數(shù)據(jù)存在頻率差為Δf，幅值為A的噪聲項(xiàng)，其對(duì)解調(diào)值的影響為

式中：k為解調(diào)值序號(hào)。由式(5)可得解調(diào)頻點(diǎn)附近的噪聲會(huì)使解調(diào)結(jié)果出現(xiàn)等幅的余弦震蕩，增加梯度信號(hào)噪聲，因此需在解調(diào)前加入帶通濾波環(huán)節(jié)以提高解調(diào)精度。并且由式(4)可知，頻率差為Δf噪聲項(xiàng)會(huì)引起解調(diào)結(jié)果呈周期為1(2Δf)的余弦振蕩，因此解調(diào)后對(duì)梯度信號(hào)低通濾波器截止頻率的選擇與原始信號(hào)帶通濾波器的帶寬參數(shù)息息相關(guān)。如果帶通濾波器的帶寬參數(shù)為f，解調(diào)后對(duì)梯度信號(hào)的低通濾波器截止頻率必須低于f才有降噪效果。在實(shí)際工作中可結(jié)合梯度儀自身性能和外部環(huán)境噪聲水平綜合考慮是否對(duì)梯度信號(hào)進(jìn)行低通濾波。

3 實(shí)驗(yàn)數(shù)據(jù)處理與結(jié)果分析

將同樣的搖擺試驗(yàn)數(shù)據(jù)進(jìn)行正弦解調(diào)，梯度信號(hào)如圖6所示，相關(guān)數(shù)據(jù)見表3，從中可以看出，梯度信號(hào)在水平搖擺試驗(yàn)時(shí)的測(cè)量偏差可在一定程度上得到解決，偏差幅度從600 E降至20 E。

圖6 搖擺試驗(yàn)梯度儀Γxx-Γyy解調(diào)信號(hào)(正弦解調(diào)）Fig.6 Demodulation signal of gravity gradiometer Γxx-Γyyin swing test (sinusoidal demodulation)

表3 搖擺試驗(yàn)梯度信號(hào)Γxx-Γyy的均值和標(biāo)準(zhǔn)差(正弦解調(diào))Tab.3 Mean values and standard deviations of gravity gradient signals in swing test (sinusoidal demodulation)

為進(jìn)一步驗(yàn)證正弦解調(diào)方式的合理性，分別將兩種解調(diào)方式的解調(diào)值與原始數(shù)據(jù)在2ω(Hz)處的FFT（快速傅式變換）值作比較，其結(jié)果見表4。從表4中可以看出，方波解調(diào)得到的梯度信號(hào)與原始數(shù)據(jù)在2ω(Hz)處的FFT值相差較大，是由于引入的高頻諧波噪聲耦合項(xiàng)不穩(wěn)定導(dǎo)致的，而正弦解調(diào)的信號(hào)與FFT值一致，因此使用正弦解調(diào)方式得到的梯度信號(hào)更加接近真實(shí)值。

再在解調(diào)前對(duì)原始數(shù)據(jù)增加帶通濾波環(huán)節(jié)，根據(jù)梯度儀輸出信號(hào)特征并在考慮儀器測(cè)量信息空間分辨率的基礎(chǔ)上，在采樣頻率800 Hz的條件下設(shè)計(jì)一組合理的帶通濾波器參數(shù)：中心頻率為2ω(Hz)，窗口函數(shù)為凱塞窗，帶通帶寬為0.01 Hz，過渡帶帶寬為0.005 Hz，濾波器階數(shù)為64 000。將濾波后數(shù)據(jù)解調(diào)得到梯度信號(hào)如圖7所示，相關(guān)數(shù)據(jù)見表5?？梢妿V波后梯度信號(hào)降噪效果明顯，梯度信號(hào)噪聲峰峰值降為1000E左右，標(biāo)準(zhǔn)差降為原來的20%左右。

表4 搖擺試驗(yàn)梯度儀輸出信號(hào)兩種解調(diào)方式在2ω(Hz)處FFT值對(duì)比Tab.4 Comparison on FFT values of gradient signals at 2ω(Hz) by the two demodulation modes in swing test

圖7 帶通濾波后梯度張量Γxx-Γyy解調(diào)信號(hào)Fig.7 Band-pass filtered gravity gradient tensor Γxx-Γyydemodulation signal

表5 使用帶通濾波器后解調(diào)效果對(duì)比(未平滑處理)Tab.5 Comparison on demodulation effects after band-pass filtering (no smoothing)

在梯度儀實(shí)際工作中，即使在解調(diào)前增加帶通濾波環(huán)節(jié)，梯度信號(hào)噪聲仍然很大，仍需對(duì)其進(jìn)行平滑濾波處理。將通過帶通濾波器后解調(diào)得到的梯度信號(hào)與原始梯度信號(hào)分別使用 N秒平滑濾波并計(jì)算標(biāo)準(zhǔn)差，對(duì)比結(jié)果如圖8所示，相關(guān)數(shù)據(jù)見表6。

從圖8中可以看出，在使用帶通濾波處理后對(duì)解調(diào)出來的梯度信號(hào)使用100 s以下的平滑濾波處理，其對(duì)噪聲的抑制效果不明顯，如需進(jìn)一步降噪需進(jìn)行更長時(shí)間的平滑濾波處理，這是由于帶通濾波器的通帶帶寬參數(shù)為0.01 Hz造成的。從表6中還可以看出，僅使用帶通濾波處理后，解調(diào)的梯度信號(hào)均值會(huì)發(fā)生幅值為4E的數(shù)據(jù)偏移現(xiàn)象，但再經(jīng)過平滑濾波后偏移現(xiàn)象消失，這是由于在2ω(Hz)上存在信號(hào)串?dāng)_，這里可以對(duì)梯度信號(hào)使用平滑濾波的手段消除。

圖8 梯度信號(hào)N秒平滑曲線對(duì)比Fig.8 Comparison on curves of gradient signals before and after smoothing for N seconds

表6 重力梯度不同數(shù)據(jù)處理方法結(jié)果對(duì)比Tab.6 Comparison on different data processing methods for gravity gradient data

改進(jìn)后的重力梯度信號(hào)數(shù)據(jù)處理過程如圖 9所示。通過增加帶通濾波環(huán)節(jié)的手段可將梯度信號(hào)噪聲降至原來的20%左右，達(dá)到提高梯度儀測(cè)量信息空間分辨率的目的。將解調(diào)方式改為正弦解調(diào)可以有效地消除測(cè)量偏差，提高儀器測(cè)量分辨率。

圖9 改進(jìn)后重力梯度數(shù)據(jù)處理框圖Fig.9 Flow chart for processing the improved gravity gradiometer data

4 總 結(jié)

本文針對(duì)梯度儀輸出信號(hào)特征，提出一種重力梯度信號(hào)提取的有效方法：

1）使用帶通濾波器對(duì)梯度儀原始輸出信號(hào)進(jìn)行濾波，可將重力梯度信號(hào)頻點(diǎn)外的干擾噪聲濾掉；

2）對(duì)濾波后的信號(hào)用正弦解調(diào)代替方波解調(diào)，以消除諧波頻點(diǎn)處的噪聲引起的測(cè)量偏差；

3）將解調(diào)后的梯度信號(hào)進(jìn)行平滑濾波處理，進(jìn)一步來提高測(cè)量分辨率。

通過采用上述方法對(duì)梯度儀原始測(cè)量數(shù)據(jù)進(jìn)行離線處理，測(cè)量偏差可在一定程度上得到解決，同時(shí)信號(hào)噪聲標(biāo)準(zhǔn)差可有效降低到原來的20%，噪聲峰峰值降低到1000E，提高了梯度儀的測(cè)量精度和測(cè)量信息的空間分辨率，為實(shí)現(xiàn)高精度動(dòng)態(tài)重力梯度測(cè)量提供了一種數(shù)據(jù)處理的方法。

（References）：

[1] Roberts D R, Chowdhury P J, Lowe S, et al. Airborne gravity gradiometer surveying of petroleum systems under Lake Tanganyika, Tanzania[C]//ASEG-PESA2015. 2015: 1-4.

[2] DiFrancesco D, Meyer T, Christensen A, et al. Gravity gradiometry—today and tomorrow[C]//11thSAGA Biennial Technical Meeting and Exhibition. Swaziland, 2009: 80-83.

[3] Dransfield M H, Christensen A N. Performance of airbornegravity gadiometers[J]. The Leading Edge, 2013, 32(8): 908-922.

[4] Lee J B. FALCON gravity gradiometer technology[J]. Exploration Geophysics, 2001, 32(4): 247-250.

[5] Christensen A N, Dransfield M H, Galder C V. Noise and repeatability of airborne gravity gradiometry[J]. First Break, 2015, 33: 55-63.

[6] 蔡體菁, 錢學(xué)武, 丁昊. 旋轉(zhuǎn)加速度計(jì)重力梯度儀重力梯度信號(hào)仿真[J]. 物探與化探, 2015, 39(S1): 76-79. Cai Ti-jing, Qian Xue-wu, Ding Hao. Signal simulation of gravity gradiometer of rotating accelerometer[J]. Geophysical and Geochemical Exploration, 2015, 39(S1): 76-79.

[7] 劉昱, 張宇, 潘軍, 等. 重力梯度測(cè)量中光柵角編碼器安裝誤差的影響及分析[J]. 中國慣性技術(shù)學(xué)報(bào), 2015, 23(1): 6-8. Liu Yu, Zhang Yu, Pan Jun, et al. Influence and analysis of optical angular encoder’s installation error in gravity gradiometry[J]. Journal of Chinese Inertial Technology, 2015, 23(1): 6-8.

[8] O’keefe G J, Lee J B, Turner R J, et al. Gravity gradiometer. US: No.5922951[P]. 1999-10-05.

[9] Van Leeuwen E H, McCracken K G, Lee J B. Airborne gravity gradiometers. US: No.6883372 B2[P]. 2005.

[10] 錢學(xué)武, 蔡體菁, 尹航. 一種旋轉(zhuǎn)加速度計(jì)重力梯度儀重力梯度解調(diào)方法[J]. 物探與化探, 2015, 39(S1): 80-83. Qian Xue-wu, Cai Ti-jing, Yin Hang. A gravity gradient demodulation method for rotation accelerometer gravity gradiometer[J]. Geophysical and Geochemical Exploration, 2015, 39(S1): 80-83.

[11] Moody M V, Paik H J. A superconducting gravity gradiometer for inertial navigation[C]//2004 IEEE Position Location and Navigation Symposium. 2004: 775-781.

[12] 錢學(xué)武, 蔡體菁. 旋轉(zhuǎn)加速度計(jì)重力梯度儀加速度計(jì)標(biāo)度因數(shù)實(shí)時(shí)反饋調(diào)整方法[J]. 中國慣性技術(shù)學(xué)報(bào), 2016, 24(2): 148-153. Qian Xue-wu, Cai Ti-jing. Real-time feedback adjusting methods for accelerometer scale factor in gravity gradiometer for rotating accelerometer[J]. Journal of Chinese Inertial Technology, 2016, 24(2): 148-153.

Gradient signal demodulation and filtering for rotating accelerometer gravity gradiometer

YANG Ye, LI Da, GAO Wei
(Tianjin Navigation Instrument Research Institute, Tianjin 200131, China)

In rotating accelerometer gravity gradiometer, due to the limitation of the accelerometer’s process and performance level and the effects of various installation errors, there are lots of noises in the system output signal, and the SNRs are extremely low. In order to effectively extract the real gravity gradient signals in strong noises, the measurement errors caused by harmonic interference must be reduced in the process of signal demodulation. Combined with the error generation mechanism, the influences of different demodulation methods on the gravity gradient signal demodulation are analyzed and compared, and a reasonable method of signal demodulation is chosen. Meanwhile, a band-pass filter is added in the data processing process to further reduce the influence of the dynamic noise. Combined with dynamic swing test data, the simulations show that the above measures have obvious effect on the noise suppression of the original gravity data, and the noise amplitude of the demodulated signal is reduced to 20%, which improves the spatial resolution of the system, and has the engineering application value.

gravity gradiometer instrument; signal extraction; band-pass filter; synchronous demodulation

U666.1

：A

2016-10-02；

：2016-11-20

國家高技術(shù)研究發(fā)展計(jì)劃（863計(jì)劃）（2011AA060501）；總裝預(yù)研項(xiàng)目（51309040702）

楊曄（1968—），男，研究員，碩士生導(dǎo)師，從事導(dǎo)航、制導(dǎo)與控制方面研究工作。E-mail: Liuyuyangye@eyou.com

1005-6734(2016)06-0701-05

10.13695/j.cnki.12-1222/o3.2016.06.001