繞射波逆時偏移成像方法研究

王曉毅，陳璽,楊振,張江杰

1 中國科學院地質(zhì)與地球物理研究所，中國科學院油氣資源研究重點實驗室，北京 100029 2 中國科學院地球科學研究院，北京 100029 3 南方海洋科學與工程廣東省實驗室(廣州)，廣州 511458 4 廣州海洋地質(zhì)調(diào)查局，廣州 510075 5 中國科學院大學，北京 100049

0 引言

地下介質(zhì)中廣泛發(fā)育的斷層、尖滅、溶洞等中小尺度不連續(xù)體，對油氣生成、運移和存儲具有重要意義.由于繞射波的形成與裂縫等構造相關聯(lián)(黃建平等，2015；欒錫武和李繼光，2018)，其可以作為偏移處理對象來對目標地質(zhì)體進行精細成像.

目前，相關研究主要集中在如何利用繞射波與反射波在運動學上的區(qū)別進行波場分離與成像.常規(guī)的繞射波成像技術以射線類偏移方法為主.其中，針對繞射波較強的情形，對地震數(shù)據(jù)采用繞射增強的疊加算法是一種有效的策略(Faccipieri et al.,2016)，但其難以處理能量較弱的繞射波.利用反射波和繞射波不同的運動學特征，在地震資料中對二者進行分離，同樣能夠達到成像的目的(Merzlikin and Fomel,2017;Rad et al.,2018).然而，當?shù)叵陆橘|(zhì)構造復雜時，反射波在變換空間的形態(tài)不再簡單，極大地提高了分離的難度.此外，還可以在偏移道集中區(qū)分繞射波和反射波(Zhang and Zhang,2014)，該類方法對復雜反射具有很好的適應性，但對于深度域三維成像而言，其運算成本較高.Li和Zhang(2019)提出基于反穩(wěn)相偏移產(chǎn)生的時差道集來甄別和增強繞射波信號，該方法要求地層傾角信息精確，不需要過大內(nèi)存就可實現(xiàn)繞射波成像.

射線類方法在提取地震波的傳播方向和走時上具有靈活簡便的優(yōu)點，但在復雜構造區(qū)域的成像上往往效果有限.基于此類方法的繞射波成像雖然提高了繞射體的成像分辨率，但其成像位置的準確性難以滿足要求.利用波動類方法進行繞射波成像精度較高，需要的計算成本也較大.例如，Liu等(2016)借助相移插值(PSPI)法提取傾角道集，并通過中值濾波實現(xiàn)了繞射點成像.逆時偏移作為一種經(jīng)典的波動方程類成像方法，具有諸多優(yōu)點，如不受地層傾角的限制，能夠靈活處理回轉(zhuǎn)波和棱柱波等.根據(jù)波傳播射線的幾何特征，選擇不同的波場組分進行互相關運算，可以實現(xiàn)對某一類構造的單獨成像(王一博等，2016;Zhang et al.,2019a,b).當入射波沿右下方傳播而出射波沿左上方傳播時，反射界面是正傾角的地層；入射波沿左下方傳播并且出射波沿右上方傳播時，反射界面是負傾角地層.由于在不連續(xù)點處產(chǎn)生的繞射是沿各個方向的，因此在多種成像條件下繞射波的能量都可以得到聚焦.Zhang等(2019a,b)通過一種乘法成像條件(the multiplication imaging condition)，利用正傾角地層剖面和負傾角地層剖面相乘得到單獨的繞射點位置剖面.這一方法較為簡單，但偏移剖面只包括位置信息，忽略了同相軸的幅值和相位信息.

本文基于逆時偏移對繞射波能量進行聚焦，運用多種方法得到單獨的繞射點剖面.首先，本文詳細分析了繞射波和反射波的運動學特征以及各種成像方式的特殊性.借助Hilbert變換波場和Fourier變換，實現(xiàn)波場分離.基于此，我們獲取了單一傾向構造的成像剖面以及更高質(zhì)量的傾角域共成像點道集.隨后，我們通過乘法成像條件獲取了單獨的繞射點剖面并在傾角域拾取了繞射波能量.此外，我們還對無需波場方向分解的正負傾角構造成像公式進行了推導，并提出了一種在成像剖面上快速獲取繞射點信息的方案.

1 基本理論

目前，基于逆時偏移對繞射點和間斷點實現(xiàn)定位的策略分為以下兩種：選擇不同的震源波場和檢波器波場組分互相關以得到正傾角和負傾角地層構造，相乘獲取單獨的繞射點剖面(Zhang et al,2019a,b,2020)；逆時偏移產(chǎn)生共成像點傾角道集，然后在角道集中甄別、拾取繞射波能量(汪天池等,2020).以上兩種方案中，選擇合適組分和提高角道集質(zhì)量都可以利用Hilbert變換實現(xiàn)波場分離.

1.1 Hilbert變換與波場分解

在頻率-波數(shù)域中，根據(jù)圓頻率ω和波矢量k的符號，可以將原始波場按照傳播方向分解為若干個分量.在二維情形下，將時間-空間域的波場轉(zhuǎn)化到頻率-波數(shù)域需要三維傅里葉變換.其中，沿時間軸方向的傅里葉變換需要用到所有時刻的波場快照，這極大地增加了存儲成本并制約了運行效率.

通過Hilbert變換對原始信號進行改造，可以只保留頻譜和振幅譜中的正頻率部分.圖1a是均值為0.25，方差為0.02的高斯子波信號，其振幅譜和相位譜分別如圖1b藍色線和橙色線所示.將原始信號進行Hilbert變換之后，得到以原始信號為實部，圖1c中橙色線段為虛部的復數(shù)信號.圖1d是復數(shù)信號所對應的振幅譜(藍色線)和相位譜(橙色線)，可見復數(shù)信號在負頻率部分振幅為0，正頻率部分振幅擴大為原來2倍.換言之，Hilbert變換使信號只保留正頻率部分.將這一性質(zhì)應用于波場分解，可以避免圓頻率ω符號的判斷，從而減少沿時間軸方向波場的存儲.

圖1 原始信號與其經(jīng)Hilbert變換之后得到的復數(shù)信號(a)原始信號：均值為0.25，方差為0.02的高斯子波；(b)與信號(a)對應的振幅(藍色線)和相位譜(橙色線)；(c)復數(shù)信號：以原始信號為實數(shù)部，以Hilbert變換后信號為虛數(shù)部(如橙色線段所示)；(d)與復數(shù)信號(c)對應的振幅(藍色線)和相位譜(橙色線).注意振幅在負頻率部分為0.Fig.1 The original signal and the complex signal obtained after the Hilbert transform to it(a)The original signal:a Gaussian wavelet with mean value of 0.25 and variance of 0.02;(b)The amplitude and phase spectrum (blue and orange lines)corresponding to the signal (a);(c)The complex signal:taking the original signal as the real part and the Hilbert transformed signal (the orange line)as the imaginary part;(d)The amplitude and phase spectrum (blue and orange lines)corresponding to the signal (c).Note that the amplitude is zero in the negative frequency part.

對于聲波而言，標量波動方程可以表述為

(1)

(2)

(3)

在這里，i=lu,ru,ld或者rd分別表示左-上行波，右-上行波，左-下行波和右-下行波.四者對應的積分區(qū)域(kmin,kmax)分別為((0,-∞),(+∞,0))，((-∞,0),(0,-∞))，((+∞,0),(0,+∞))和((0,+∞),(-∞,0)).Re{·}表征取實數(shù)部算子.顯然，公式(3)也可以通過空間域Hilbert變換實現(xiàn).需要特別指出的是，盡管我們在論述的過程中選用的是標量聲波方程，使用一階速度-應力方程同樣可以先對震源進行Hilbert變換后延拓得到伴隨波場，再實現(xiàn)波場分解.

1.2 針對不同地層結(jié)構的成像條件

逆時偏移作為一種基于雙程波動方程延拓的方法，成像結(jié)果往往受低頻噪音以及偏移假象的影響.解決這一問題的途徑之一是在成像條件中只保留下行震源波場和上行檢波器波場的互相關(Claerbout,1971;Fei et al.,2015)，即

(4)

其中，s和r分別表示震源波場和檢波器波場，下標u和d則分別代表上行和下行方向，tmax是最大記錄時間.按照波的傳播方向，震源波場和檢波器波場除了可以分為上下組分之外，還可以進一步分為左右組分,即sd=sld+srd,ru=rlu+rru.因此，(4)式又可以進一步寫為

(5)

利用震源波場和檢波器波場不同組分之間的互相關，可以實現(xiàn)對單一傾向的地質(zhì)構造的成像.圖2a是一正傾角地層，地震波從震源激發(fā)向右下方向傳播至界面，沿左上方向發(fā)生反射并被地表檢波器接收.與之對應，對于負傾角地層，地震波是向左下方傳播而反射波則沿右上方向傳播至檢波器，如圖2b所示.不同于傾斜的連續(xù)反射界面，繞射點可以同時滿足兩種情況下的入射波和反射波的方向特征.換言之，如果我們利用(5)式的第一項I1(x)，則只有正傾角地層和繞射點被成像；如果利用(5)式的第二項I2(x)，則只有負傾角地層和繞射點被成像.進一步地，針對繞射點定位的乘法成像條件被提出(Zhang et al.,2019a,b)，即

圖2 針對正傾角(a)和負傾角(b)的地層構造，入射波與出射波之間的關系Fig.2 The relationship between the incident wave and the outgoing wave for the formation structure with positive dip angle (a)and negative dip angle (b)

(6)

在這里，ns表示總炮數(shù).

不同于前兩者，采用I3(x)和I4(x)得到的構造剖面沒有特定性，無論正、負傾角地層、平界面或者繞射點，都可能被成像.并且由互易性定理可知，當觀測系統(tǒng)對地下介質(zhì)能夠充分覆蓋時，I3(x)和I4(x)的成像結(jié)果具有很好的相似性.

為了更好地利用多次波、回轉(zhuǎn)波、棱柱波等波型，我們還可以采用Causal成像條件(Liu et al.,2011;Li et al.,2019)，即

(7)

與之對應，有

(8)

不同于(5)式，上式中的成像條件可以在不顯式對波場方向分解的情況下實現(xiàn).經(jīng)過推導(詳見附錄)，(8)式可以寫為

(9)

其中

(10)

(11)

1.3 高質(zhì)量傾角道集的提取

在地震資料處理中，角度域共成像點道集具有十分重要的作用.在張角域，反射波和繞射波信號之間沒有明顯的區(qū)別，因此區(qū)分二者能量比較困難.但是在傾角域，反射波能量只集中在地層的真實傾角附近，而繞射波的能量則連續(xù)分布在一定的范圍內(nèi).所以，通過先提取傾角道集再拾取繞射波能量的方式，可以凸顯繞射點的正確位置.Jin等(2014)將基于逆時偏移計算角道集的方法總結(jié)為三類：方向矢量法(DVB)，局部平面波分解法(LPWD)和局部相移成像條件法(LSIC).其中，方向矢量法運算量少，方便快捷，并且結(jié)果分辨率高.為此，我們采用Poynting矢量法計算共成像點傾角道集.

本文選擇聲波方程的一階速度-應力形式對波場進行延拓，因此Poynting矢量可以直接由速度和應力之積得到，即

p=-τv,

(12)

在這里，p表示Poynting矢量，τ表示應力張量.v表示質(zhì)點的運動速度矢量.Poynting矢量的幅值代表能流密度的大小，其方向指向波的傳播方向.由(12)式可知，對于任意一個空間位置，Poynting矢量只能給出一個傳播方向.所以在波發(fā)生重疊的區(qū)域，例如強反射界面附近，計算得到的角度準確性不足.

利用Hilbert變換伴隨波場，首先將原始聲壓場和質(zhì)點速度分量分離為沿不同方向(左上，左下，右上和右下)的組分，再使用(12)式分別獲取不同波精確的傳播方向.雖然波場分解增加了方案的復雜度和運算量，但能夠避免由波重疊引起的Poynting矢量不準確的問題，從而提升傾角道集的質(zhì)量.

圖3描述了波由震源出發(fā)，在界面處發(fā)生反射最終被檢波器接收的過程.假設在成像點處發(fā)生反射，利用震源波場Poynting矢量ps和檢波器波場Poynting矢量pg可以計算出地層傾角α，即

圖3 傾角計算的幾何示意Fig.3 Geometric diagram of the dip-angle calculation

(13)

在偏移過程中，將每個時刻震源波場和檢波器波場的互相關值疊加到成像點對應的α處即可得到共成像點傾角道集.

在傾角道集中，甄別、拾取繞射點能量既可以人工操作，也可以利用中值濾波等方法來自動實現(xiàn).值得注意的是，繞射能量也會受到偏移孔徑的影響，對于有限的空間采樣，繞射點埋深越大，其在傾角域分布的范圍越窄.而窄化的繞射點能量會弱化和反射點之間的區(qū)別，增大繞射點成像的難度.

1.4 成像后繞射點定位

前兩種方案都是以逆時偏移為基礎，根據(jù)地震波遇到不同地質(zhì)構造時射線路徑的特征，實現(xiàn)繞射點定位和能量聚集的.因此，它們的過程相對復雜，存儲成本較高，運算量成倍增加.如果能夠直接在成像剖面上拾取繞射點和間斷點，方法將不局限于逆時偏移技術，且效率也會有效提升.目前，已有相關研究展開.例如，奚先和黃江清(2020)通過卷積神經(jīng)網(wǎng)絡(CNN)實現(xiàn)了在地震剖面上散射體的定位和成像.本文提出另一種簡單直接的繞射點定位方案.

(14)

與式(6)相同，我們采用乘法成像條件獲取單獨的繞射點成像剖面，即

Idif(x)=I+(x)·I-(x).

(15)

顯然，本方案與方案一相比所需的運算量較少、可以對工區(qū)剖面做局部化運算、也不局限于逆時偏移所產(chǎn)生的地震剖面.

這里我們給出一個簡單的實例.圖4a是輸入圖像，包括左傾、右傾、水平和豎直四條線段和一個孤點.圖4b是與4a相對應的振幅譜，按照kx·kz≥0和kx·kz<0，我們分別得到了圖4c和4d.進一步，經(jīng)過反傅里葉變換，我們得到了如圖4e和4f的構造.需要注意的是，在圖4e中，除了傾角為正值的線段和孤點，還包括連續(xù)的水平線段和垂直線段，這是因為我們在其波數(shù)域中保留了kx·kz=0的區(qū)域.對于實際的地震剖面，為了避免在較平緩界面或者接近垂直的構造處成像，我們會在波數(shù)域kx和kz接近0的位置將值設置為零.最后采用(15)式乘法運算，我們得到了孤立點以及唯一的交叉點的位置，如圖4g中紅圈所示.

圖4 原始輸入圖像，不同構造分離和繞射點成像圖(a)是原始圖像，包括不同傾向的線段和繞射點.圖(b)是與圖(a)相對應的幅度譜，按照(14)式可分為圖(c)和(d).經(jīng)過反傅里葉變換，可得圖(e)和(f).進一步，按照乘法成像條件可得同時包含交叉點和孤點的最終剖面(g).Fig.4 Original input image,separation of different structures,and imaging of diffraction pointsPanel (a)is the original image,including line segments with different inclinations and an isolated point.Panel (b)is the amplitude spectrum corresponding to figure (a),which can be divided into panels (c)and (d)according to equation (14).After the inverse Fourier transform,panels (e)and (f)can be obtained.Furthermore,by using the multiplication imaging condition,we have a final profile (g)containing both an intersection point and an isolated point.

2 數(shù)值算例

為了驗證上述理論的正確性和方案的可行性，我們在簡單模型和Marmousi2模型上分別進行了實驗并對結(jié)果進行了展示.

2.1 簡單模型

在此模型中，我們設置了正、負傾角反射面、平界面以及三個繞射點，如圖5所示.40個炮點在地表均勻分布，間隔設為75 m，檢波器則分布在整個地表，間隔設為5 m.震源選擇為主頻30 Hz，時間采樣間隔0.5 ms的雷克子波.地震數(shù)據(jù)由時間二階、空間八階的有限差分算法模擬得到.該實驗選擇3000 m·s-1的均勻模型作為偏移背景速度.

圖5 簡單模型，其中包括正負傾角界面，平界面以及三個繞射點Fig.5 A simple model.It includes positive-and negative-dip interfaces,a flat interface,and three diffraction points

在震源波場和檢波器波場中分別選取不同的組分采用互相關成像條件，可以得到不同的成像結(jié)果.圖6a和6b分別對應式(5)中的I2(x)和I1(x)，即分別對負和正傾角構造進行了成像.除了連續(xù)界面，繞射點或者不連續(xù)點在剖面上也有所顯示.通過I1(x)和I2(x)剖面相乘，可以獲取只包含繞射點位置的成像結(jié)果，即圖6f.另外，我們在圖6c和6d中展示了I3(x)和I4(x)所對應的剖面，由于震源點和檢波器點之間的互易性，兩者有一定的相似性.特別地，我們計算了以(4)式作為成像條件的剖面Isum(x)，如圖6e所示.按照第三種方案，我們采用獲取到的完整剖面作為輸入，分別提取了負傾角構造和正傾角構造，如圖7a和7b所示.應當注意的是，我們在波數(shù)域的分離避開了波數(shù)分量為0的區(qū)域，從而避免了單一傾角構造剖面上的平界面的出現(xiàn).最后，通過乘法成像條件，獲取了圖7c中的繞射點剖面.從圖中可見，方案一和方案三在簡單的模型上的繞射點定位結(jié)果基本一致.圖8a和8b是按照成像條件(11)得到的負傾角和正傾角剖面，在求取過程中需對震源波場和檢波器波場進行Hilbert變換，而無需伴隨波場的延拓和原始波場的方向分解.圖8c是由圖8a和8b得到的繞射點剖面.

圖6 簡單模型不同成像條件下的成像結(jié)果(a)是負傾角構造，與成像條件I2(x)對應；(b)是正傾角構造，與成像條件I1(x)對應；(c)和(d)分別是在成像條件I3(x)和I4(x)下的剖面；(e)完整的成像剖面；(f)繞射點剖面,為(a)和(b)的乘積.Fig.6 Imaging results of the simple model under different imaging conditions(a)is the negative-dip structure,which corresponds to the imaging condition I2(x);(b)is the positive-dip structure,which corresponds to the imaging condition I1(x);(c)and (d)are the profiles under imaging conditions I3(x)and I4(x),respectively;(e)is the complete profile;(f)is the profile of diffraction points,which is the product of (a)and (b).

圖7 簡單模型按照方案三得到的成像剖面(a)負傾角構造；(b)正傾角構造；(c)繞射點剖面,為(a)和(b)的乘積.Fig.7 Imaging profiles of the simple model obtained from scheme 3(a)The negative-dip structure;(b)The positive-dip structure;(c)The profile of diffraction points,which is the product of (a)and (b).

圖8 簡單模型按照式(11)中的成像條件得到的成像剖面(a)負傾角構造；(b)正傾角構造；(c)繞射點剖面,為(a)和(b)的乘積.Fig.8 Imaging profiles of the simple model according to image conditions in equation (11)(a)The negative-dip structure;(b)The positive-dip structure;(c)The profile of diffraction points,which is the product of (a)and (b).

在計算疊加剖面的同時，我們提取了共成像點傾角道集，并在圖9中展示了分別與I1、I2、I3、I4和Isum所對應的角道集在橫向1500 m處的片段.從圖中可以看出，平界面在以I1和I2為成像條件提取的道集上沒有能量聚集，但是繞射點在四種成像條件下均可被呈現(xiàn).在疊加后的角道集9e上，繞射波能量連續(xù)分布在一段區(qū)域內(nèi)，而反射波的能量則聚焦在0°附近.以這些特征為依據(jù)，我們在傾角域手動拾取了繞射波能量.圖10展示了拾取的結(jié)果，繞射點和界面不連續(xù)點都得到了很好的呈現(xiàn)，且保留了同相軸的振幅和相位特征.至此，三種繞射點成像方案在簡單模型上都得到了良好的效果.

圖9 簡單模型在不同成像條件下的傾角道集(a),(b),(c)和(d)分別與成像條件 I1,I2,I3和I4對應；(e)是(a),(b),(c)和(d)之和.Fig.9 Dip gathers for the simple model under different imaging conditions(a),(b),(c),and (d)correspond to the imaging conditions I1,I2,I3,and I4,respectively.(e)is the sum of (a),(b),(c),and (d).

圖10 簡單模型在傾角域手動拾取能量得到的繞射點剖面Fig.10 Diffraction point profile of the simple model obtained by manually picking up energy in the dip-angle domain

2.2 Marmousi2模型

相對復雜的Marmousi2模型的速度變化如圖11a所示.該模型在水平方向上有1325個采樣點，豎直方向上有467個采樣點.震源埋深為5 m,總炮數(shù)為40，炮間距設置為165 m.而檢波器間隔設為5 m,在整個地表均有分布.地震數(shù)據(jù)同樣由時間二階、空間八階的有限差分算法模擬得到.偏移速度由真實模型經(jīng)過高斯光滑之后得到，如圖11b所示.

圖11 Marmousi2模型真實速度(a)與偏移所用的背景速度(b)Fig.11 True velocity of the Marmousi2 model (a)and background velocity used for migration (b)

利用震源波場的左下組分和檢波器波場的右上組分之間的互相關運算I2，我們得到了如圖12a所示的負傾角構造；與之相反，利用震源波場的右下組分和檢波器波場的左上組分之間的互相關運算I1，我們得到了如圖12b所示的正傾角構造.然后，借助正負傾角剖面的乘積，實現(xiàn)了繞射點和間斷點的定位，如圖12c所示.同樣，我們計算了采用I3和I4成像條件時的結(jié)果，即13a和13b，二者具有很好的相似性.

圖12 Marmousi2模型的負傾角剖面(a)和正傾角剖面(b)，以及最終的繞射點剖面(c)Fig.12 The negative-and positive-dip profiles of the Marmousi2 model,and the final diffraction point profile

圖13 Marmousi2模型在成像條件I3(a)和I4(b)下的結(jié)果，以及完整的剖面(c)Fig.13 (a)and (b)are the imaging profiles of the Marmousi2 model under imaging conditions I3 and I4,respectively.(c)is the complete profile

圖14 Marmousi2模型按照方案三得到的成像剖面(a)負傾角構造；(b)正傾角構造；(c)繞射點剖面,為(a)和(b)的乘積.Fig.14 Imaging profiles of the Marmousi2 model obtained from scheme 3(a)The negative-dip structure;(b)The positive-dip structure;(c)The profile of diffraction points,which is the product of (a)and (b).

圖15 Marmousi2模型按照式(11)中的成像條件得到的成像剖面(a)負傾角構造；(b)正傾角構造；(c)繞射點剖面,為(a)和(b)的乘積.Fig.15 Imaging profiles of the Marmousi2 model according to imaging conditions in equation (11)(a)The negative-dip structure;(b)The positive-dip structure;(c)The profile of diffraction points,which is the product of (a) and (b).

在圖16中，我們分別展示了與I1，I2，I3和I4相對應的傾角道集的片段，其位于地表2500 m處.相較于簡單模型，Marmmousi2模型提取到的角道集噪音更加發(fā)育，繞射點的識別特征也不夠強.在此，我們采用中值濾波對繞射波信號進行了甄別和判定，得到了如圖17的結(jié)果.當繞射點埋深較大且觀測系統(tǒng)孔徑有限時，繞射波信號與反射波信號并無顯著區(qū)別，為能量拾取帶來了一定的困難.因此，與圖12c相比，成像結(jié)果除在正確位置附近(如斷層處)比較清晰和準確外，在其他區(qū)域也有部分能量殘余.

圖16 Marmousi2模型按照不同成像條件得到的傾角道集(a),(b),(c)和(d)分別與成像條件 I1,I2,I3和I4對應；(e)是(a),(b),(c)和(d)之和.Fig.16 Dip-angle gathers for the Marmousi2 model under different imaging conditions(a),(b),(c),and (d)correspond to the imaging conditions I1,I2,I3,and I4,respectively.(e)is the sum of (a),(b),(c),and (d).

圖17 Marmousi2模型在傾角道集通過中值濾波得到的繞射點剖面Fig.17 Diffraction point profile of the Marmousi2 model is obtained by median filtering the dip-angle gathers

3 結(jié)論

本文以繞射點為成像目標，分析了基于逆時偏移的各種繞射波成像方案.既利用了逆時偏移方法對地下構造歸位準確的優(yōu)點，又盡量減少計算成本.通過理論分析和數(shù)值算例驗證，我們將三種方案總結(jié)如下：

(1)借助震源波場和檢波器波場的不同組分之間的互相關，可以實現(xiàn)不同構造的成像.利用繞射點在正傾角剖面和負傾角剖面均能呈現(xiàn)的特點，通過乘法條件可以實現(xiàn)繞射點的定位.這一方法經(jīng)過與Hilbert變換和伴隨波場延拓相結(jié)合，存儲成本大大降低，計算效率得以提高；

(2)從Causal成像條件出發(fā)，推導正負傾角構造成像公式，從而可以利用空間域的Hilbert變換避免震源波場和檢波器波場的方向分解；

(3)Poynting矢量法提取傾角道集進而拾取繞射能量.為了提升角道集的質(zhì)量，我們采用先波場分離再計算傾角的策略，避免了Poynting矢量在波場重疊區(qū)域求取不準的問題.這一方法的運算量和存儲量都較大，并且當偏移孔徑有限時，繞射波信號和反射波信號之間的區(qū)別會被弱化，從而在拾取上造成困難；

(4)通過地震剖面在波數(shù)域的分解，實現(xiàn)正負構造剖面的提取.然后采取乘法運算，確定繞射點和間斷點的位置.其中，為了避免在比較平緩或接近垂直界面上的錯誤成像，我們犧牲少許繞射點信息，將波數(shù)分量接近于0的區(qū)域值設置為0.相較于前兩種方案，這種圖像處理的方法具有以下優(yōu)點：運算和存儲成本極低，不受偏移方法的限制，針對大工區(qū)可以只對目標區(qū)域做局部化成像.當然，本方法也需要信噪比較高的剖面作為輸入.

附錄A

由正文可知

(A1)

在頻率-波數(shù)域，上式可表示為

(A2)

因為

(A3)

在此，我們定義中間變量

(A4)

顯然，有

(A5)

同理，有

(A6)

(A7)

(A8)

將上式變換回空間-時間域，可得

(A9)

其中，

(A10)

上標“±”表示取“+”或者取“-”.顯然，我們可以通過對空間進行Hilbert變換實現(xiàn)式(A9),即

(A11)

其中，H算子定義為Hx(s)=s+iHx(s)，Hz(s)=s+iHz(s),Re{·}表示取實數(shù)部運算.

與上述推導類似，我們得出

(A12)

這與Zhang等(2020)的結(jié)論一致.