劉 備, 邊少鋒, 紀(jì) 兵,*, 賢鵬飛

(1. 海軍工程大學(xué)電氣工程學(xué)院, 湖北 武漢 430033; 2. 廣西空間信息與測繪重點(diǎn)實(shí)驗(yàn)室, 廣西 桂林 530001)

0 引 言

海洋重力測量作為地球重力場研究的一個重要組成部分，為海洋資源勘探、軍事應(yīng)用、基礎(chǔ)海洋研究等領(lǐng)域提供了重要的數(shù)據(jù)支撐[1-3]。目前，海洋重力測量的主要方法包括衛(wèi)星重力測量、航空重力測量和船載重力測量，船載海洋重力測量是獲取高精度以及高分辨率海洋重力場信息的重要方法之一[4-8]。常規(guī)的重力測量包括靜態(tài)重力測量和動態(tài)重力測量。但由于海洋環(huán)境的特殊性,靜態(tài)重力測量無法滿足當(dāng)前海洋開發(fā)和未來海戰(zhàn)的需求。動態(tài)重力測量以其能夠快速獲取高精度、高分辨率的優(yōu)勢成為獲取海洋重力場信息的主要手段之一[9-10]。動態(tài)重力測量獲取的高精度、高分辨率海洋重力場信息能夠?yàn)楹Ｑ箝_發(fā)、海底地形反演、未來海戰(zhàn)等提供一些關(guān)鍵信息[11-12]。衛(wèi)星重力測量作為動態(tài)重力測量的主要手段之一,其空間分辨率目前可以達(dá)到100 km,測量精度可達(dá)到2 mGal,若需要更高精度的海洋重力場測量成果,則需要通過船載重力測量方式獲取重力場數(shù)據(jù)對衛(wèi)星重力場數(shù)據(jù)進(jìn)行補(bǔ)充[13-15]。

船載重力測量由于其測量環(huán)境的特殊性,其數(shù)據(jù)處理成為制約測量成果精度的關(guān)鍵問題之一。對于高精度海洋重力測量數(shù)據(jù)處理與分析,文獻(xiàn)[2-3]將Chekan-AM重力儀搭載在不同測量船只上進(jìn)行測量任務(wù),并對測量結(jié)果進(jìn)行處理分析,搭載在兩艘測量船上的測量精度分別為0.2 mGal、0.6 mGal,并利用測量結(jié)果給出了波羅的海地區(qū)局部大地水準(zhǔn)面構(gòu)建的初步結(jié)果。Hunegnaw等[16]對北大西洋地區(qū)的重力測量數(shù)據(jù)進(jìn)行了處理分析,結(jié)果為交叉點(diǎn)處的標(biāo)準(zhǔn)差為4.03 mGal,經(jīng)平差處理后提高到1.5 mGal。蔡劭琨等[17]將經(jīng)驗(yàn)?zāi)B(tài)分解用于動態(tài)重力測量數(shù)據(jù)處理當(dāng)中,將重力異常的內(nèi)符合精度從0.88 mGal提高到了0.611 mGal。董慶亮等[18]對丹麥技術(shù)大學(xué)重力數(shù)據(jù)用于檢查海洋重力測量成果的可行性進(jìn)行了研究,對海洋重力測量數(shù)據(jù)的檢核方法起到了補(bǔ)充作用。吳燕雄等[8]以激光干涉絕對重力儀為基礎(chǔ)設(shè)計了一套船載絕對重力儀測量系統(tǒng),并給出了該系統(tǒng)工作的誤差修正方法、動態(tài)限制條件和修正精度,經(jīng)修正后該系統(tǒng)測量精度優(yōu)于±1.1 mGal。黃謨濤等[19]對?？罩亓y量誤差進(jìn)行深入分析,提出了一種適用于補(bǔ)償各類海空重力儀動態(tài)效應(yīng)剩余影響的通用模型,經(jīng)驗(yàn)證,該模型將海洋重力測量的內(nèi)符合精度從±9.35 mGal提高到±1.01 mGal。

對于船載海洋重力測量來說,外界環(huán)境會產(chǎn)生大量的噪聲。如何減弱噪聲、消除噪聲對測量成果的影響,提高測量精度,具有非常重要的意義[20-22]。目前常用的降噪方法是采用有限沖激響應(yīng)(finite impulse response, FIR)低通濾波器,FIR可以降低高頻噪聲,但濾波器最優(yōu)參數(shù)的普適性不高[23-24]。孫鶴泉等[25]將最大重疊離散小波變換(maximum overlap discrete wavelet transform, MODWT)用于海洋重力觀測航行數(shù)據(jù)的降噪中,并對該方法的可靠程度進(jìn)行了驗(yàn)證。基于小波變換具有適用性高的優(yōu)點(diǎn),本文將小波變換方法引入到船載重力場數(shù)據(jù)處理分析中,并利用實(shí)測數(shù)據(jù)驗(yàn)證小波變換方法用于船載重力場數(shù)據(jù)降噪中的精度和可靠性。

1 船載重力測量

本文用于驗(yàn)證實(shí)驗(yàn)的船載海洋重力場數(shù)據(jù)是將dgM1型船載重力儀搭載到廣海局“海洋四號”測量船上在某海域開展重力實(shí)驗(yàn)測量任務(wù)所測得的重力場數(shù)據(jù)。

dgM1型船載重力儀(見圖1)是一款針對船載平臺的動態(tài)重力儀,系統(tǒng)由重力儀本體、不間斷電源(uninterrupted power supply, UPS)、全球?qū)Ш叫l(wèi)星系統(tǒng)(global navigation satellite system, GNSS)接收機(jī)、顯控設(shè)備/實(shí)時顯控軟件等組成。其中,重力傳感器采用高精度石英撓性加速度計,姿態(tài)測量傳感器采用高精度光學(xué)陀螺。重力儀在線實(shí)時完成數(shù)據(jù)采集與處理,通過可視化軟件實(shí)時顯示重力異常、位置、速度、姿態(tài)等信息。

圖1 dgM1重力儀Fig.1 Gravimeter dgM1

該次重力實(shí)驗(yàn)測量任務(wù)開始前,測量人員在珠江口外海進(jìn)行了儀器精度測試,完成了重復(fù)線和交叉點(diǎn)測量,在確保儀器精度滿足要求之后,測量船再去測區(qū)開展作業(yè),測線圖如圖2所示。

圖2 測線圖Fig.2 Tracks of gravimetry

測線包括東西測線、南北測線、45°和135°斜線。測線間距5′。其中L1、L3、Z1、Z2、B1、D1測線各有兩條重復(fù)線。

2 數(shù)據(jù)處理方法及評價指標(biāo)

2.1 小波變換理論

小波變換的出現(xiàn)彌補(bǔ)了傅里葉變換中的不足。小波是一種長度有限、平均值為0的波形,主要特點(diǎn)有:在時域具有緊支性或近似緊支性;直流分量為0。小波變換的數(shù)學(xué)公式可以表示如下:

(1)

常用的小波基函數(shù)有Haar小波、多貝西(Daubechies,Db)小波、Symlets小波、墨西哥小波(Mexican hat wavelet, Mexh)和雙正交樣條小波(biorthogonal spline wavelet, BSW)等。本文選擇了Db小波和Haar小波兩種基函數(shù)用于數(shù)據(jù)的變換處理,并對變換結(jié)果進(jìn)行對比分析。

2.2 誤差比率

誤差比率是指先將相鄰重力異常值作差,然后再將兩個相鄰的差值相比得到的數(shù)值,具體計算如下:

(2)

誤差比率可以顯示重力異常差值的變化趨勢及規(guī)律。本文選取誤差比率作為指標(biāo)旨在驗(yàn)證信噪比較低的數(shù)據(jù)經(jīng)過小波降噪后是否可以保留原始數(shù)據(jù)的變化規(guī)律及趨勢。

3 數(shù)據(jù)處理

3.1 漂移檢查

重力儀固有的一大缺點(diǎn)就是存在零點(diǎn)漂移問題。但是,只要漂移變化的幅度不大,并且存在一定的規(guī)律性,漂移就可以通過一定的方法被檢測出來并加以改正。通過對數(shù)據(jù)處理分析得出該航次重力儀漂移量,如表1所示,儀器的日漂移為0.14 mGal,月漂移為4.2 mGal,該船次重力測量總漂移量為5.71 mGal。

3.2 測線內(nèi)符合精度

對重復(fù)測線計算其內(nèi)符合精度,計算公式如下:

(3)

式中:Mrms為重復(fù)測線測量值均方根誤差;δgi1為重復(fù)測線的第一次測量的第i個點(diǎn)的值;δgi2為重復(fù)測線的第二次測量第i個點(diǎn)的值。具體各測線的內(nèi)符合精度如表2所示。

表2 各測線內(nèi)符合精度表

續(xù)表2

通過各測線的標(biāo)準(zhǔn)差統(tǒng)計結(jié)果可以看出dgM1重力儀的測量精度在海況良好的情況下可以滿足科考和其他海洋勘探任務(wù)的測量精度要求。

3.3 數(shù)據(jù)降噪

本文進(jìn)行數(shù)據(jù)降噪處理是對所有測線總體數(shù)據(jù)進(jìn)行整體處理,小波變換方法選取了Db和Haar兩種小波基函數(shù)。Db小波選擇的是Db8小波基,分解層數(shù)為8層;Haar小波分解層數(shù)同樣為8層。降噪處理結(jié)果如圖3～圖10以及表3、表4所示。圖3中，數(shù)據(jù)中添加的5倍隨機(jī)噪聲服從均勻分布。圖4中添加的正態(tài)分布隨機(jī)噪聲屬于一般正態(tài)分布。圖7中，Db-R1為5倍隨機(jī)噪聲與原始數(shù)據(jù)之差；Db-R11為經(jīng)Db小波變換后的5倍隨機(jī)噪聲與原始數(shù)據(jù)之差；Db-R2為正態(tài)分布隨機(jī)噪聲與原始數(shù)據(jù)之差；Db-R22為經(jīng)Db小波變換后的正態(tài)分布隨機(jī)噪聲與原始數(shù)據(jù)之差。圖8中，H-R1為5倍隨機(jī)噪聲與原始數(shù)據(jù)之差；H-R11為經(jīng)Haar小波變換后的5倍隨機(jī)噪聲與原始數(shù)據(jù)之差；H-R2為正態(tài)分布隨機(jī)噪聲與原始數(shù)據(jù)之差；H-R22為經(jīng)Haar小波變換后的正態(tài)分布隨機(jī)噪聲與原始數(shù)據(jù)之差。

圖3 原始數(shù)據(jù)與隨機(jī)誤差Fig.3 Original data and random error

圖4 正態(tài)分布隨機(jī)誤差Fig.4 Normal distribution random error

圖5 Db8小波變換Fig.5 Db8 wavelet transform

圖6 Haar小波變換Fig.6 Haar wavelet transform

圖7 Db8小波變換后誤差對比圖Fig.7 Error diagram of after Db8 wavelet transform

圖8 Haar小波變換后誤差對比圖Fig.8 Error diagram of after Haar wavelet transform

圖9 Db8小波變換后誤差比率圖Fig.9 Error ratio diagram of Db8 wavelet transform

圖10 Haar小波變換后誤差比率圖Fig.10 Error ratio diagram of Haar wavelet transform

表3 標(biāo)準(zhǔn)差統(tǒng)計表

表4 相對精度統(tǒng)計表

表5 誤差比率統(tǒng)計表

對比原始數(shù)據(jù)、噪聲數(shù)據(jù)以及經(jīng)過降噪處理的數(shù)據(jù)可以發(fā)現(xiàn),加入5倍隨機(jī)噪聲和正態(tài)分布隨機(jī)噪聲的數(shù)據(jù)經(jīng)過Db8小波變換處理后數(shù)據(jù)絕對精度分別提高了1.024 2 mGal和9.274 7 mGal;采用Haar小波變換后,分別提高了0.804 4 mGal和9.109 5 mGal。通過圖5、圖6以及表3可知,經(jīng)過小波變換處理后的數(shù)據(jù)相對精度至少提高56%,隨著噪聲的增大,最高可提高92%。對圖9、圖10以及表4中的誤差比率變化情況進(jìn)行分析,可以看出經(jīng)過降噪后的重力異常數(shù)據(jù)能夠很好地保持?jǐn)?shù)據(jù)的基本變化趨勢不變,并且經(jīng)過小波變換后的數(shù)據(jù)視圖更加簡潔,能夠較為輕松地判別測量過程中是否存在異常體導(dǎo)致的測量結(jié)果突變,為利用重力異常進(jìn)行探測提供判別依據(jù)。

4 結(jié) 論

本文研究了小波變換方法在船載重力測量數(shù)據(jù)降噪處理與分析中的應(yīng)用,將Db小波和Haar小波兩種小波基函數(shù)引入到重力場數(shù)據(jù)降噪中,并利用dgM1重力儀實(shí)測數(shù)據(jù)驗(yàn)證降噪方法的精度與可靠性,結(jié)果表明小波變換方法可以很好地提高數(shù)據(jù)精度,加入5倍隨機(jī)誤差和正態(tài)分布隨機(jī)誤差的數(shù)據(jù)經(jīng)過Db8小波變換處理后數(shù)據(jù)絕對精度分別提高了1.024 2 mGal和9.274 7 mGal,相對精度分別提高了71%和93%;采用Haar小波變換后,分別提高了0.804 4 mGal和9.109 5 mGal,相對精度分別提高了56%和91%。對誤差比率變化情況進(jìn)行分析可以推斷,經(jīng)過小波變換方法降噪后的數(shù)據(jù)可以保留重力異常變化趨勢,并且更加簡潔,能夠?yàn)橹亓μ綔y判別是否存在異常體提供判別依據(jù)。