速度誤差和地震噪聲對最小二乘逆時偏移的影響分析

馬振，孫成禹，彭鵬鵬，姚振岸

(1.中國石油大學(xué)(華東) 地球科學(xué)與技術(shù)學(xué)院，山東 青島 266580； 2.雅安市雨城區(qū)自然資源和規(guī)劃局，四川 雅安 625000； 3.東華理工大學(xué) 地球物理與測控技術(shù)學(xué)院，江西 南昌 330013)

0 引言

隨著計算機技術(shù)的快速發(fā)展以及勘探開發(fā)目標的日益精細化和復(fù)雜化，對地震疊前深度偏移技術(shù)的精確性、成像分辨率和保幅性均提出了更高要求，因此在常規(guī)逆時偏移的基礎(chǔ)上基于反演理論的最小二乘逆時偏移應(yīng)運而生。與常規(guī)逆時偏移相比，最小二乘逆時偏移具有高分辨率、高保幅性以及噪聲少的優(yōu)勢，對于復(fù)雜的地質(zhì)構(gòu)造成像更加精確。然而，最小二乘逆時偏移對于偏移速度場的精度要求較高，其精確成像依賴準確的偏移速度，速度估算不準確時點繞射在偏移剖面表現(xiàn)出上翹或下拉現(xiàn)象，即通常說的偏移過量或偏移不足。同時在實際的勘探中地震資料中往往混雜著大量的噪聲，噪聲對速度分析的精確性有直接影響，同時也影響地震偏移成像的品質(zhì)，分析噪聲對最小二乘逆時偏移的影響程度也是有必要的。因此研究最小二乘逆時偏移在不同條件下的成像特征，有助于處理人員做出適當?shù)恼{(diào)整，使最小二乘逆時偏移更好的為復(fù)雜構(gòu)造、巖性和地層圈閉油氣藏的勘探與開發(fā)服務(wù)。

Lebras等[1]最先提出最小二乘偏移反演的思想，Lambaré等[2]進行了完善，Nemeth等[3]在地震偏移中引入最小二乘算法，進行Kirchhoff最小二乘偏移。逆時偏移[4-7]成像技術(shù)出現(xiàn)之后，國內(nèi)外學(xué)者在逆時偏移的基礎(chǔ)上發(fā)展了基于反演理論的最小二乘逆時偏移[8-10]。在偏移技術(shù)[11-13]的不斷發(fā)展過程中，國內(nèi)外學(xué)者對其不同條件下的成像特征做了相關(guān)的研究，Versteeg等[14]分析了疊前深度偏移對速度模型的敏感性，認為平滑后的真實速度模型仍可用作疊前深度偏移。Lines等[15]在雙曲時差方程的基礎(chǔ)上，研究了速度誤差、炮間距等對成像深度的影響。Grubb等[16]分析了AVO和偏移結(jié)果同相軸深度受速度不確定性的影響。葉月明等[17]分析了幾種典型疊前深度偏移方法對偏移速度的敏感性，為實際情況選取合適的偏移方法提供依據(jù)。Zhu等[18]研究了速度誤差對疊前和疊后偏移成像的影響，從數(shù)學(xué)和幾何方面分析了偏移剖面中出現(xiàn)的點繞射上翹和下拉現(xiàn)象。Parkes等[19]提出偏移速度場對成像結(jié)果影響較大，偏移成像過程中重點應(yīng)放在提高輸入速度場的精確性方面。熊曉軍等[20]分析了速度誤差在偏移剖面中引起的假象，指出偏移速度存在誤差不僅會帶來深度誤差，同時對構(gòu)造體的形態(tài)產(chǎn)生影響。Shan等[21]研究了逆時偏移對偏移速度的敏感性，指出RTM在強波阻抗界面處更容易產(chǎn)生虛假反射。勾麗敏等[22]研究了隨機噪聲和規(guī)則噪聲對疊前深度偏移層速度估計精度的影響。鄒少峰等[23]分析了地震信號中的噪聲對偏移結(jié)果的影響，得出疊前深度偏移比疊加具有更強的壓制噪聲能力。

為了分析速度誤差和地震噪聲對最小二乘逆時偏移的影響，本文基于Marmousi模型進行最小二乘逆時偏移成像，將其成像結(jié)果與逆時偏移成像結(jié)果進行了對比，分析最小二乘逆時偏移在高分辨率成像中的優(yōu)勢;同時針對偏移過程中常見的速度系統(tǒng)誤差、速度平滑以及地震資料中存在噪聲的情況，本文利用不同偏移速度及不同含噪的地震數(shù)據(jù)進行最小二乘逆時偏移成像，對比分析不同條件下的成像結(jié)果同時與逆時偏移成像結(jié)果進行對比，得到幾點結(jié)論，可為實際生產(chǎn)提供基礎(chǔ)指導(dǎo)。

1 方法原理

1.1 最小二乘逆時偏移

在常密度介質(zhì)的前提下，波動方程在頻率域內(nèi)表示為：

[2+ω2s2(x)]u(x,xs,ω)=-f(ω)δ(x-xs) ,

(1)

其中:ω為頻率，s(x)表示介質(zhì)的慢度，u表示總波場，f(ω)為震源函數(shù)，δ(x-xs)為脈沖函數(shù)。

由波動方程可知，在研究區(qū)域內(nèi)速度擾動與波場的擾動存在如下關(guān)系：

(2)

相應(yīng)波場為：

u(x,xs,ω)=u0(x,xs,ω)+us(x,xs,ω) 。

(3)

其中:s0(x)為背景參考慢度，Δs(x)為介質(zhì)慢度的擾動，u0(x,xs,ω)為背景慢度產(chǎn)生的波場，us(x,xs,ω)為擾動慢度產(chǎn)生的擾動波場。

同時,背景波場u0(x,xs,ω)滿足下式：

[2+ω2s2(x)]u0(x,xs,ω)=-f(ω)δ(x-xs)。

(4)

將式(2)(3)(4)代入式(1)，基于線性Born近似，可以得到：

(5)

引入背景介質(zhì)中Green函數(shù)：

(6)

可以得到背景波場u0(x,xs,ω)及擾動波場us(x,xs,ω)：

u0(x,xs,ω)=f(ω)G(x,xs,ω) ,

(7)

G(xr,x,ω)dx。

(8)

其中:xr為檢波點位置，xs為震源位置，G(x,xs,ω)表示從震源位置xs傳播到散射點x的Green函數(shù)，G(xr,x,ω)表示從散射點x傳播到檢波點xr的Green函數(shù)。

式(8)僅表示波傳播過程中，忽略了二階以上的散射波的一階散射波的傳播現(xiàn)象。Born近似描述的物理含義為：在背景波場u0(x,xs,ω)中沒有散射場；在擾動波場中僅包含一次散射場。目前對最小二乘逆時偏移的討論主要基于式(8)，稱之為線性化最小二乘逆時偏移。為了使后續(xù)表達式更加簡潔，將式(8)改寫為：

us=Lm,

(9)

其中:L(ω,xs,xr,x)=ω2f(ω)G(x,xs,ω)G(xr,x,ω)，其中L通常被稱為反偏移算子。為求取模型參數(shù)擾動m，基于最小二乘方法求解(9)式。最小二乘意義下的誤差泛函可表示為：

F(m)=‖us-Lm‖2。

(10)

求解目標函數(shù)，其最小二乘解為：

m=(LTL)-1LTus=H-1LTus,

(11)

式(11)中,H-1表示Hessian逆矩陣，LTus即為逆時偏移的結(jié)果。逆時偏移過程中采用PML吸收邊界，2N階有限差分的穩(wěn)定性條件為：

(12)

其中:Δt為時間步長，Δx和Δz分別為沿x和z軸的空間采樣間隔，v為相應(yīng)坐標點上的速度，a是對應(yīng)的差分系數(shù)。

通過式(11)可以得出， 最小二乘逆時偏移的核心算法是通過判斷偏移結(jié)果與地震數(shù)據(jù)的匹配程度來判斷偏移成像的準確性，同時根據(jù)數(shù)據(jù)殘差對成像結(jié)果進行修正，通過不斷的迭代實現(xiàn)二者的最佳匹配，當數(shù)據(jù)殘差滿足期望值時則可以停止迭代過程，得到更加真實的成像結(jié)果。

1.2 速度誤差對最小二乘逆時偏移的影響

最小二乘逆時偏移是一種高精度的疊前深度偏移成像方法，其對于偏移速度場的精度要求較高，速度估算不準確時，點繞射在偏移剖面表現(xiàn)出上翹或下拉現(xiàn)象，即通常說的偏移過量或偏移不足。在幾何地震學(xué)中，偏移速度與地層深度可以表示為：

(13)

(14)

其中:Δvm表示速度誤差，Δz為深度誤差。從式(14)可以直觀看出，當偏移速度偏大即Δvm>0時，此時Δz>0，成像位置變深；當偏移速度偏小即Δvm<0時，則Δz<0，成像位置變淺。當?shù)刭|(zhì)構(gòu)造復(fù)雜，速度局部變化劇烈時，深度誤差與速度誤差之間的關(guān)系趨于復(fù)雜，但成像位置與實際地層位置相對關(guān)系不變。

1.3 地震噪聲對最小二乘逆時偏移的影響

以積分的角度分析，最小二乘逆時偏移過程可以視為疊加過程，當?shù)卣鹦盘栔写嬖谠肼晻r，對于N次覆蓋的某成像點，從信號和噪聲兩方面可將成像結(jié)果表示為：

(15)

2 數(shù)據(jù)測試

2.1 精確速度成像結(jié)果

本文選取如圖1所示的Marmousi模型對最小二乘逆時偏移的適應(yīng)性進行研究。這里Marmousi模型大小為橫向4.61 km，縱向2.84 km。地震記錄用時間2階空間8階有限差分正演所得。正演模擬時，垂向和水平方向網(wǎng)格大小均為10 m，時間步長 1 ms，記錄時長5 s，子波主頻25 Hz。炮記錄共92炮(起始炮點位置0.05 km，終止炮點位置4.6 km)，每炮共461道接收，每道5 000個采樣點，時間采樣率1 ms，記錄5 s，道間距10 m，炮間距均為50 m。圖2分別為逆時偏移和最小二乘逆時偏移迭代20次的精確成像結(jié)果。當用精確的速度作為偏移速度時，逆時偏移和最小二乘逆時偏移都能正確成像。在中部背斜構(gòu)造處，由于背斜右翼較為陡峭，最小二乘逆時偏移構(gòu)造分層更加明顯，分辨率更高，且最小二乘逆時偏移對斷層構(gòu)造及小尺度構(gòu)造展布更加清晰。

圖1 Marmousi模型Fig.1 Marmousi model

a—逆時偏移;b—最小二乘逆時偏移a—reverse time migration;b—least squares reverse time migration圖2 精確偏移速度成像結(jié)果Fig.2 Accurate migration velocity imaging results

2.2 速度誤差的影響

2.2.1 速度系統(tǒng)誤差的影響

分別利用Marmousi模型速度整體增大5%和減小5%作為偏移速度進行最小二乘逆時偏移，圖3為最小二乘逆時偏移迭代20次后的成像結(jié)果。對于偏移速度整體減小5%的成像結(jié)果，剖面中同相軸的位置變淺；斷點處繞射波不能收斂，出現(xiàn)向下彎曲的同相軸，影響斷層的邊界識別；剖面1.5 km以下的部分同相軸較為雜亂，構(gòu)造形態(tài)失真，難以反映真實的構(gòu)造形態(tài)，同時對地質(zhì)層位追蹤解釋造成困難。對于偏移速度整體增大5%的成像結(jié)果，剖面中同相軸位置變深；斷點處繞射波不能收斂，出現(xiàn)向上彎曲的同相軸，對斷層的邊界識別產(chǎn)生不利影響；剖面1.5 km以下的部分同相軸較為雜亂，同樣很難反映真實的構(gòu)造形態(tài)。中深部左側(cè)的高速楔狀體上界面精確位置位于2.25 km處，偏移速度減小時，上界面成像位置為2.1 km，偏移速度增大時，上界面成像位置為2.55 km，偏移速度的增大顯然比偏移速度減小引起的深度誤差更大。但是對于該高速楔狀體的形態(tài)而言，偏移速度減小時，頂?shù)捉缑鎶A角明顯變小，相比偏移速度增大的成像結(jié)果形態(tài)失真程度更大。圖4為相同誤差偏移速度下逆時偏移的成像結(jié)果，通過與最小二乘逆時偏移成像結(jié)果對比可以發(fā)現(xiàn)，逆時偏移和最小二乘逆時偏移成像結(jié)果都受速度誤差影響較大，兩者都難以反映地下實際構(gòu)造形態(tài)。

a—偏移速度減小5%;b—偏移速度增加5%a—migration velocity decreases by 5%;b—migration velocity increases by 5%圖3 不同系統(tǒng)誤差偏移速度的最小二乘逆時偏移成像結(jié)果Fig.3 Least squares reverse-time migration imaging results with different system error migration velocities

a—偏移速度減小5%;b—偏移速度增加5%a—migration velocity decreases by 5%;b—migration velocity increases by 5%圖4 不同系統(tǒng)誤差偏移速度的逆時偏移成像結(jié)果Fig.4 Reverse time migration imaging results with different system error migration velocities

圖5為不同偏移速度進行偏移成像時，最小二乘逆時偏移迭代收斂的情況。從圖中可以明顯看出，利用精確偏移速度進行偏移時，偏移迭代過程泛函收斂快，偏移速度存在誤差時，偏移迭代過程泛函收斂很慢，得到的偏移結(jié)果的反偏移合成數(shù)據(jù)和原始地震記錄之間存在的誤差大，偏移結(jié)果難以逼近地下的真實構(gòu)造形態(tài)。

圖5 不同偏移速度誤差下歸一化數(shù)據(jù)殘差收斂曲線對比Fig.5 Comparison of normalization data residual convergence curves under different migration velocity errors

2.2.2 速度平滑的影響

偏移速度的平滑處理相當于對其進行濾波，常規(guī)的空間域平滑處理為簡單平均法，即求取臨近速度值的平均值，平滑效果與鄰域的大小直接相關(guān)。本文采用二維高斯函數(shù)

(16)

進行速度的平滑處理，同時用平滑后的速度相對于真實速度的整體差異來衡量平滑程度，定義平滑后的速度的誤差率e：

(17)

其中，vture為真實的速度，vsmooth為平滑后的速度，nx,nz為二維速度模型的大小。

分別使用參數(shù)σ為2、5、10，利用式(16)對Marmousi模型速度進行平滑處理得到圖6的速度模型，此時誤差率e的值分別為2.7%、5.24%、7.28%。

分別利用不同平滑程度的偏移速度進行最小二乘逆時偏移，圖7為不同偏移速度的迭代20次的最小二乘逆時偏移成像結(jié)果。通過成像結(jié)果的對比可知，對偏移速度進行較小的平滑處理之后，即e=2.7%時，最小二乘逆時偏移成像質(zhì)量得到提高，剖面更加清晰。當e=5.24%時，在剖面的背斜及深部潛山處出現(xiàn)錯誤畫弧同相軸。當e=7.28%時，在斷層處，斷面出現(xiàn)微弱的扭曲，在深部潛山錯誤畫弧同相軸現(xiàn)象加劇。整體而言，偏移速度場進行較小的平滑處理，對最小二乘逆時偏移成像剖面質(zhì)量影響不大；當速度平滑程度較大時，會產(chǎn)生錯誤的同相軸，將會對后續(xù)的構(gòu)造解釋，層位分析產(chǎn)生不利影響。圖8為不同偏移速度的逆時偏移成像結(jié)果，可以發(fā)現(xiàn)較小的平滑處理時，即e=2.7%時，逆時偏移成像質(zhì)量同樣得到提高；在e=5.24%時，逆時偏移斷層斷面出現(xiàn)微弱扭曲現(xiàn)象，e=7.28%時，斷層斷面扭曲現(xiàn)象加劇，剖面中出現(xiàn)畫弧現(xiàn)象；相比逆時偏移而言，最小二乘逆時偏移對平滑程度較大的速度誤差更加敏感。

a—真實速度(e=0%);b—偏移速度模型1(e=2.7%);c—偏移速度模型2(e=5.24%);d—偏移速度模型3(e=7.28%)a—accurate speed(e=0%);b—migration velocity model 1(e=2.7%);c—migration velocity model 2(e=5.24%);d—migration velocity model 3(e=7.28%)圖6 不同平滑程度的Marmousi模型Fig.6 Marmousi model with different smoothness

a—速度未做平滑;b—e=2.7%的平滑速度;c—e=5.24%的平滑速度;d—e=7.28%的平滑速度a—velocity without smoothness;b—e=2.7% smoothing speed;c—e=5.24% smoothing speed;d—e=7.28% smoothing speed圖7 不同平滑程度的偏移速度最小二乘逆時偏移成像結(jié)果Fig.7 Least squares reverse-time migration imaging results with different smooth migration velocity

a—速度未做平滑;b—e=2.7%的平滑速度;c—e=5.24%的平滑速度;d—e=7.28%的平滑速度a—velocity without smoothness;b—e=2.7% smoothing speed;c—e=5.24% smoothing speed;d—e=7.28% smoothing speed圖8 不同平滑程度的偏移速度逆時偏移成像結(jié)果Fig.8 Reverse time migration imaging results with different smooth migration velocity

圖9為不同平滑程度的偏移速度進行偏移成像時，最小二乘逆時偏移迭代收斂的情況?？梢钥闯?隨著偏移速度平滑程度的增大，偏移結(jié)果的反偏移合成數(shù)據(jù)和原始地震記錄之間的匹配程度也隨之變差。由于速度平滑是速度的隨機變化，速度場誤差的數(shù)學(xué)期望值較小，因此平滑后的速度相比速度系統(tǒng)誤差對最小二乘逆時偏移的影響小。

圖9 不同平滑程度偏移速度下歸一化數(shù)據(jù)殘差收斂曲線對比Fig.9 Comparison of normalized data residuals convergence curves with different smooth migration velocity

2.3 隨機噪聲的影響

為與實際勘探過程中采集地震數(shù)據(jù)噪聲規(guī)律相符，本文對炮集中的全部地震道加入隨機噪聲，其中隨機噪聲的構(gòu)造步驟如下。

首先,針對炮記錄生成區(qū)間[-1,1]的二維隨機數(shù)：

noise[m,n]=rand[m,n] ,

(18)

其中m表示炮記錄中單道樣點數(shù)，n表示炮集中的道數(shù)。

然后,計算振幅均衡因子scale:

(19)

其中ns為信噪比，S為原始炮集信號中最大振幅，N為噪聲noise的均方根振幅。

最后,輸出信噪比為ns的含隨機噪聲的地震記錄output：

output=signal+scale×noise ,

(20)

其中，signal為原始地震記錄信號。

分別在地震記錄加入信噪比為5和信噪比為1的隨機噪聲，圖10為去除直達波之后加入隨機噪聲的地震記錄。

圖11為利用精確的偏移速度對上述地震記錄進行最小二乘逆時偏移迭代20次得到的成像結(jié)果。當信噪比為5時，最小二乘逆時偏移仍然能夠較好地壓制噪聲，噪聲對偏移剖面影響不大。當信噪比是1時，偏移剖面受到極大的影響，只有強反射軸能夠勉強辨認，成像結(jié)果較差。整體而言，當?shù)卣鹳Y料的信噪比較小時，隨著噪聲的增強，偏移剖面中的弱小反射被淹沒在噪聲之中，弱同相軸難以分辨，尤其是斷層斷面、深部潛山嚴重受噪聲影響，難以反映地下構(gòu)造形態(tài)。對比圖12逆時偏移的成像結(jié)果可以發(fā)現(xiàn)，在信噪比為5時，逆時偏移和最小二乘逆時偏移成像品質(zhì)都較高，抗噪能力相當。當信噪比為1時，可以發(fā)現(xiàn)逆時偏移的抗噪能力優(yōu)于最小二乘逆時偏移，這是因為逆時偏移是利用雙程波進行成像，減弱了噪聲的影響。而最小二乘逆時偏移在進行迭代時，對偏移結(jié)果做反偏移與地震記錄匹配計算誤差泛函，地震噪聲的存在，使得迭代時誤差泛函難以收斂，因而噪聲較大時，最小二乘偏移較逆時偏移更為敏感。

a—原始地震記錄;b—ns=5的地震記錄;c—ns=1的地震記錄a—seismic record without noise;b—seismic records of ns=5;c—seismic records of ns=1圖10 加入不同程度噪聲的地震記錄Fig.10 Seismic records with adding different degrees noise

圖11 不同含噪地震記錄最小二乘逆時偏移成像結(jié)果Fig.11 Least squares reverse-time migration imaging results of different noisy seismic records

圖12 不同含噪地震記錄逆時偏移成像結(jié)果Fig.12 Reverse time migration imaging results of different noisy seismic records

3 結(jié)論與討論

本文分析了速度誤差和地震噪聲對最小二乘逆時偏移成像的影響，基于Marmousi模型進行最小二乘逆時偏移及逆時偏移成像，針對實際資料處理中經(jīng)常出現(xiàn)的偏移速度系統(tǒng)誤差、速度平滑、地震資料噪聲情況分別利用最小二乘逆時偏移進行了偏移成像，并且與其精確成像和對應(yīng)逆時偏移結(jié)果進行了對比，得到了以下幾點認識：

1)在偏移速度存在系統(tǒng)誤差的情況下，最小二乘逆時偏移成像深度會出現(xiàn)誤差，剖面中斷點處繞射波無法收斂，影響斷層的邊界識別；剖面中深部出現(xiàn)雜亂、不連續(xù)的同相軸，形態(tài)失真，很難反映真實的構(gòu)造形態(tài)；逆時偏移和最小二乘逆時偏移成像結(jié)果都受速度誤差影響較大，兩者都難以反映地下實際構(gòu)造形態(tài)。相比速度系統(tǒng)減小的情況，偏移速度增大時引起的深度誤差更大，而對構(gòu)造形態(tài)影響相對較小。

2)適當?shù)钠扑俣绕交梢愿纳谱钚《四鏁r偏移剖面質(zhì)量，一般來說，速度誤差率e小于3%可以適當改善偏移剖面質(zhì)量，大于3%后，則會影響成像品質(zhì)，較大的平滑會產(chǎn)生錯誤的同相軸，相比逆時偏移而言，最小二乘逆時偏移對平滑程度較大的速度誤差更加敏感。

3)最小二乘逆時偏移的成像質(zhì)量會隨著地震資料中噪聲增加而變差，當?shù)卣鹳Y料的信噪比較小時，隨著噪聲的增強，偏移剖面中的弱小反射被淹沒在噪聲之中，弱同相軸難以分辨，尤其是斷層斷面、深部潛山嚴重受噪聲影響，難以反映地下構(gòu)造形態(tài)，在噪聲較大時，最小二乘偏移較逆時偏移更為敏感。

通過本文分析，最小二乘逆時偏移的正確成像依賴精確的偏移速度，最小二乘逆時偏移及逆時偏移受速度系統(tǒng)誤差影響較大，在實際資料的處理過程中，應(yīng)盡量提高偏移速度的精度，同時可以對偏移速度進行適當?shù)钠交?，有利于最小二乘逆時偏移的成像質(zhì)量；對地震數(shù)據(jù)進行最小二乘逆時偏移成像之前，應(yīng)降低其背景噪聲，從而減小其對成像剖面的影響。