繞射波疊前共虛震源道集分離方法

楊城增 張宣堂 盛同杰 黎小偉 馬運利 劉 韜

(①中國石化華北油氣分公司勘探開發研究院，河南鄭州 466400；②中國礦業大學(北京)地球科學與測繪工程學院，北京 100083；③中國石化石油勘探開發研究院，北京 100083)

0 引言

繞射波作為地下局部不連續體的地震響應，攜帶了斷裂、溶洞、尖滅、陷落柱等小尺度非均質體的重要信息，是提高地震勘探分辨率的重要工具[1]。然而，在地震記錄中，繞射波能量比反射波能量低1～2個數量級，在偏移成像時，與繞射相關的小尺度異常體常被強能量的反射體掩蓋，難以識別[2-3]；同時，在以反射波為有效信號的傳統地震數據處理流程中，繞射波也通常作為噪聲被壓制[4]。因此，對繞射波與反射波進行分離，并對繞射波單獨成像，成為地下小尺度非均質體高分辨率成像的關鍵。

20世紀50年代，地球物理學家Krey[5]在研究含有斷層或尖滅的反射地震記錄時發現了繞射波。近年來，為提高小尺度地質體的成像分辨率，各種繞射波分離方法相繼被提出。根據分離繞射波時所處的地震數據處理階段，可以將繞射波分離方法分為以下兩類。

(1)在偏移前進行繞射波分離的方法。Nowak等[6]根據共炮點道集中繞射波與反射波同相軸橫向位置差異，使用加權拉東變換，在變換域中提取繞射波。平面波破壞濾波器(PWD)在繞射波分離中應用廣泛，Fomel[7]改進了該方法，只需局部平面波場斜率一個參數。孔雪等[8]根據繞射波與反射波在平面波記錄上的幾何形態差異，應用PWD壓制反射波。朱生旺等[9]對共炮點道集記錄進行平面波分解，首先利用局部傾角濾波方法分離繞射波大傾角信息成分，再利用預測反演方法分離繞射波低傾角信息成分，將兩次分離結果相加得到最終分離的繞射波。魏巍等[10]在共炮點道集上通過動校正拉平反射波，增加其相干性，而繞射波無法被拉平，相干性弱；然后根據繞射波與反射波相干性差異，利用奇異值分解技術提取繞射波。Rad等[11]將共反射波面元方法擴展至疊前域，提出了一種基于波前屬性的疊前繞射分離流程。Gong等[12]結合Radon變換、Stolt偏移和反偏移，在炮域分離繞射波。此外，多道奇異譜分析[13]、相干波場疊加方法[14]、在線字典學習方法[15]等也被用于繞射波分離。

(2)在偏移后進行繞射波分離的方法。Khai-dukov等[16]提出了聚焦—反聚焦方法分離繞射波。Moser等[17]利用反穩相濾波器修正Kirchhoff偏移核函數，實現繞射波分離與成像。Silvestrov 等[18]利用逆時偏移生成共成像點道集，根據反射波與繞射波在共成像點道集上的傾角差異分離反射波。Zhao等[19]根據馬氏距離構建了振幅衰減函數，據此對Kirchhoff核函數進行修正，有效分離出繞射波。Yu等[20]基于局部平面波方程提出了數據規則化方法，實現了繞射波精確、穩定提取。在方位傾角道集中，反射波具有稀疏特性，而繞射波則表現為低秩分量，根據這種差異，Zhao等[21]提出一種稀疏低秩矩陣分解方法提取繞射波。Li等[22]構建了垂直旅行時差道集，在該道集中，繞射波是平坦的，而反射波表現為上凸形態，根據形態差異分離繞射波。汪天池等[23]證明了逆時偏移產生的傾角道集比Kirchhoff積分偏移產生的傾角道集具有更強的反射聚焦能力，更適合繞射波分離。羅騰騰等[24]在傾角域共成像點道集中，引入迭代收縮高分辨率Radon變換分離繞射波，具有良好的抗噪能力。劉培君等[25]通過構建反穩相濾波算子，在深度域進行繞射波分離與成像。

傳統疊前繞射波分離方法大多在共炮檢距道集開展，但由于深層繞射波能量分布范圍大，因此造成繞射波能量損失。針對此問題，本文提出一種疊前繞射波分離方法，通過共虛震源變換將共炮點道集轉換成共虛震源道集，并在共虛震源道集中進行偏移—反射波去除—反偏移處理，提取繞射波信號，解決在地震波場交叉和相切情況下繞射波弱信號難以分離的問題，最后通過共虛震源逆變換得到共炮點道集繞射波場。

1 方法原理

1.1 共虛震源變換

以均勻介質模型描述反射波和繞射波的旅行時特征(圖1)。根據反射定律，反射波旅行時可表示為

圖1 反射波和繞射波射線路徑及聚焦特性示意圖

(1)

式中：v為地震波傳播速度；h為炮檢距；d為炮點到反射界面法向距離；θ為界面傾角。在共炮點道集中，炮點坐標固定，h變化，上式可寫為

(2)

該旅行時曲線的頂點坐標為(-2dsinθ，2dcosθ/v)，漸近線為±1/v。

繞射波旅行時可以表示為

(3)

式中：d1表示炮點至繞射體D的距離；φ為炮點至繞射體的射線角度。繞射波旅行時曲線的頂點坐標為(-d1sinφ，(d1+d1cosφ)/v)，漸近線為±1/v。

根據式(2)和式(3)，繞射波(圖2的藍色箭頭所示)和反射波(圖2的紅色箭頭所示)都表現為雙曲線形態，在共炮點道集中分離繞射波具有一定難度。

圖2 共炮點道集中繞射波和反射波形態示意圖[26]

事實上，繞射波和反射波具有不同的聚焦性。如圖1a所示，反射波在虛震源SV處聚焦，其旅行時等于從虛震源至接收點的時間。繞射波則在真正的繞射體位置聚焦，其旅行時等于從炮點到繞射體、再到接收點的時間(圖1b)。在一個共炮點道集中，同一個反射界面產生的反射波可以視為來自于一個公共虛震源，而繞射波則不能，因此可以利用反射波旅行時，將反射波映射到對應的虛震源位置上。通過這種方式，共炮點道集被轉換成一個新的道集——共虛震源道集(CVSG)。

假設地下某一位置點是虛震源，將共炮點道集轉換為CVSG可以表示為

δ{t-tR[S,SV(x,z),R]}dt

(4)

式中：x、z分別代表虛震源的水平位置和深度；R表示接收點；U表示記錄時間為t的共炮點道集數據；V是對應的CVSG數據；在各向同性或等效介質中，由震源S激發、R點接收、對應(x,z)處虛震源的反射波走時可以簡單地表示為

tR[S,SV(x,z),R]=

(5)

其中：t0是零炮檢距雙程旅行時，等于從虛震源SV到震源S的時間；W0為波前半徑；α為射線的出射角；vrms為均方根速度。在直射線假設條件下，有

(6)

式中xS表示震源的水平坐標。根據式(6)，可求得式(4)中的反射波旅行時。值得注意的是，SV是一個假想的震源，并沒有真實的速度信息，使用均方根速度vrms計算從虛震源SV到接收點的旅行時。在簡單情況下，反射波的出射角等于界面的真傾角，因此，可以用已知的實際傾角范圍限制SV的位置。

在共炮點道集向CVSG轉換的過程中存在一個問題：如圖1a所示，當接收點位于界面延伸至地表位置Rm左側時，可以計算從虛震源SV至接收點的反射旅行時；但是當接收點位于Rm右側時，無法觀察到反射波。一般情況下，可以用S和SV計算Rm的近似位置處理這一問題。轉換過程中存在的另一個問題是：當虛震源和接收點相距較遠時，會產生拉伸畸變現象。如圖3和式(4)所示，CSVG中的數據點和共炮點道集中的數據點具有一一對應的映射關系。當時差Δt小于時間采樣間隔dt時，將會出現拉伸變形的現象。在均勻介質中，時差Δt可以表示為

圖3 炮集U和共虛震源道集V的映射關系示意圖[23]j是時間t的離散序號

Δt=tR[S,SV(x,zi+1),R]-tR[S,SV(x,zi),R]

(7)

式中：i為深度z離散序號；LRSV(x,zi)表示接收點R到虛震源SV(x,zi)的距離。

當虛震源和接收點相距越遠，Δt的值越小，拉伸變形就越嚴重，在界面陡峭和淺層情況下尤為突出。參考動校正的處理方法，當拉伸畸變無法接受時，切除變形嚴重的部分。

共炮點道集和CVSG具有精確的映射關系，可以通過逆變換由CVSG恢復炮集數據

U(S,R,t)=

(8)

假設地下每一個點都是虛震源點，則可逐點給定坐標對共炮點道集數據進行共虛震源變換，來自同一反射面的反射事件將被校正到同一個虛震源。本文并不在虛震源位置進行反射波的聚焦疊加，而是依次排列成為線性同相軸(圖4紅色箭頭所示)。繞射波無虛震源，所以不會在共虛震源道集中表現為水平線性，而是表現為彎曲狀(圖4藍色箭頭所示)。

圖4 CVSG中繞射波和反射波形態示意圖[26]

1.2 繞射波分離方法

利用Stolt偏移、正則化PWD方法和反偏移在CVSG中提取繞射波弱信號。

Stolt偏移方法是目前最快速的偏移算法，廣泛應用于地震數據處理。為了便于Stolt偏移，對共虛震源道集進行深時轉換，轉換后的CVSG在三維情況下可以寫為V(x,y,z=0,t)，其傅里葉變換為

exp(ikxx+ikyy-iωt)dxdydt

(9)

則成像數據M(x,y,z)的傅里葉變換為

(10)

式中vm為偏移速度。

式(10)即為Stolt偏移公式。圖5為圖4的偏移結果，反射波表現為平面波，繞射波收斂為繞射點。

圖5 共虛震源道集偏移結果

在CVSG經Stolt偏移后，反射波表現出線性連續特征具有平面波的性質，符合平面波假設，可以通過PWD濾波器預測[7]。然而，在偏移數據中，繞射波收斂為繞射點，繞射波不符合平面波假設，難以通過PWD濾波器預測，在預測過程中，繞射波會作為預測誤差得到。因此，可以使用PWD方法得到偏移剖面中的繞射波

C(σ)M=Q(σ)

(11)

式中：M=(m1,m2,…,mN)，其中mi為第i道偏移數據，N為總道數；Q(σ)為預測誤差，σ為斜率矢量；C(σ)是PWD算子，定義為

(12)

其中：I為單位矢量；Bl,n(σl)表示由第l道偏移數據預測第n道偏移數據的預測誤差算子[7]。

從式(12)中可以看出，PWD算子所需要的唯一參數就是局部斜率σ，只需要求解局部斜率，再將其代入PWD算子，就可以利用式(11)得到繞射波。因此，獲得繞射波的關鍵就是估計局部斜率。

求解局部斜率是一個迭代過程，可以使用高斯—牛頓法[7]，在Q(σ)近似為零條件下，有

C′[σ(K)]Δσ(K)M+C[σ(K)]M=0

(13)

式中：K為迭代次數；C′[σ(K)]為C[σ(K)]關于σ(K)的微分；Δσ(K)為每次迭代中傾角斜率矢量的更新步長，局部斜率矢量的更新公式為σ(K+1)=σ(K)+Δσ(K)。

上述是傳統的求解局部斜率的過程，在求解過程中，并未考慮與斷層等不連續地質體有關的地震數據對反演過程干擾作用和地震反射信號平滑特征，反演求解不穩定，局部斜率估計不準確，應用不準確的局部斜率分離繞射波，可能會導致反射波殘留、分離出的繞射波不完整等問題。因此，為了避免上述問題，本文采用Yu等[20]提出的正則化局部斜率估計模型

ψ〈L1σ,σ〉+β〈L2σ,σ〉→min

(14)

式中：ψ、β均為正則化因子，用于平衡預測地震道和觀測地震道的擬合情況；〈·,·〉表示兩個向量的內積；L1為負拉普拉斯算子，即

(15)

L2為差分算子，即

(16)

上述正則化傾角估計模型保證了反演求解穩定性，也確保三維地震強反射信號估計準確性[20]。

為快速求解上述局部斜率估計模型，目標函數JΨ,β(σ)的梯度可表示為[20]

(17)

模型更新方案為

σ(K+1)=A{σ(K)-ξKg[σ(K)]}

(18)

其中：ξK為更新步長；對于任意ξK>0，A(·)表示斜率矢量的可行集投影。

在得到局部斜率后，可由式(11)得到預測誤差Q(σ)，即偏移數據中的繞射波。

在利用上述方法去除偏移數據中的反射波后，再進行反偏移

(19)

通過逆傅里葉變換可得CVSG繞射波

exp(-ikxx-ikyy+iωt)dkxdkydω

(20)

最后，可以通過式(8)將繞射波CVSG變換為繞射波共炮點道集。

本文方法的處理流程如圖6所示。

圖6 疊前繞射波分離流程

2 模型和實際數據應用

構建巖溶—沖溝模型驗證本文繞射波分離方法的效果，并與常規的PWD方法對比，闡明本文方法的優勢；通過實際資料處理驗證本文方法的適用性。

2.1 巖溶—沖溝模型

巖溶-沖溝模型(圖7)尺寸為9000m×4000m，1000～2300m深度范圍有兩個薄層，2000～3500m深度范圍內隨機分布著一些形狀、尺寸不同的巖溶體，最小的直徑為4m，最大的直徑為100m，填充速度為4200m/s，圍巖的速度為4700m/s。圖8為模型的反射系數剖面，其中小尺度巖溶構造清晰，可用于對比分析不同分離方法的效果。對該模型使用Kirchoff正演算法生成合成數據，震源使用主頻為25Hz的Ricker子波；共901炮，炮點位置從0到9000m變化，炮間距為10m；每炮901道，道間距為10m，采樣間隔為4ms，其共炮檢距剖面如圖9a所示。

圖7 巖溶—沖溝速度模型

圖8 巖溶—沖溝模型的反射系數剖面

應用PWD法對圖9a所示的共炮檢距剖面進行繞射波分離，結果如圖9b所示。PWD方法有效去除了平滑、連續的反射波，但繞射波也遭到嚴重的損傷，尤其在繞射波頂點，即繞射波與反射波相干處。

對模擬數據應用本文方法分離繞射波，并將分離結果(共炮點道集)抽成共炮檢剖面，如圖9c所示。本文方法有效去除了剖面中平滑、連續的反射波，比PWD方法更完整地保留了繞射波信息，為繞射波的清晰成像奠定了基礎。

圖9 巖溶—沖溝模型兩種繞射波分離方法的共炮檢距剖面對比

圖10a是PWD法分離的繞射波偏移成像結果，雖然繞射波弱信號在一定程度上突顯，但由于分離的繞射波不完整，繞射點成像不清晰，信噪比較低。圖10b是使用本文方法分離的繞射波偏移成像結果，得益于反射波去除干凈、繞射波保留完整，繞射點成像清晰完整，小尺度巖溶能夠突顯。

圖10 巖溶—沖溝模型兩種方法分離的繞射波偏移剖面對比

2.2 實際資料應用

實際資料采集于內蒙古杭錦旗。該區奧陶系縫洞型儲集體發育，具有埋藏深、尺度小、非均質性強等特點，繞射波響應較弱。在常規反射成像剖面(圖11a)上縫洞型儲集體反映不清晰，精細識別與刻畫難度大，目的層的沖溝和巖溶體有顯示(紅色圓圈所示)。由于小尺度巖溶體的存在，平滑線性的反射層呈現向下彎曲的形態，但由于背景強反射的影響，分辨率較低，很難識別其中的小尺度巖溶體。對實際資料使用本文方法分離繞射波，并對分離結果進行成像(圖11b)。由于強反射已被移除，原來被強反射掩蓋的小尺度巖溶突顯(紅色圓圈內)，分辨率和聚焦性顯著提高。相干體屬性可用于斷裂識別，圖12a是實際資料的反射波相干屬性切片，雖然可以看到巖溶體和斷裂，但小尺度巖溶體和斷裂難以識別。圖12b為繞射波能量屬性切片，可識別的小尺度巖溶體和斷裂明顯增多，小尺度巖溶體和斷裂的邊界、形態也更清晰(紅色圓圈內)。

圖11 實際數據成像剖面對比

圖12 反射波相干屬性(a)與繞射波能量屬性(b)的切片對比

3 結論

本文通過共虛震源變換利用共炮點道集生成CVSG，對CVSG應用Stolt偏移使雙曲形態的繞射波能量收斂，使用正則化PWD去除反射波，再利用反偏移重新得到CVSG，最后利用共虛震源反變換得到共炮點道集繞射波場。本文方法解決了反射與繞射波相切或相交時繞射波難以分離的問題，避免了傳統方法在共炮檢距道集上進行分離造成的深層繞射波能量損失。數值模型測試驗證了本文提出的疊前繞射波分離方法可以提升小尺度地質體成像分辨率，實際資料的應用進一步證明了本文方法對探明縫洞型儲集體具有積極意義。本文方法中共虛震源正、反變換雖然數學上是可逆的，但實現過程中有能量損失，如何減小這種損失是下一步的研究方向。