馬國(guó)慶,孟令順,杜曉娟,張鳳旭

(吉林大學(xué)地球探測(cè)科學(xué)與技術(shù)學(xué)院,吉林長(zhǎng)春 130026)

0 前言

在重、磁數(shù)據(jù)的處理中,往往在處理前,首先要對(duì)原始數(shù)據(jù)進(jìn)行擴(kuò)邊處理。常用的擴(kuò)邊方法有余弦衰減至零和對(duì)折擴(kuò)邊方法,這二種方法都存在明顯的邊界效應(yīng)[1],而邊界效應(yīng)不利于對(duì)異常進(jìn)行解釋。為了減小邊界效應(yīng),作者在本文提出了三方向擴(kuò)邊方法。

磁測(cè)數(shù)據(jù)受外界因素影響較大,即使一些小的磁性干擾,都會(huì)給實(shí)際數(shù)據(jù)的處理和解釋帶來困難,所以必須利用相應(yīng)的處理方法進(jìn)行消除。中值濾波具有消除突跳點(diǎn)和保護(hù)邊緣的特征[2、3],作者在本文對(duì)中值濾波方法進(jìn)行了改進(jìn),提出了優(yōu)化中值濾波方法。實(shí)踐證明,利用該方法對(duì)磁測(cè)數(shù)據(jù)進(jìn)行處理,能夠取得較好的效果。

1 三方向擴(kuò)邊方法

1.1 三方向擴(kuò)邊方法的原理

實(shí)際測(cè)量中,某一測(cè)點(diǎn)的數(shù)據(jù)不僅受到目標(biāo)地質(zhì)體的影響,還受到來自周圍一定范圍內(nèi)的其它地質(zhì)體的影響,因此該點(diǎn)的數(shù)據(jù)是一個(gè)綜合反映[4]。三方向擴(kuò)邊方法的原理是對(duì)一個(gè)點(diǎn)的三方向同時(shí)進(jìn)行擴(kuò)邊,取擴(kuò)邊值的平均值作為最后的擴(kuò)邊結(jié)果。這樣做可綜合周圍地質(zhì)體的影響,減小邊界效應(yīng)。

1.2 三方向擴(kuò)邊方法計(jì)算步驟

進(jìn)行縱向擴(kuò)邊,對(duì)于D點(diǎn)(見下頁圖1),首先分別計(jì)算出Y2,6→D,Y2,9→D,Y2,12→D三個(gè)方向上相鄰點(diǎn)的異常差平均值;然后根據(jù)邊界值和異常差,分別在三個(gè)方向上計(jì)算出D點(diǎn)的值;最后取三個(gè)方向上D值的平均值,做為D點(diǎn)的擴(kuò)邊值。B、C二點(diǎn)因?yàn)樵跀U(kuò)邊數(shù)據(jù)的邊角處,因此采用特殊的擴(kuò)邊方法。對(duì)于B點(diǎn),由于只在三個(gè)方向當(dāng)中的二個(gè)方向上存在原始數(shù)據(jù),所以只利用這二個(gè)方向的邊界值和差值來計(jì)算出最后結(jié)果。對(duì)于C點(diǎn),將采用“二長(zhǎng)一短”的方法來進(jìn)行計(jì)算,因?yàn)樵谌齻€(gè)方向當(dāng)中有二個(gè)方向上的數(shù)據(jù)量充足,而另一個(gè)方向的數(shù)據(jù)量不足,最后的結(jié)果也是三個(gè)方向的結(jié)果的平均。在橫向也采用上述方法進(jìn)行擴(kuò)邊。對(duì)于A點(diǎn),因?yàn)槭墙稽c(diǎn),所以在橫向、縱向擴(kuò)邊計(jì)算完畢以后再對(duì)其進(jìn)行擴(kuò)邊。

圖1 三方向擴(kuò)邊方法Fig.1 Three directions of expansion methods

三方向擴(kuò)邊方法的計(jì)算公式為:

當(dāng)i=0時(shí)為D點(diǎn)的計(jì)算公式:

其余點(diǎn)的計(jì)算公式可以依據(jù)公式(3)來獲得。

1.3 模型數(shù)據(jù)驗(yàn)算

理論模型:一個(gè)平坦的地區(qū)存在二個(gè)位置不同,深度不同的垂直磁化的球體,所產(chǎn)生的疊加磁異常見圖2。

圖2 理論模型計(jì)算結(jié)果Fig.2 The results of theoreticalmodel calculations

采用對(duì)折擴(kuò)邊方法和三方向擴(kuò)邊方法,分別對(duì)理論數(shù)據(jù)進(jìn)行擴(kuò)邊,如圖3所示。

模型試算表明,作者在本文提出的三方向擴(kuò)邊法,基本上不會(huì)引起邊界效應(yīng),并且能使邊界信息得到很好的保留。

圖3 不同擴(kuò)邊方法的結(jié)果對(duì)比Fig.3 The comparison of differentmethods of edge expansion

2 優(yōu)化中值濾波在磁法數(shù)據(jù)處理中的應(yīng)用

在實(shí)際測(cè)量中,當(dāng)外界干擾較大,對(duì)異常的處理和解釋產(chǎn)生較大的影響時(shí),必須想辦法將干擾消除掉。

中值濾波是常用的一種非線性平滑濾波。它是一種領(lǐng)域運(yùn)算,類似于卷積,但不是加權(quán)求和計(jì)算,而是把領(lǐng)域中的像素按灰度等級(jí)進(jìn)行排序,然后選擇該組的中間值作為輸出像素值。它能減弱或消除傅里葉空間的高頻分量,并且影響低頻分量。因?yàn)楦哳l分量對(duì)應(yīng)圖像中的區(qū)域邊緣的灰度值,具有較大較快變化的部份,所以該濾波可將這些分量濾除,使圖像平滑[5]。

作者在本文以中值濾波為基礎(chǔ),結(jié)合干擾的特征以及異常的屬性,提出了優(yōu)化中值濾波的概念。面積性優(yōu)化中值濾波是以環(huán)帶的形式來進(jìn)行的(見圖4),環(huán)帶的多少由干擾異常的范圍來決定。

圖4中的矩形區(qū)域?yàn)楦蓴_區(qū),數(shù)字“1”代表的是第一個(gè)環(huán)帶(包括邊界的數(shù)據(jù))內(nèi)的數(shù)據(jù);數(shù)字“2”環(huán)指的是第二個(gè)環(huán)帶(包括邊界數(shù)據(jù))與第一個(gè)環(huán)帶(不包括邊界數(shù)據(jù))之間的數(shù)據(jù);數(shù)字“3”環(huán)及數(shù)字“4”環(huán)的數(shù)據(jù)獲得方式,與數(shù)字“2”環(huán)是一致的。

中值濾波的點(diǎn)數(shù)不同,處理后的結(jié)果反映的信息也不同。利用這種形式來進(jìn)行中值濾波的優(yōu)點(diǎn)是:

(1)由于有三個(gè)點(diǎn)以上的異常才稱為異常,因此選擇一環(huán)的目的是去除突跳點(diǎn),在無干擾的區(qū)域不會(huì)使異常幅值變小。

(2)二環(huán)、三環(huán)是一個(gè)過渡區(qū)域,因?yàn)樵谀承﹨^(qū)域干擾的范圍可能會(huì)小一些,所以加這二個(gè)環(huán)帶的目的,是為了去除范圍小的干擾。

圖4 優(yōu)化中值濾波示意圖Fig.4 Schematic diagram of filter optimization

(3)四環(huán)是主要的濾波區(qū)域,在這個(gè)區(qū)域內(nèi)無干擾的數(shù)據(jù)多于干擾數(shù)據(jù),對(duì)這個(gè)環(huán)帶內(nèi)的數(shù)據(jù)進(jìn)行中值濾波,能夠較好地去除干擾。

對(duì)于剖面濾波是以不同大小的窗口對(duì)曲線進(jìn)行濾波。

對(duì)數(shù)據(jù)進(jìn)行濾波之后,根據(jù)其異常和干擾的特點(diǎn),加上相應(yīng)的限制條件會(huì)使處理結(jié)果更可靠。

3 實(shí)際數(shù)據(jù)驗(yàn)算

以四川省綿陽市射洪縣地區(qū)的一條日變曲線為例,來驗(yàn)證優(yōu)化中值濾波方法的實(shí)用性。

從圖5中可以看出,在“1”處、“2”處由于受到人為的影響而出現(xiàn)蹦跳。在理論上,一天的日變曲線在沒有大小磁暴的干擾時(shí),應(yīng)為一條比較平滑的曲線,因此需要將日變干擾造成的蹦跳改正掉。圖6(見下頁)是在不同窗口下,對(duì)實(shí)測(cè)日變曲線進(jìn)行中值濾波后的結(jié)果。

圖5 四川省綿陽市射洪縣地區(qū)日變曲線(2009年4月27日)Fig.5 Diurnal variation curve on Shehong county area,Mianyang,Sichuan province

從圖6中可以看出,不同窗口的中值濾波,對(duì)“2”號(hào)干擾都起到了作用。但是在“1”號(hào)干擾處,三個(gè)窗口和五個(gè)窗口的中值濾波效果明顯不如七個(gè)窗口的中值濾波,但是七個(gè)窗口的中值濾波在極值處存在一定壓制。

見圖7(見下頁),結(jié)合不同窗口對(duì)數(shù)據(jù)的處理效果及干擾為負(fù)值的特點(diǎn),將以上三個(gè)窗口濾波后的數(shù)據(jù)進(jìn)行比較,選取最大值來作為最終的結(jié)果。

從圖7中可以看出,優(yōu)化中值濾波可以很好地去掉磁測(cè)中的偶然誤差,這種方法可以應(yīng)用到實(shí)測(cè)數(shù)據(jù)處理中。

圖8(見下頁)是內(nèi)蒙古呼圖格地區(qū)的實(shí)測(cè)平面磁法數(shù)據(jù),以此為例來驗(yàn)證擴(kuò)邊方法和優(yōu)化中值濾波的聯(lián)合應(yīng)用效果。

圖8中左側(cè)的負(fù)異常條帶,是由高壓線引起的。為了對(duì)資料進(jìn)行準(zhǔn)確的地質(zhì)解釋,必須將這種干擾去除掉。利用優(yōu)化中值濾波對(duì)其進(jìn)行處理的步驟[6]如下頁圖9所示。

由原始圖像可知,干擾異常為負(fù)值且異常寬度為1點(diǎn)～4個(gè)點(diǎn),一環(huán)、二環(huán)、三環(huán)帶的作用主要是去掉突跳點(diǎn)和去除范圍小的干擾,四環(huán)帶對(duì)磁干擾條帶有很好的去除作用,但是對(duì)極值有壓制。由于本地區(qū)進(jìn)行磁法測(cè)量的主要目的是尋找鐵磁性物質(zhì),因此在中值濾波之后,取四個(gè)值的最大值做為最終的計(jì)算結(jié)果(見下頁圖10)。

圖11(見下頁)是采用與優(yōu)化中值濾波同樣大小窗口的方域中值濾波。

從圖10、圖11二幅圖的對(duì)比中可以看出,利用三方向擴(kuò)邊方法,在邊界區(qū)域未出現(xiàn)邊界效應(yīng),所以這種擴(kuò)邊方法實(shí)用性較好。

從圖10與圖11的對(duì)比中可以看出,優(yōu)化中值濾波的效果比普通中值濾波的效果較好,在將干擾很好的去除的前提下,并沒有使異常的極值降低很多。從原始圖像中我們可以發(fā)現(xiàn),在負(fù)值區(qū)域內(nèi),坐標(biāo)為(16,42)的位置存在一個(gè)小的正異常,普通的方域中值濾波將這個(gè)異常給模糊掉了,而優(yōu)化中值濾波將其完好的保留。在優(yōu)化中值濾波使異常的處理與解釋更可靠,因此這種方法實(shí)用性較高。

圖6 不同窗口所做的中值濾波結(jié)果對(duì)比Fig.6 Comparative results of median filter in different windows

圖7 經(jīng)過濾波后的日變曲線與實(shí)測(cè)結(jié)果的對(duì)比Fig.7 Comparison of diurnal variation curve and measured results after filter

圖8 內(nèi)蒙古呼圖格地區(qū)的實(shí)測(cè)磁異常平面圖Fig.8 The measured magnetic map of Hutuge area,Neimenggu

圖9 優(yōu)化中值濾波的計(jì)算步驟Fig.9 Calculation steps of optimization median filter

4 結(jié)論

從實(shí)際資料的處理中可以看出,作者在本文提出的三方向擴(kuò)邊方法不會(huì)引起邊界效應(yīng),并且還可以很好地保留邊界信息。從圖像的對(duì)比中我們會(huì)看出,優(yōu)化中值濾波比普通中值濾波效果好,不僅能將干擾去除,而且還不會(huì)降低幅值,對(duì)信息有很好的保留效果。因此以上方法可以在實(shí)際數(shù)據(jù)的處理中得以廣泛應(yīng)用。

圖10 優(yōu)化中值濾波后的結(jié)果Fig.10 Results of optimization Median filter

圖11 方域中值濾波后的結(jié)果Fig.11 Results of domain median filter

[1] 段本春,徐世浙.磁(重力)異常局部場(chǎng)與區(qū)域場(chǎng)分離處理中的擴(kuò)邊方法研究[J].物探化探計(jì)算技術(shù),1997,19(4):299.

[2] 劉財(cái),王典,劉洋,等.二維多級(jí)中值濾波技術(shù)在隨機(jī)噪聲消除中的應(yīng)用初探[J].石油地球物理勘探,2005,40(2):163.

[3] 邵 敏,何展翔,張立恩.二維中值濾波在電磁測(cè)深曲線去噪中的應(yīng)用[J].物探化探計(jì)算技術(shù),2001,23(3):226.

[4] 劉東甲,程方道.劃分重力區(qū)域場(chǎng)與局部場(chǎng)的多次切割法[J].物探化探計(jì)算技術(shù),1997,19(1):32.

[5] SPECTOR A.Comments and errata to classic papers[J].Geophysics,1985,50(11):2278.

[6] ARCE G R,FOSTER R E.Detail-preserving rankedor-der based filter for image processing[J].IEEE Transon ASSP,1987,37(1):83.