沈恒旭，李才明

?

Tikhonov迭代正則化向下延拓應(yīng)用探討

沈恒旭，李才明

(成都理工大學(xué) 地球物理學(xué)院，四川 成都 610059)

重磁位場(chǎng)延拓廣泛應(yīng)用于地質(zhì)勘探、能源勘探、軍事導(dǎo)航等領(lǐng)域的重磁資料數(shù)據(jù)處理中。Tikhonov迭代正則化向下延拓方法實(shí)現(xiàn)快速高精度向下延拓，能夠有效壓制高頻振蕩，實(shí)現(xiàn)過場(chǎng)源體而不奇異。但該方法對(duì)地質(zhì)體的響應(yīng)程度及對(duì)于地質(zhì)體物性的反映尚需研究，為了探討該方法對(duì)不同地質(zhì)體的響應(yīng)情況，以Tikhonov迭代正則化向下延拓算法為基礎(chǔ)，深入討論該方法對(duì)于不同地質(zhì)體向下延拓的結(jié)果。并以多個(gè)不同深度的球體以及板狀體等模型模擬不同地質(zhì)條件，借此來討論延拓結(jié)果對(duì)地質(zhì)體的響應(yīng)情況。理論場(chǎng)和實(shí)測(cè)資料的向下延拓成果表明，該方法能夠有效反映地質(zhì)體埋深，也能校準(zhǔn)確區(qū)分水平和垂直疊加異常體。

位場(chǎng)；Tikhonov迭代正則化；向下延拓；場(chǎng)源體

(InstituteofGeophysics,ChengduUniversityofTechnology,ChengduSichuan610059,China)

1 引 言

重磁位場(chǎng)向下延拓是重磁數(shù)據(jù)處理以及反演解釋中一種重要的數(shù)據(jù)處理方法。通過位場(chǎng)的向下延拓可以獲得更多的關(guān)于重磁位場(chǎng)的數(shù)據(jù)，從而更加有利于對(duì)地下介質(zhì)進(jìn)行定性定量分析。位場(chǎng)向下延拓可以突出疊加在區(qū)域背景上的局部異常，使之盡量少受區(qū)域場(chǎng)的影響，而且可以“放大”某些在低緩異常中不夠明顯的異常特征(如拐點(diǎn)、極值點(diǎn)、零值點(diǎn)等)，有利于進(jìn)一步的解釋推斷[1]。位場(chǎng)向下延拓?cái)U(kuò)充了已有重磁位場(chǎng)的數(shù)據(jù)量，從而可以利用位場(chǎng)進(jìn)行輔助導(dǎo)航[2-5]以及進(jìn)行磁異常邊緣識(shí)別[6]。重磁位場(chǎng)的空間分布特征，具有較強(qiáng)的直觀性，是研究場(chǎng)源體的重要依據(jù)[7]。然而，從數(shù)學(xué)上我們知道，位場(chǎng)向下延拓屬于數(shù)學(xué)物理反問題，是一類不適定問題，邊界值ΔT(x,0)有極微小的變化時(shí)，所得的解可以變化很大[8]。為此，眾多學(xué)者致力于提高向下延拓精度和延拓速度。徐世浙院士提出了位場(chǎng)向下延拓的積分迭代法，能夠穩(wěn)定向下延拓20倍點(diǎn)距，提高了向下延拓精度[9,10]；陳龍偉等實(shí)現(xiàn)了Landweber迭代法位場(chǎng)向下延拓[11]；多位學(xué)者基于維納濾波[12]、窗函數(shù)[13]、Tikhonov正則化[14,15]、自適應(yīng)濾波[16]、有限元法[17]等進(jìn)行位場(chǎng)延拓或應(yīng)用；另外劉東甲等研究了位場(chǎng)向下延拓的波數(shù)域迭代法[18]。但是，迭代積分法雖然延拓精度較高，但是存在多次迭代而導(dǎo)致計(jì)算數(shù)據(jù)量較大，不適合處理數(shù)據(jù)量較大的數(shù)據(jù)。波數(shù)域延拓雖然計(jì)算數(shù)據(jù)量小，延拓速度快，但是延拓精度受限，并不能很好地實(shí)現(xiàn)位場(chǎng)向下延拓。曾小牛等[19]將迭代積分法的高精度延拓和波數(shù)域延拓的快速計(jì)算能力相結(jié)合，并根據(jù)傅初黎等提出的不適定問題的迭代Tikhonov正則化方法[20]，提出了迭代Tikhonov正則化向下延拓法，并對(duì)其波數(shù)域算子進(jìn)行了推導(dǎo)。本文以曾小牛等提出的迭代Tikhonov正則化向下延拓法為基礎(chǔ)，分別對(duì)不同的正演模型數(shù)據(jù)進(jìn)行向下延拓，從而探討該方法對(duì)不同地質(zhì)體、不同磁化方向的應(yīng)用效果，進(jìn)而為該方法的應(yīng)用提供參考。

2 Tikhonov迭代正則化向下延拓方法原理

如圖1所示，設(shè)z軸向下為正，場(chǎng)源位于z=h平面以下(h>0)，若z=0觀測(cè)平面上的位場(chǎng)T(x，0)為已知，則由T(x，0)求z=h平面上的位場(chǎng)數(shù)據(jù)T(x,h)稱為位場(chǎng)的向下延拓[21]。

圖1 向下延拓示意圖Fig.1 Schematic diagram of downward continuation

(1)

由式(1)可知，向下延拓可以表示為

T(x,z)=T(x,0)*φ(x)

(2)

即觀測(cè)平面上的位場(chǎng)值與一個(gè)核函數(shù)的褶積，*代表褶積運(yùn)算，由富氏變換的性質(zhì)可知空間域的褶積與波數(shù)域的乘積相對(duì)應(yīng)。對(duì)(2)式兩端做富氏變換可得

(3)

根據(jù)劉天佑[21]和曾小牛等[19]的研究，可以得到迭代Tikhonov正則化方法的波數(shù)域算子的推導(dǎo)過程為

gn(x)=atikT(x)+btikgn-1(x)

(4)

其中g(shù)n(x)為迭代n次結(jié)果；T(x)為原始數(shù)據(jù)；atik和btik為延拓系數(shù)。

對(duì)(1)式兩端做富氏變換可得

Φ(w)Gn-1(w)]

(5)

由數(shù)學(xué)歸納法可得

BtikGn-2(w)]

(6)

由此可以得到迭代Tikhonov正則化法的波數(shù)域算子為

(7)

3 正演模型數(shù)據(jù)延拓

本文借鑒李才明、管志寧等表述的球體正演公式以及有限延深厚板的正演公式進(jìn)行正演計(jì)算。

(2x2-R2)cosiscosIcosA′-(x2+R2)·

siniscosIcot-1Isin2A′-6RxsiniscosIcosA′]

(8)

(9)

式(1)為球體磁異常正演公式，式(2)為有限延深厚板磁異常正演公式。式(1)中ms為有效磁矩(A·m)，x為水平坐標(biāo)(m)，R球體中心埋深(m)，is、I、A′分別為有效磁化傾角、地磁傾角和測(cè)線方向與地磁在水平面投影的夾角(°)。式(2)中Ms為有效磁化強(qiáng)度(A/m)，u0真空磁導(dǎo)系數(shù)4π×10-7(H/m)，is、I同式(1)，α為板傾角(°)，rA、rB、rC、rD分別為板四個(gè)角點(diǎn)到測(cè)點(diǎn)的距離(m)，φA、φB、φC、φD分別為板四個(gè)角點(diǎn)與測(cè)點(diǎn)法線的夾角(°)。

根據(jù)曾小牛等[19]可知，Tikhonov迭代正則化向下延拓主要通過正則化因子和迭代次數(shù)實(shí)現(xiàn)高精度、快速的向下延拓，正則化因子和迭代次數(shù)的選取對(duì)于延拓的精度和速度有著至關(guān)重要的影響，通過對(duì)Tikhonov正則化向下延拓波數(shù)域算子的研究，可以為正則化因子的選取以及迭代次數(shù)的選取提供指導(dǎo)，總體來說正則化因子取0～1之間的數(shù)值，迭代次數(shù)取1～100較為合適。對(duì)于一般情況，淺部磁性介質(zhì)可以選取較大正則化因子和較少的迭代次數(shù)，而深部磁性介質(zhì)則需要較小的正則化因子和較多的迭代次數(shù)。

圖2和圖3分別為球體異常剖面圖、球體異常下延拓等值線圖以及球體正演模型圖(三個(gè)球體分別位于：100 m、500 m、900 m，埋深均為100 m，垂直磁化Ms=5 000A·m2，dh=10 m，正則化因子α=0.1，迭代次數(shù)n=5，點(diǎn)距x=10 m)?？梢钥闯?，通過Tikhonov迭代正則化向下延拓后，隨著延拓深度的增加，磁異常并沒有出現(xiàn)明顯的畸變。磁異常明顯收斂為三個(gè)正異常區(qū)和兩個(gè)負(fù)異常區(qū)，而且可以發(fā)現(xiàn)正異常區(qū)兩側(cè)均存在一個(gè)負(fù)異常區(qū)，這與球體的空間等值線圖相似，也就是說通過Tikhonov迭代正則化法向下延拓可以近似的得到磁性地質(zhì)體的近似空間等值線圖，通過比較正演模型圖和延拓成果圖，可以清楚地發(fā)現(xiàn)三個(gè)正異常所對(duì)應(yīng)的位置恰恰是正演模型中磁性球狀地質(zhì)體所對(duì)應(yīng)的位置。從而通過Tikhonov迭代正則化向下延拓方法可以近似獲得地下磁性介質(zhì)的位場(chǎng)空間分布情況，進(jìn)而可以較為清楚地反映磁性地質(zhì)體的埋深。

圖2 三個(gè)同深度球體正演磁異常剖面Fig.2 Forward magnetic abnormal profile of the three similar deep spheres

圖3 三個(gè)同深度球體向下延拓結(jié)果以及正演模型Fig.3 Forward model and the downward continuation result of the three similar deep spheres

圖4和圖5分別為球體異常剖面圖、球體異常下延拓等值線圖以及球體正演模型圖(三個(gè)球體分別位于：300 m、500 m、800 m，埋深為50 m、150 m、250 m，垂直磁化M1=20 000 A·m2，M2=2×105A·m2，M1=2×106A·m2，dh=10 m，正則化因子α=0.05，迭代次數(shù)n=10，點(diǎn)距x=10 m)。可以看出，Tikhonov迭代正則化向下延拓后，磁異常同樣沒有出現(xiàn)明顯的畸變現(xiàn)象，延拓后磁異常收斂為三個(gè)正異常區(qū)域和三個(gè)負(fù)異常區(qū)域，而且可以較清楚地看出，正異常區(qū)域兩側(cè)伴隨著負(fù)異常，形成類似雙扭線的形狀。另外可以看出，正異常在收斂的時(shí)候，淺部的異常收斂較為集中，而且異常值較大；隨著深度的增加，深部的正異常收斂相對(duì)較緩，而且異常值相對(duì)較小。盡管由于深度不同而造成的收斂程度出現(xiàn)不同，但是通過與正演模型的對(duì)比，依然可以通過Tikhonov迭代正則化向下延拓得到的等值線圖中的正異常收斂區(qū)判斷出磁性地質(zhì)體的埋深。但是可以注意到當(dāng)深部異常體埋深正確時(shí)，淺部的異常體埋深往往會(huì)偏大。

圖4 三個(gè)不同深度球體正演異常剖面Fig.4 Forward magnetic abnormal profile of the three different deep spheres

圖5 不同深度球體異常向下延拓結(jié)果圖以及正演模型Fig.5 Forward model and the downward continuation result of the three different deep spheres

圖6和圖7分別為球體異常剖面圖、球體異常下延拓等值線圖以及球體正演模型圖(兩個(gè)球體均位于：500 m處，埋深分別為100 m、500 m，垂直磁化M1=2×104A·m2、M2=2×105A·m2，dh=10 m，正則化因子α=0.5，迭代次數(shù)n=175，點(diǎn)距x=10 m)?？梢钥闯?，通過Tikhonov迭代正則化向下延拓后，磁異常并沒有出現(xiàn)明顯的畸變現(xiàn)象，而且延拓后形成兩個(gè)正異常收斂區(qū)和四個(gè)負(fù)異常收斂區(qū)，可以注意到在正異常收斂區(qū)的兩側(cè)都伴隨有負(fù)異常的收斂區(qū)。觀察正異常收斂區(qū)，可以發(fā)現(xiàn)淺部的正異常收斂區(qū)較深部正異常收斂區(qū)收斂較為集中，而且正異常值較大。通過與正演模型對(duì)比，Tikhonov迭代正則化向下延拓依然可以較為清楚地反映磁性地質(zhì)體的埋深，而且可以很好地區(qū)分疊加磁性介質(zhì)。

圖6 疊加球體正演數(shù)據(jù)Fig.6 Forward data of overlay spheres

圖7 疊加球體向下延拓成果圖以及正演模型Fig.7 Forward model and downward continuation result of overlay spheres

圖8和圖9分別為有限延深厚板ΔT、Ha、Za異常剖面圖、有限延深厚板ΔT、Ha、Za異常下延拓等值線圖以及有限延深厚板正演模型圖(厚板頂界面位于：100 m處，寬度為100 m，傾角45°，順層磁化Ms=10 A·m2，dh=10 m，正則化因子α=0.1，迭代次數(shù)n=10，點(diǎn)距x=10 m)。可以看出，Tikhnovo迭代正則化向下延拓對(duì)ΔT磁異常進(jìn)行處理后，磁異常沒有出現(xiàn)明顯的畸變現(xiàn)象，磁異常逐漸收斂為一個(gè)較為明顯的正異常和一個(gè)較為明顯的負(fù)異常，但是在正異常的左側(cè)存在一個(gè)微弱的中低異常收斂區(qū)，在負(fù)異常的右側(cè)存在一個(gè)中高異常的收斂區(qū)，通過與正演模型相比較，雖然在一定程度上可以反映磁性介質(zhì)的諸如埋深寬度等物性參數(shù)，但是效果并不明顯。Tikhonov迭代正則化向下延拓對(duì)Za磁異常處理后，磁異常沒有出現(xiàn)明顯的畸變現(xiàn)象，磁異常逐漸收斂為一個(gè)正異常區(qū)域和兩個(gè)負(fù)異常區(qū)域，兩個(gè)負(fù)異常區(qū)域伴隨于正異常區(qū)域的兩側(cè)。通過與正演模型比較可以發(fā)現(xiàn)，正異常收斂區(qū)對(duì)應(yīng)的位置與厚板的頂界面埋深相近，而且正異常區(qū)域的寬度也與厚板的寬度相近，正異常區(qū)的傾向與正演模型的磁化方向也存在近似鏡像的關(guān)系。Tikhonov迭代正則化向下延拓處理Ha磁異常，可以發(fā)現(xiàn)磁異常沒有出現(xiàn)明顯的畸變現(xiàn)象，而且磁異常逐漸收斂為兩個(gè)正異常區(qū)和一個(gè)負(fù)異常區(qū)，負(fù)異常區(qū)兩側(cè)伴隨存在正異常區(qū)。與正演模型比較，可以看出雖然負(fù)異常區(qū)域的傾向與厚板的磁化傾角方向一致，但是對(duì)于磁性介質(zhì)的埋深和寬度并不能得到很好的反映。綜合比較ΔT、Ha、Za延拓結(jié)果，可以清楚地發(fā)現(xiàn)Za向下延拓等值線圖能夠較為明顯地反映地下磁性介質(zhì)的埋深、寬度、磁化方向等參數(shù)。

圖8 有限延深厚板ΔT、Ha、Za異常剖面Fig.8 The ΔT、Ha、Za exceptions profile of the finite deep thick plate

圖9 有限延深厚板ΔT、Ha、Za異常延拓圖以及正演模型Fig.9 Forward model and downward continuation results of the finite deep thick plate

4 實(shí)際資料數(shù)據(jù)處理

現(xiàn)以某工區(qū)的高精度實(shí)測(cè)資料為例，運(yùn)用Tikhonov迭代正則化向下延拓法對(duì)其進(jìn)行處理。

圖10和圖11分別為某工區(qū)內(nèi)呈NW-SE向的高精度實(shí)測(cè)磁異常剖面圖和對(duì)實(shí)測(cè)異常進(jìn)行Tikhonov正則化向下延拓得到的成果圖。由于測(cè)區(qū)地勢(shì)存在一定程度的起伏，所以對(duì)原始磁測(cè)資料進(jìn)行曲化平等處理后得到如圖10所示的磁異常剖面圖。通過圖11，可以看出在淺部存在一些小的正異常區(qū)， 在深度200m水平坐標(biāo)500m處存在一處明顯的正異常，而且在正異常兩側(cè)存在明顯的兩個(gè)伴隨負(fù)異常，同樣在水平坐標(biāo)1 500 m、深度550 m處和水平坐標(biāo)2 650 m、深度500 m處存在兩個(gè)較為明顯的正異常區(qū)域，其兩側(cè)也存在較為明顯的伴隨負(fù)異常。結(jié)合圖4和圖5推斷這三處存在三個(gè)磁性介質(zhì),深度分別對(duì)應(yīng)為200 m、550 m和500 m，后經(jīng)過鉆探小組鉆孔驗(yàn)證，其位于水平坐標(biāo)500 m、深度184.5 m處可見鈦磁鐵礦，與向下延拓推斷相近。

圖10 某工區(qū)高精度磁異常剖面Fig.10 A work-area high-precision abnormal profile

圖11 某工區(qū)高精度磁異常剖面向下延拓成果Fig.11 Downward continuation result of a work-area high-precision abnormal profile

5 結(jié) 論

通過上述討論可以看出，Tikhonov迭代正則化向下延拓法結(jié)合了迭代法和波數(shù)域延拓法的優(yōu)點(diǎn)，在提高向下延拓精度的同時(shí)也提高了延拓的速度。該方法具有向下延拓距離大、精度高、過場(chǎng)源不奇異等優(yōu)點(diǎn)，但也存在延拓精度過多依賴于迭代次數(shù)和正則化因子的選取等缺點(diǎn)。通過以上簡(jiǎn)述與比較，可以得出如下結(jié)論：

1)Tikhonov迭代正則化向下延拓處理后，會(huì)出現(xiàn)多個(gè)正異常和負(fù)異常的收斂區(qū)，通過研究可以看出，磁性體對(duì)應(yīng)的正異常收斂區(qū)兩側(cè)一般會(huì)伴隨出現(xiàn)負(fù)異常區(qū)，這種特征可以作為實(shí)際應(yīng)用中推測(cè)異常體的依據(jù)。

2)正負(fù)異常收斂區(qū)呈近似雙扭線狀，這與球體磁位場(chǎng)的空間等值線圖相似，正異常收斂區(qū)往往對(duì)應(yīng)地下磁性介質(zhì)。

3)不同深度的異常體引起的磁異常在經(jīng)過Tikhonov迭代正則化向下延拓處理后，其對(duì)應(yīng)的正異常區(qū)域收斂程度會(huì)存在差異，而且對(duì)于同一剖面中不同深度的磁性體在延拓等值線圖中反映的深度會(huì)存在一定差異。

4)Tikhonov迭代正則化向下延拓法能夠很好地區(qū)分疊加異常。

5)垂直磁化的磁性體或者對(duì)Za進(jìn)行延拓后得到的成果圖能夠更加直觀地反映地下磁性介質(zhì)的物性參數(shù)，因此在進(jìn)行向下延拓時(shí)可以通過化極或者參數(shù)變換等方式來提高向下延拓的成果。

[1]李才明.重磁勘探原理與方法[M].北京：科學(xué)出版社,2013.

[2]許大欣.利用重力異常匹配技術(shù)實(shí)現(xiàn)潛艇導(dǎo)航[J].地球物理學(xué)報(bào),2005,48(4):812-816.

[3]郭才發(fā),胡正東,張士峰,等.地磁導(dǎo)航綜述[J].宇航學(xué)報(bào),2009, 30(4):1 314-1 319.

[4]楊云濤,石志勇,關(guān)貞珍,等.地磁場(chǎng)在導(dǎo)航定位系統(tǒng)中的應(yīng)用[J].中國(guó)慣性技術(shù)學(xué)報(bào),2007,15(6):686-692.

[5]彭富清.地磁模型與地磁導(dǎo)航[J].海洋測(cè)繪,2006,26(2):73-75.

[6]雷文敏,吳健生.基于Canny算子深源磁異常圖像的邊緣檢測(cè)[J].工程地球物理學(xué)報(bào),2013,10(2):127-131.

[7]王邦華,王理.重磁位場(chǎng)的正則化向下延拓[J].物探化探計(jì)算技術(shù),1998(1):31-36.

[8]管志寧.地磁場(chǎng)與磁力勘探[M].北京：地質(zhì)出版社，2005.

[9]徐世浙.位場(chǎng)延拓的積分—迭代法[J].地球物理學(xué)報(bào),2006,49(4):1 176-1 182.

[10]徐世浙.迭代法與FFT法位場(chǎng)向下延拓效果的比較[J].地球物理學(xué)報(bào),2007,50(1):285-289.

[11]陳龍偉,張輝,鄭志強(qiáng),等.水下地磁輔助導(dǎo)航中地磁場(chǎng)延拓方法[J].中國(guó)慣性技術(shù)學(xué)報(bào),2007,15(6):693-697.

[12]Clarke G K C. Optimum Second-Derivative and Downward-Continuation Filters[J]. Geophysics, 1969, 34(3):424-437.[13]Ku C C. The Use of Linear Filtering in Gravity Problems[J]. Geophysics, 2012,36(36).

[14]安玉林,管志寧.濾除高頻干擾的正則化穩(wěn)定因子[J].物化探計(jì)算技術(shù),1985(1):13-23.

[15]梁錦文. 位場(chǎng)向下延拓的正則化方法[J].地球物理學(xué)報(bào),1989,32(5):600-608.

[16]王萬銀，劉展，肖慶疆.位場(chǎng)向下延拓的自適應(yīng)濾波方法[C]//周鳳桐.勘查地球物理勘查地球化學(xué)文集.北京：地質(zhì)出版社,2004:30-35.

[17]宋滔,羅勇,王宇航.有限單元法在二維位場(chǎng)延拓中的應(yīng)用[J].工程地球物理學(xué)報(bào),2012,9(2):134-138.

[18]劉東甲,洪天求,賈志海,等.位場(chǎng)向下延拓的波數(shù)域迭代法及其收斂性[J].地球物理學(xué)報(bào),2009,52(6):1 599-1 605.

[19]曾小牛,李夕海,韓紹卿,等.位場(chǎng)向下延拓三種迭代方法之比較[J].地球物理學(xué)進(jìn)展,2011,26(3):908-915.

[20]傅初黎,李洪芳,熊向團(tuán).不適定問題的迭代Tikhonov正則化方法[J].計(jì)算數(shù)學(xué),2006,28(3):237-246.

[21]劉天佑.磁法勘探[M].北京：地質(zhì)出版社，2013.

The Application of Tikhonov Iteration Regularization Downward Continuation

Shen Hengxu,Li Caiming

The potential field downward continuation of gravity magnetic is applied to the gravity magnetic data interpretation in the geo-exploration, energy exploration and military navigation. The method of Tikhonov iteration regularization downward continuation realizes speedy and high-accuracy downward continuation. The method can suppress high-frequency shock through source body without strange. However, the respond degree of the method for the geologic body and the reflect of the method for the physical parameters need to be studied. The paper based on the Tikhonov iteration regularization downward continuation method discusses the downward continuation result for different geological bodies. By dealing with the magnetic data of different geological factors simulated by several balls of different deeps and plate bodies, we discuss the response result of the downward continuation result for the geological body. Some tests for the theoretical field data and observed field data indicate that the new normalized method will lead to a real representation of potential field in the whole half-space under the observable surface and can be used in determining the vertical distributing of source bodies.

potential field; Tikhonov iteration regularization; downward continuation; source body

1672—7940(2016)05—0621—06

10.3969/j.issn.1672-7940.2016.05.011

沈恒旭(1980-)，男，碩士研究生，主要從事礦產(chǎn)勘查與工程勘測(cè)方向的學(xué)習(xí)與研究。E-mail:1218886808@qq.com

P631

A

2016-05-10