應(yīng)用迭代收縮高分辨率Radon變換的繞射波分離與成像方法

羅騰騰 徐基祥 孫夕平

(中國(guó)石油勘探開(kāi)發(fā)研究院，北京 100083)

0 引言

地層中的斷層、裂隙、尖滅和孔洞等小尺度地質(zhì)體對(duì)油氣勘探具有重要意義，這些小尺度地質(zhì)體在地震記錄上的響應(yīng)主要為能量較弱的繞射波。Klem等[1]指出地震數(shù)據(jù)中的繞射能量通常比反射能量弱1～2個(gè)數(shù)量級(jí)。在常規(guī)的地震數(shù)據(jù)處理流程中，來(lái)自反射界面以外的繞射波同相軸常被視作干擾噪聲而被濾除。即使同時(shí)將繞射波和反射波準(zhǔn)確偏移歸位，弱繞射能量仍會(huì)被大尺度地質(zhì)體反射能量淹沒(méi)。為了充分利用地震數(shù)據(jù)中的繞射信息，可從全波場(chǎng)中提取繞射波，實(shí)現(xiàn)繞射波單獨(dú)成像，進(jìn)而提高繞射目標(biāo)的成像分辨率。

早在二十世紀(jì)五十年代，Krey[2]就利用繞射波信息研究斷層、斷點(diǎn)等，但直到二十世紀(jì)八十年代，學(xué)者們才真正意識(shí)到繞射波分離成像的重要性，并相繼提出不同的繞射波場(chǎng)分離方法。根據(jù)地震數(shù)據(jù)處理階段的不同，可以將繞射波分離與成像方法歸納為以下三類(lèi)：

(1)基于疊前數(shù)據(jù)的繞射波分離與成像。該類(lèi)方法是利用疊前道集中繞射波與反射波的運(yùn)動(dòng)學(xué)特征差異，從全波場(chǎng)中提取繞射波，進(jìn)而實(shí)現(xiàn)繞射波單獨(dú)成像。Landa等[3]提出在共炮檢距道集上利用相位校正和定性分析檢測(cè)小尺度繞射體。Khaidu-kov等[4]在疊前共炮點(diǎn)道集中利用反射波與繞射波的運(yùn)動(dòng)學(xué)特征差異，即反射波聚焦于震源點(diǎn)的鏡像位置而繞射波仍保持發(fā)散的狀態(tài)，通過(guò)反射波聚焦、切除、再反聚焦方法成功分離出繞射波場(chǎng)。Nowak等[5]根據(jù)共炮點(diǎn)道集中反射波與繞射波同相軸橫向位置差異，借助于加權(quán)Radon變換分離反射繞射波場(chǎng)。Bansal等[6]討論了多種利用疊前數(shù)據(jù)壓制反射信號(hào)、增強(qiáng)繞射能量的方法，主要包括共炮檢距傾角濾波、正常時(shí)差校正傾角濾波、特征向量濾波、Radon濾波等，并且對(duì)比、分析了不同方法的優(yōu)缺點(diǎn)，認(rèn)為特征向量法是較為有效的繞射波場(chǎng)分離方法。Taner等[7]利用平面波分解(Plan Wave Decomposition，PWD)濾波器生成平面波記錄，依據(jù)繞射波與反射波在平面波記錄中的同相軸差異，在疊前道集中分離繞射波場(chǎng)。黃建平等[8]對(duì)PWD解構(gòu)濾波器的使用方法和在疊前、疊后道集上分離成像繞射波的流程進(jìn)行了系統(tǒng)總結(jié)?？籽┑萚9]和劉玉金等[10]證實(shí)了平面波記錄具有區(qū)分繞射波與反射波的特性，并利用平面波解構(gòu)濾波技術(shù)分離獲得繞射波場(chǎng)。朱生旺等[11]解決了PWD濾波器中繞射波場(chǎng)的低傾角信息估計(jì)不足問(wèn)題，提高了繞射波成像結(jié)果的橫向分辨率。

(2)在疊加過(guò)程中進(jìn)行繞射波分離與成像。Kanasewich等[12]直接對(duì)地震信號(hào)沿著繞射波走時(shí)曲線進(jìn)行疊加，獲得了共炮檢距道集和共斷層點(diǎn)道集上的繞射波時(shí)間剖面。Berkovich等[13]、Dell等[14]分別利用繞射波多聚焦和共反射面的繞射走時(shí)計(jì)算方法，將選取的繞射波同相軸進(jìn)行最優(yōu)疊加，得到包含繞射波的疊加剖面。Asgedom等[15]利用修正的共反射面(Common Reflection Surface，CRS)技術(shù)討論了相似度和多信號(hào)分類(lèi)兩種方法在分離反射波與繞射波時(shí)的優(yōu)勢(shì)。

(3)在偏移過(guò)程中進(jìn)行繞射波分離與成像。Landa等[16]推導(dǎo)了深度偏移傾角域共成像點(diǎn)道集(Common Imaging Gathers，CIGs)中反射波和繞射波的形態(tài)，并分析了偏移速度對(duì)兩者的影響。Zhu等[17]在小波束偏移過(guò)程中，利用角度域局部成像矩陣中反射與繞射能量的分布差異，即反射能量沿著傾角方向呈線性分布而繞射能量彌散分布在整個(gè)矩陣中，分離出繞射能量并對(duì)其單獨(dú)成像。Klokov等[18]將相似度掃描反射頂點(diǎn)去除與混合Radon變換方法相結(jié)合，在傾角域CIGs上進(jìn)行繞射波場(chǎng)分離。Decker 等[19]在偏移傾角道集中，對(duì)每個(gè)固定傾角的成像剖面應(yīng)用PWD解構(gòu)濾波技術(shù)，獲得單獨(dú)的繞射波成像剖面。Zhang等[20]提出在單炮的反射角域CIGs中通過(guò)切除菲涅耳帶提取繞射波場(chǎng)。黃建平等[21]和李振春等[22]將逆時(shí)偏移與最小二乘偏移的優(yōu)勢(shì)相結(jié)合，分別發(fā)展了適用于復(fù)雜構(gòu)造的保幅成像處理方法。劉夢(mèng)麗等[23]提出一種滿足能量守恒的逆散射成像條件，并將其與最小二乘逆時(shí)偏移(Reverse Time Migration，RTM)相結(jié)合，實(shí)現(xiàn)了角度濾波成像。汪天池等[24]在傾角域中利用逆時(shí)偏移提取傾角道集，然后采用中值濾波方法實(shí)現(xiàn)反射波與繞射波的分離。

近些年來(lái)，大量學(xué)者嘗試在傾角域CIGs上進(jìn)行繞射波與反射波的分離。根據(jù)傾角域CIGs中反射波響應(yīng)曲線呈雙曲規(guī)律而繞射波呈現(xiàn)擬線性規(guī)律的差異，發(fā)展了平面波解構(gòu)濾波、掃描相似度頂點(diǎn)去除與混合Radon變換相結(jié)合[18]、相似譜能量分析等繞射波分離方法。同時(shí)，由于在傾角域CIGs中繞射波同相軸對(duì)速度誤差的敏感性，Landa等[16]、Reshef等[25]、Li等[26]將其用于偏移速度分析以提高速度提取精度。

鑒于傳統(tǒng)混合域Radon變換采用迭代重加權(quán)算法進(jìn)行求解，每一步迭代常常涉及到大規(guī)模矩陣的求逆問(wèn)題，計(jì)算量巨大。為此，本文引入迭代收縮閾值算法(Iterative Shrinkage Thresholding Algorithm，ISTA)，發(fā)展了一種基于迭代收縮高分辨率Radon變換的繞射波分離與成像方法。首先采用迭代收縮閾值算法，通過(guò)時(shí)間域簡(jiǎn)單的模型收縮步驟獲得Radon域模型的稀疏性，以解決常規(guī)Radon變換中反射波與繞射波在Radon域中能量團(tuán)聚焦性差、繞射波場(chǎng)分離效果不理想的問(wèn)題； 然后在編程實(shí)現(xiàn)算法的基礎(chǔ)上，通過(guò)模型數(shù)據(jù)驗(yàn)證了該方法的有效性和適用性。

1 方法原理

1.1 傾角域共成像點(diǎn)道集波場(chǎng)特征

常速介質(zhì)疊后偏移(圖1)反射面的深度z(x)可以表示為[17]

z(x)=z1+xtanα0

(1)

式中：z1為傾斜反射界面與z軸的交點(diǎn)深度；α0為反射界面的傾角；x為橫坐標(biāo)。在地表y處的自激自收地震響應(yīng)可表示為

(2)

式中：v是介質(zhì)速度；w為旅行時(shí)間。

由模型坐標(biāo)(x，z)和數(shù)據(jù)坐標(biāo)(y，w)之間的映射關(guān)系可得偏移公式，即

(3)

τ=wcosα

(4)

圖1 零炮檢距反射示意圖

式中：vm表示偏移速度；α為偏移傾角；τ表示截距時(shí)間。將式(2)代入式(3)和式(4)中，消去y可得傾斜反射界面在傾角域的成像表達(dá)式，即

(5)

式中γ=vm/v。根據(jù)式(1)所定義的橫坐標(biāo)xm處的反射面深度為z0，可得

xmsinα0+z1cosα0=z0cosα0

(6)

將式(6)代入式(5)，可得傾角域CIGs中固定橫向位置的反射波響應(yīng)表達(dá)式，即

(7)

式中τ0=2z0/v為反射點(diǎn)處的雙程垂直旅行時(shí)。

將式(7)對(duì)α求偏導(dǎo)可得到

(8)

從式(8)可以看出，當(dāng)偏移速度正確，即γ=1時(shí)，反射響應(yīng)為具有穩(wěn)相頂點(diǎn)的“笑臉”狀曲線，且穩(wěn)相頂點(diǎn)位于α=α0處。

類(lèi)似地，可以得到地下繞射點(diǎn)的傾角域響應(yīng)表達(dá)式

(9)

式中ξ=(x-x0)/z0，ξ是繞射點(diǎn)x0和成像點(diǎn)x之間的橫向距離與繞射點(diǎn)深度z0的比值。當(dāng)偏移速度準(zhǔn)確(vm=v)且成像點(diǎn)在繞射點(diǎn)處(ξ=0)，式(9)簡(jiǎn)化為

τ(α)=τ0

(10)

由式(10)可以看出，此時(shí)繞射響應(yīng)為一水平直線。

對(duì)于式(7)和式(9)，利用τ(α)=2z(α)/vm和τ0=2z0/v，即可從時(shí)間域轉(zhuǎn)換到偽深度域。深度域反射界面在傾角域的響應(yīng)表達(dá)式為

(11)

地下繞射點(diǎn)在傾角域的響應(yīng)表達(dá)式為

(12)

式中ρ=x-x0。

通過(guò)一個(gè)簡(jiǎn)單的常速介質(zhì)模型(圖2)可直觀地闡明傾角域CIGs響應(yīng)曲線的特點(diǎn)，并且可以進(jìn)一步分析兩者之間的差異。該模型由一個(gè)垂直斷層和一個(gè)繞射體組成，記錄共150炮，最大炮檢距為500m，時(shí)間采樣間隔為1ms，檢波點(diǎn)間距為10m。

圖2 常速介質(zhì)模型垂直斷層和繞射體分別位于0.8km和1.2km

圖3和圖4分別顯示了0.8km處斷層和1.2km處繞射體在不同偏移速度下產(chǎn)生的傾角道集?？梢钥闯?，在傾角域CIGs上，位于繞射點(diǎn)正上方的繞射響應(yīng)可分為三種情況：當(dāng)偏移速度小于準(zhǔn)確的介質(zhì)速度時(shí)，繞射波表現(xiàn)為向上彎曲的“笑臉”狀(凹形)(圖3a、圖4a)；當(dāng)偏移速度準(zhǔn)確時(shí)，繞射波同相軸表現(xiàn)為線性形態(tài)(圖3b、圖4b)；當(dāng)偏移速度大于準(zhǔn)確的介質(zhì)速度時(shí)，繞射波同相軸則為向下彎曲的“哭臉”狀(凸形)(圖3c、圖4c)。而反射界面響應(yīng)曲線無(wú)論偏移速度準(zhǔn)確與否，始終保持為具有穩(wěn)相頂點(diǎn)的“笑臉”狀，且穩(wěn)相頂點(diǎn)橫坐標(biāo)指示界面的傾角。

綜上所述，在傾角域CIGs中，反射響應(yīng)與繞射響應(yīng)曲線之間存在著明顯的幾何形態(tài)差異，因而可以利用兩者之間形態(tài)差異作為波場(chǎng)分離的依據(jù)。

圖3 不同偏移速度誤差時(shí)斷層處的傾角域CIGs(a)偏移速度為介質(zhì)速度的90%； (b)偏移速度等于介質(zhì)速度； (c)偏移速度為介質(zhì)速度的110%藍(lán)色和紅色箭頭分別指示反射波和繞射波同相軸，圖4同

圖4 不同偏移速度誤差時(shí)繞射體處的傾角域CIGs(a)偏移速度為介質(zhì)速度的90%； (b)偏移速度等于介質(zhì)速度； (c)偏移速度為介質(zhì)速度的110%

1.2 常規(guī)Radon變換

二維拋物線Radon變換的定義式可以表示為

(13)

式中：d(t，u)為時(shí)間—空間域地震數(shù)據(jù)，其中t為時(shí)間，u為炮檢距；m(τ，q)為Radon域數(shù)據(jù)，其中q表示曲率參數(shù)。

因此，拋物線Radon正變換的離散形式可以表示為

(14)

式中：D(ω，u)和M(ω，q)分別為d(t，u)和m(τ，q)的傅里葉變換結(jié)果，其中ω為頻率；n為Radon域的曲率個(gè)數(shù)。

由式(14)可知，Radon變換可以表示為一個(gè)線性方程組的求解問(wèn)題，即

d=Lm

(15)

式中：d為地震數(shù)據(jù)；m為模型數(shù)據(jù)；L是Radon變換算子，由地震數(shù)據(jù)的采集參數(shù)和Radon變換參數(shù)共同確定。

求式(15)的反問(wèn)題，不能完全滿足存在性、唯一性和穩(wěn)定性，故該反問(wèn)題是不適定的。對(duì)于不適定問(wèn)題，為獲得唯一且穩(wěn)定的解，最直接方法就是在反問(wèn)題中引入與模型信息相關(guān)的正則化項(xiàng)以求解稀疏正則解，即高分辨率Radon模型解。

對(duì)于式(15)，可構(gòu)造如下目標(biāo)函數(shù)，即

(16)

(17)

對(duì)式(17)求極小值，可以得到m的常規(guī)最小二乘解，即

m=(LTL+βI)-1LTd

(18)

式中I為單位矩陣。由式(18)可以看出，在傳統(tǒng)最小二乘意義下求得的Radon模型解分辨率較低。由于對(duì)每一個(gè)頻率所采用的阻尼因子相同，且阻尼因子的實(shí)際作用是對(duì)Radon模型數(shù)據(jù)進(jìn)行平滑，導(dǎo)致模型空間的能量團(tuán)無(wú)法較好地聚焦，降低了Radon變換的分辨率。因此，常規(guī)最小二乘Radon變換不能滿足繞射波場(chǎng)分離的高分辨率需求。

1.3 迭代收縮高分辨率Radon變換

為了獲得稀疏的高分辨率Radon模型解，一般為正則項(xiàng)選用L1范數(shù)、數(shù)據(jù)誤差項(xiàng)選用L2范數(shù)，即p=1、s=2，則式(18)可轉(zhuǎn)換為

(19)

在求解稀疏反演問(wèn)題(式(19))時(shí)，本文利用Lu[27]提出的混合頻率—時(shí)間域迭代收縮閾值算法，即在時(shí)間域通過(guò)模型收縮步驟獲得Radon模型的稀疏性，矩陣與向量的乘積操作在頻率域?qū)崿F(xiàn)。具體來(lái)說(shuō)，在第k步迭代過(guò)程中，更新的模型解為

mk=Tσ{mk-1+2t0F-1{(LTL)-1LT[F(d)-

LF(mk-1)]}

(20)

式中：k為當(dāng)前迭代次數(shù)；t0是迭代步長(zhǎng)；F、F-1分別表示正、反傅里葉變換；Tσ是一個(gè)收縮算子，定義為

(21)

(22)

本文基于迭代收縮算法的高分率Radon變換處理步驟如下。

(1)給定閾值σ、最大迭代次數(shù)K和二維中值濾波器的大小以及初始迭代步長(zhǎng)t1，將當(dāng)前迭代次數(shù)k設(shè)置為0。

(2)選取Radon模型的最小二乘解作為模型迭代的初始值，m0=F-1[(LTL)-1LTd]，同時(shí)設(shè)置k=1。

(3)在迭代變量k=1，2，...，K時(shí)，執(zhí)行下述迭代過(guò)程：

①將上次的模型值mk-1傅里葉變換到頻率域得到F(mk-1)；

②計(jì)算Radon算子L與頻率域模型值的乘積即LF(mk-1)；

③計(jì)算頻率域中數(shù)據(jù)擬合誤差F(d)-LF(mk-1)；

④計(jì)算算子的偽逆矩陣(LTL)-1LT與數(shù)據(jù)擬合誤差的乘積，即(LTL)-1LT[F(d)-LF(mk-1)]；

⑤利用傅里葉反變換將頻率域數(shù)據(jù)變換到時(shí)間域，即F-1{(LTL)-1LT[F(d)-LF(mk-1)]}；

⑥利用帶有收縮運(yùn)算的方程Tσ計(jì)算模型更新值mk。

(4)循環(huán)增加迭代次數(shù)k，檢查是否達(dá)到了最大迭代次數(shù)，如果“否”，則轉(zhuǎn)到第(3)步。

(5)輸出Radon模型的高分辨率解。

利用上述迭代收縮高分辨Radon變換算法，在求解過(guò)程中僅需計(jì)算一次偽逆矩陣，相比于傳統(tǒng)的迭代重加權(quán)最小二乘(Iterative Reweighted Least Squares，IRLS)求解算法，極大地降低了計(jì)算量，更適用于實(shí)際地震資料處理。

2 模型試算

2.1 有效性

下面用模型試算驗(yàn)證本文算法的有效性。模擬CIGs中包含三條反射波和兩條繞射波同相軸(圖5)，反射波頂點(diǎn)分別位于深度600、1000、1500m處，兩條繞射波同相軸則分別處于深度1300、2000m的位置。

圖6為模型數(shù)據(jù)在不同迭代次數(shù)下的試算結(jié)果。測(cè)試過(guò)程中，設(shè)置最大迭代次數(shù)K=20、步長(zhǎng)t0=0.5、閾值算子σ=0.05，二維均值濾波器的大小與模型相同。

圖5 模擬數(shù)據(jù)的CIG記錄

從圖6可以看出，隨著迭代次數(shù)的增加，Radon域模型的分辨率逐漸提高，即反射波和繞射波能量團(tuán)聚焦性逐漸變好，剪刀狀拖尾能量逐漸變?nèi)?，分離后的繞射波場(chǎng)也越來(lái)越干凈。在此模擬數(shù)據(jù)的測(cè)試過(guò)程中，當(dāng)?shù)螖?shù)K=20時(shí)，Radon域中的反射波和繞射波能量聚焦性很好，分辨率高，切除反射后的繞射波場(chǎng)中幾乎不存在殘余反射能量。

2.2 抗噪能力

為了測(cè)試本文方法的抗噪能力，在圖5的模擬數(shù)據(jù)中加入了高斯白噪，得到如圖7所示數(shù)據(jù)。圖8左為圖7的常規(guī)最小二乘Radon變換結(jié)果，包含大量隨機(jī)噪聲，圖8右為分離后得到的時(shí)空域繞射波波場(chǎng)，由于在Radon域?qū)υ肼暤膲褐撇蛔?，?yán)重影響分離后的繞射波場(chǎng)質(zhì)量。圖9為迭代收縮高分辨率Radon變換結(jié)果?？梢钥闯觯旱湛s高分辨率Radon變換方法能夠壓制原始數(shù)據(jù)中的隨機(jī)噪聲，有效提高原始數(shù)據(jù)的信噪比(圖9左)；本文方法從低信噪比的數(shù)據(jù)中提取繞射波的能力較強(qiáng)，可很好地壓制模擬數(shù)據(jù)中的隨機(jī)噪聲(圖9右)。

圖6 不同迭代次數(shù)下迭代收縮高分辨率Radon變換結(jié)果

圖8 常規(guī)最小二乘Radon變換結(jié)果左為Radon域全波場(chǎng)，右為時(shí)空域繞射波波場(chǎng)

圖9 迭代收縮高分辨率Radon變換結(jié)果左為Radon域全波場(chǎng)，右為時(shí)空域繞射波波場(chǎng)

2.3 計(jì)算效率

為了說(shuō)明本文方法的計(jì)算效率相較于IRLS算法的優(yōu)勢(shì)，不同迭代次數(shù)下本文方法和IRLS算法在模擬數(shù)據(jù)上的運(yùn)算時(shí)間如表1所示。從表中可以看出，IRLS算法的運(yùn)算時(shí)間與迭代次數(shù)具有緊密的聯(lián)系。當(dāng)?shù)?次時(shí)，本文方法運(yùn)算時(shí)間比IRLS算法略大，這是由于基于迭代收縮算法的高分辨率Radon變換需要在迭代前對(duì)變換算子L進(jìn)行一次廣義逆(LTL)-1LT的計(jì)算，并將其儲(chǔ)存于內(nèi)存中，使后續(xù)迭代只需要計(jì)算矩陣向量相乘以及向量之間的運(yùn)算即可，這一特性也是本文方法能夠?qū)崿F(xiàn)快速運(yùn)算的核心所在。隨著迭代次數(shù)的增加，本文方法運(yùn)算所需時(shí)間變化平緩，而IRLS算法的運(yùn)算時(shí)間隨迭代次數(shù)的增加而陡然變化。通過(guò)對(duì)比可以看出，針對(duì)此模擬數(shù)據(jù)，當(dāng)?shù)螖?shù)超過(guò)20次后，本文方法的計(jì)算效率是常規(guī)IRLS算法的兩倍左右。

表1 兩種Radon變換方法運(yùn)算時(shí)間 s

3 實(shí)際資料應(yīng)用

選取中國(guó)西部M探區(qū)實(shí)際地震資料測(cè)試本文方法。M探區(qū)溶蝕孔洞和裂縫儲(chǔ)層十分發(fā)育，利用傳統(tǒng)的成像處理流程難以對(duì)儲(chǔ)層內(nèi)部的多種繞射目標(biāo)體精確成像。

首先，利用本文方法對(duì)某一成像點(diǎn)下Kirchhoff疊前深度偏移(PSDM)產(chǎn)生的傾角域CIGs進(jìn)行反射波與繞射波的分離測(cè)試，如圖10所示。

圖10a和10b分別為分離前實(shí)際地震資料的傾角域CIGs和經(jīng)過(guò)本文方法分離后的繞射波CIGs。從圖中可以看出，全波場(chǎng)CIGs中的強(qiáng)反射能量被明顯壓制，盡管在分離后的繞射波CIGs中(圖10b)無(wú)法清晰地觀察到分離出的弱繞射能量，但是強(qiáng)能量的反射同相軸被明顯削弱，說(shuō)明了本文方法分離反射波與繞射波是有效的。

圖10 基于迭代高收縮分辨率Radon變換的CIGs波場(chǎng)分離(a)實(shí)際地震資料的全波場(chǎng)CIGs； (b)波場(chǎng)分離后的繞射波CIGs

進(jìn)一步將本文方法應(yīng)用于該地區(qū)Inline1320測(cè)線的傾角域CIGs，圖11顯示了測(cè)線Inline1320的0～60°傾角域CIGs。由圖可見(jiàn)，強(qiáng)反射波能量呈雙曲線形態(tài)，繞射能量則呈水平直線或者曲率較小的擬線性同相軸。

圖12為經(jīng)過(guò)Radon變換壓制反射波后分離得到的繞射波傾角域CIGs。由圖可見(jiàn)，經(jīng)過(guò)本文方法處理后，大部分強(qiáng)反射能量同相軸被去除，而同一深度點(diǎn)的繞射波被較好地保留。

該測(cè)線波場(chǎng)分離前的PSDM剖面如圖13a所示，利用本文方法分離的繞射波成像剖面如圖13b所示。在原始疊加剖面(圖13a)中，由于強(qiáng)反射能量的干擾，繞射體的位置難以識(shí)別。從圖13b中可以看出，上覆地層反射能量殘余多，中深層反射殘余較弱(幾乎全部去除)，主要原因在于上覆地層采集腳印明顯，前期CDP道集預(yù)測(cè)去除反射能量不多，而下伏地層反射波能量得到準(zhǔn)確壓制。相較于原始偏移剖面，經(jīng)過(guò)繞射分離后得到的繞射波成像剖面(圖13b)中，強(qiáng)反射同相軸被濾除，淹沒(méi)的繞射體位置凸顯出來(lái)(圖13b中黃色圓圈所示)。因此，經(jīng)過(guò)迭代收縮高分辨率Radon變換方法的處理，繞射體的成像結(jié)果被凸顯出來(lái)，說(shuō)明本文方法適用于地下小尺度非均質(zhì)體的識(shí)別和精確定位。

圖11 Inline1320測(cè)線傾角域CIGs(0～60°)

圖12 Inline1320測(cè)線Radon變換后的傾角域CIGs(0～60°)

圖13 中國(guó)西部M探區(qū)實(shí)際資料(Inline1320測(cè)線)處理結(jié)果(a)波場(chǎng)分離前PSDM剖面； (b)波場(chǎng)分離后繞射波成像剖面

4 結(jié)論

本文基于迭代收縮高分辨率Radon變換的繞射波分離成像方法，利用反射波與繞射波在傾角域CIGs中的的運(yùn)動(dòng)學(xué)特征差異，并通過(guò)構(gòu)建一個(gè)常速介質(zhì)模型以說(shuō)明兩者之間的差異，進(jìn)一步利用反射波與繞射波的形態(tài)差異在傾角域CIGs中分離繞射波場(chǎng)。根據(jù)模擬數(shù)據(jù)和實(shí)際地震資料處理結(jié)果得出以下結(jié)論：

(1)在傾角域CIGs中，反射波同相軸總是呈現(xiàn)開(kāi)口向上的“笑臉”狀曲線，繞射同相軸則表現(xiàn)為擬線性，根據(jù)兩者曲線形態(tài)的差異可以實(shí)現(xiàn)波場(chǎng)分離；

(2)迭代收縮高分辨率Radon變換對(duì)反射和繞射波同相軸在Radon域的聚焦能力強(qiáng)，分辨率高，且在計(jì)算過(guò)程中采用收縮閾值算法，不僅提高了反演的收斂速度，而且較常規(guī)方法具有更好的抗噪性；

(3)本文方法雖然能夠有效分離反射波與繞射波，并利用繞射信息揭示了地下小規(guī)模繞射體，可提高地震資料解釋的精度，但由于算法中需要較多的人工干預(yù)(如迭代步長(zhǎng)、閾值算子等的選取)，存在一定的局限性。