陳生昌,魯方正,劉亞楠,周華敏

(1.浙江大學(xué)地球科學(xué)學(xué)院,浙江杭州310027;2.長(zhǎng)江水利委員會(huì)長(zhǎng)江科學(xué)院,湖北武漢430010)

一次波(一次反射波、一次散射波)地震數(shù)據(jù)的波動(dòng)方程偏移成像起源于20世紀(jì)70年代初期[1-3],是當(dāng)前地震數(shù)據(jù)處理中最重要的方法技術(shù)之一,也是獲得地下三維構(gòu)造圖像的主要技術(shù)手段[4-7]。然而,隨著油氣勘探開(kāi)發(fā)目標(biāo)的日趨復(fù)雜化和精細(xì)化(如復(fù)雜構(gòu)造巖性油氣藏和非常規(guī)油氣藏的勘探開(kāi)發(fā)),需要我們發(fā)展高分辨、高保真的反演成像方法。寬方位、寬頻帶、高密度數(shù)據(jù)采集方法技術(shù)的應(yīng)用,為我們利用波形(即振幅、相位)信息進(jìn)行高分辨、高保真的反演成像方法研究提供了數(shù)據(jù)基礎(chǔ)。地震數(shù)據(jù)全波形反演方法理論[8-10]的發(fā)展也為我們開(kāi)展高分辨、高保真的反演成像方法研究提供了借鑒。此外,地震數(shù)據(jù)的波動(dòng)方程偏移成像方法理論也需要從當(dāng)前的構(gòu)造成像走向以地層物性參數(shù)和儲(chǔ)層參數(shù)為目標(biāo)的巖性成像。

地震數(shù)據(jù)是地下介質(zhì)變化情況的客觀反映,不僅包含地震波的走時(shí)信息,還含有地震波的波形信息。地震勘探中地震波的頻帶具有特定的范圍,相應(yīng)的地震波長(zhǎng)也有一定的范圍,因此在利用地震數(shù)據(jù)研究地下介質(zhì)變化時(shí),須考慮地震波長(zhǎng)的長(zhǎng)短。地震數(shù)據(jù)的正演表達(dá)是構(gòu)建地震波反演成像方法的基礎(chǔ),不同的正演表達(dá)會(huì)導(dǎo)致不同的反演目標(biāo)和反演方法。從地震數(shù)據(jù)的全波形反演可知,要實(shí)現(xiàn)地震數(shù)據(jù)的巖性成像不僅需要利用好地震數(shù)據(jù)的走時(shí)信息,更要準(zhǔn)確利用地震數(shù)據(jù)的波形信息。為了在地震數(shù)據(jù)的偏移成像中準(zhǔn)確地利用波形信息,使當(dāng)前的構(gòu)造成像走向以地層物性參數(shù)和儲(chǔ)層參數(shù)為目標(biāo)的巖性成像,我們必須從波動(dòng)理論出發(fā)為地震數(shù)據(jù)構(gòu)建一個(gè)嚴(yán)謹(jǐn)?shù)木€性正演表達(dá)。

基于上述認(rèn)識(shí),我們從地震波的控制方程(波動(dòng)方程)出發(fā),利用擾動(dòng)理論,并根據(jù)地下非均勻體的尺寸或其物理性質(zhì)空間變化特征尺度與地震波長(zhǎng)之間的相對(duì)大小關(guān)系,將非均勻體劃分為產(chǎn)生散射波的散射體和產(chǎn)生反射波的反射體,推導(dǎo)出控制散射波的散射波方程和控制反射波的反射波方程,然后在一定的近似條件下分別得到散射地震數(shù)據(jù)和反射地震數(shù)據(jù)的線性正演表達(dá)式,再利用線性反演理論構(gòu)建基于散射地震數(shù)據(jù)和反射地震數(shù)據(jù)的線性正演表達(dá)式的地震數(shù)據(jù)波形成像方法理論。上述工作也是對(duì)我們近年相關(guān)研究工作[11-19]的補(bǔ)充、完善和系統(tǒng)化。

1 當(dāng)前地震數(shù)據(jù)波動(dòng)方程偏移的主要不足

根據(jù)CLAERBOUT[1]提出的地震數(shù)據(jù)波動(dòng)方程偏移成像理論和BERKHOUT[20-21]提出的有關(guān)地震波傳播的WRW概念模型,在給定滿足地震波運(yùn)動(dòng)學(xué)的準(zhǔn)確的光滑模型(也稱為偏移的宏觀模型)下,當(dāng)前的地震數(shù)據(jù)波動(dòng)方程逆時(shí)偏移主要由兩步組成:第一步是波場(chǎng)外推,即利用宏觀模型下的波動(dòng)方程分別對(duì)震源波場(chǎng)和記錄波場(chǎng)進(jìn)行順時(shí)外推和逆時(shí)外推;第二步是波場(chǎng)成像,即利用CLAERBOUT提出的時(shí)間一致性成像原理[22-26]建立的成像條件(公式)對(duì)外推得到的記錄波場(chǎng)和震源波場(chǎng)進(jìn)行成像,得到由地下反射系數(shù)估計(jì)的偏移成像結(jié)果,以此實(shí)現(xiàn)對(duì)反射面或散射(繞射)體空間位置的構(gòu)造成像。經(jīng)過(guò)近50年的發(fā)展,當(dāng)前的地震數(shù)據(jù)偏移成像方法主要包括:①利用描述波場(chǎng)傳播的Kirchhoff積分公式的射線Kirchhoff偏移方法和射線Beam偏移方法[27-32];②在Born近似和WKBJ近似建立的散射數(shù)據(jù)表達(dá)式或Kirchhoff近似建立的基于反射系數(shù)的反射數(shù)據(jù)表達(dá)式的基礎(chǔ)上,利用漸進(jìn)反演理論構(gòu)建的積分形式的射線Kirchhoff真振幅偏移方法[33-36];③利用單程波波動(dòng)方程和Claerbout成像條件的單程波偏移方法[37-40];④利用雙程波波動(dòng)方程和Claerbout成像條件的逆時(shí)偏移方法[41-43];⑤借助射線Kirchhoff真振幅偏移方法或利用采集照明補(bǔ)償?shù)牟▌?dòng)方程真振幅偏移方法[44-53];⑥利用Born近似的地震數(shù)據(jù)正演方程建立的波動(dòng)方程最小二乘偏移方法[54-55];⑦利用反射波方程構(gòu)建的波動(dòng)方程最小二乘偏移方法[13,17]。

由于利用了準(zhǔn)確的波動(dòng)方程進(jìn)行波場(chǎng)外推,逆時(shí)偏移被認(rèn)為是當(dāng)前理論上最為完善的偏移方法。對(duì)于標(biāo)量波方程的地震數(shù)據(jù)常規(guī)逆時(shí)偏移方法有以下計(jì)算公式。

1)光滑偏移速度模型下的地下入射波場(chǎng)計(jì)算。

(1)

2)光滑偏移速度模型下利用逆時(shí)外推重建地下反射波場(chǎng)。

(2)

3)利用成像條件(反射系數(shù)成像條件)進(jìn)行波場(chǎng)成像。

(3)

當(dāng)前地震數(shù)據(jù)波動(dòng)方程偏移成像方法存在的不足主要有以下3方面。①不區(qū)分散射數(shù)據(jù)和反射數(shù)據(jù),將散射和反射視為同義詞,用同一種公式對(duì)散射數(shù)據(jù)和反射數(shù)據(jù)進(jìn)行偏移成像[7,11]。我們知道,在波動(dòng)力學(xué)中,散射波與反射波具有不同的波場(chǎng)特征,散射波是由局部非均勻體形成的二次源產(chǎn)生的,而反射波是由空間上具有一定延續(xù)度的非均勻體形成的二次源產(chǎn)生的,可視為散射波的疊加,所以不能將基于散射概念建立的散射數(shù)據(jù)表達(dá)應(yīng)用于反射數(shù)據(jù)的偏移成像,同樣也不能將基于平面波平界面的反射概念的波場(chǎng)關(guān)系應(yīng)用于散射數(shù)據(jù)的偏移成像。②成像公式?jīng)]有正確考慮散射波、反射波的不同產(chǎn)生機(jī)制,公式(3)所表示的成像公式僅適合平面波小角度入射到無(wú)窮大平反射面的特殊情況。③缺乏與波動(dòng)方程偏移成像相適應(yīng)的地震數(shù)據(jù)嚴(yán)謹(jǐn)正演表達(dá)式,使得偏移成像不能準(zhǔn)確地利用地震數(shù)據(jù)的波形信息,得到的偏移成像結(jié)果會(huì)出現(xiàn)相位誤差。在當(dāng)前的偏移成像中,對(duì)于地震數(shù)據(jù)的正演表達(dá)有基于介質(zhì)模型參數(shù)擾動(dòng)的Born近似表達(dá)[33,56-57]和基于反射面反射率的Kirchhoff近似表達(dá)[58-60]。前者適合相對(duì)地震波長(zhǎng)為小尺度的局部非均勻體產(chǎn)生的散射數(shù)據(jù),而不適用于相對(duì)地震波長(zhǎng)為大尺度的非均勻體產(chǎn)生的反射數(shù)據(jù),后者主要用來(lái)表達(dá)反射面產(chǎn)生的反射數(shù)據(jù),但也存在問(wèn)題(見(jiàn)下文)?；诘卣饠?shù)據(jù)Born近似表達(dá)而建立的偏移成像方法得到的成像結(jié)果是有關(guān)模型參數(shù)擾動(dòng)的近似估計(jì),如果非均勻體為小尺度的局部散射體,則該近似估計(jì)可有效反映散射體的位置;如果非均勻體為大尺度的反射體,則該近似估計(jì)難以準(zhǔn)確反映反射體的邊界,因此基于模型參數(shù)擾動(dòng)Born近似的地震數(shù)據(jù)表達(dá)是不適合反射數(shù)據(jù)的構(gòu)造(反射面位置)成像。根據(jù)上述認(rèn)識(shí),我們認(rèn)為在偏移成像中對(duì)于地震散射數(shù)據(jù)和反射數(shù)據(jù)應(yīng)該有不同的正演表達(dá)公式,并由此建立其不同的偏移成像計(jì)算公式。

在應(yīng)力與應(yīng)變間滿足線性近似的假定下,地震震源激發(fā)的地震波在地下的運(yùn)動(dòng)可以用下述二階微分算子方程描述:

L(m)u=f

(4)

式中:L為與地下模型物性參數(shù)m(也稱為彈性參數(shù))有關(guān)的二階偏微分算子,也稱為波動(dòng)算子;u表示地震波在地下運(yùn)動(dòng)狀態(tài)的狀態(tài)變量,也稱為地震波場(chǎng);f為激發(fā)地震波的震源函數(shù)。方程(4)可描述震源和地下模型m變化產(chǎn)生的各種波,如直達(dá)波、散射波(繞射波)、反射波、透射波和轉(zhuǎn)換波(對(duì)于彈性波動(dòng)方程)等等,也是標(biāo)量波方程、聲波方程和彈性波方程等的抽象形式,因此被稱為一般(或通用)波動(dòng)方程,是地震波場(chǎng)模擬、地震數(shù)據(jù)全波形反演的基礎(chǔ)方程。

給定滿足地震波運(yùn)動(dòng)學(xué)的準(zhǔn)確的光滑偏移模型mb(也稱為用于偏移的宏觀模型),逆時(shí)偏移中用于波場(chǎng)傳播的通用波動(dòng)方程為:

L(mb)ui=f

(5)

由于mb的光滑性,方程(5)和其對(duì)應(yīng)的齊次方程只能描述地震波的傳播,而不能描述地震波的散射、反射和波型轉(zhuǎn)換,因此入射波場(chǎng)ui中不包含散射波場(chǎng)、反射波場(chǎng)和轉(zhuǎn)換波場(chǎng)。

逆時(shí)偏移中Claerbout成像公式(條件)(3)的反射波場(chǎng)除以入射波場(chǎng)的成像公式是平面波小角度入射到無(wú)窮大平反射面的反射系數(shù)計(jì)算公式。對(duì)于無(wú)窮大平邊界和平面波,在入射角小于臨界角的條件下,反射系數(shù)為與頻率無(wú)關(guān)的實(shí)數(shù)[61-62],有:

ur=Rui

(6)

式中:R為與入射角有關(guān)的鏡面反射系數(shù)。如果入射波不是平面波或反射面不是無(wú)窮大平邊界或入射角大于臨界角,則公式(6)不成立。因?yàn)楣?6)中的反射系數(shù)R為實(shí)數(shù),則要求式中的反射波與入射波應(yīng)有相同的波形(即有相同的相位)。我們將(6)式所表示的反射波與入射波之間的關(guān)系稱為鏡面反射方程,它也被用于反射數(shù)據(jù)的Kirchhoff近似表達(dá)和基于Kirchhoff近似表達(dá)的偏移成像中[58-59]。因此,我們認(rèn)為當(dāng)前廣泛應(yīng)用的Claerbout成像條件和Kirchhoff近似對(duì)于非無(wú)窮大平邊界和非平面波地震數(shù)據(jù)的偏移成像存在不足,它無(wú)法考慮入射波與反射波、散射波之間的相位差異。

在當(dāng)前的偏移成像方法中,將方程(5)和方程(6)作為描述地震波傳播和反射(散射)的通用方程,并將它們聯(lián)合作為偏移成像方法中地震數(shù)據(jù)的正演表達(dá)。由上述的分析可知,方程(5)和方程(6)是不同條件下相互獨(dú)立的波動(dòng)方程,它們的聯(lián)合不能作為非無(wú)窮大平邊界和非平面波情況下地震數(shù)據(jù)的正演表達(dá)。因此,我們認(rèn)為當(dāng)前的偏移成像缺乏與之適應(yīng)的地震數(shù)據(jù)嚴(yán)謹(jǐn)正演表達(dá)公式,致使其不能準(zhǔn)確地利用地震數(shù)據(jù)的波形信息。也由于方程(6)的嚴(yán)格假定條件,導(dǎo)致當(dāng)前的偏移成像方法對(duì)于散射數(shù)據(jù)和反射數(shù)據(jù)的偏移成像通常都會(huì)出現(xiàn)相位誤差。

為了使偏移成像方法更適合地下實(shí)際情況和能準(zhǔn)確地利用地震數(shù)據(jù)波形信息,彌補(bǔ)當(dāng)前地震數(shù)據(jù)偏移成像方法理論中存在的不足和修正當(dāng)前波動(dòng)方程偏移成像方法對(duì)散射數(shù)據(jù)和反射數(shù)據(jù)偏移成像存在的相位誤差,有必要為偏移成像方法提供與之相適應(yīng)的地震數(shù)據(jù)嚴(yán)謹(jǐn)正演表達(dá),以滿足地下復(fù)雜構(gòu)造油氣藏、巖性油氣藏和非常規(guī)油氣藏勘探開(kāi)發(fā)對(duì)高分辨、高保真偏移成像的要求。本文在滿足地震波運(yùn)動(dòng)學(xué)的準(zhǔn)確的光滑模型下,首先將擾動(dòng)法應(yīng)用于通用波動(dòng)方程,得到描述非均勻體產(chǎn)生的擾動(dòng)波場(chǎng)的波動(dòng)方程。依據(jù)地震波長(zhǎng)與非均勻體大小或非均勻體物理性質(zhì)空間變化特征尺度之間的相對(duì)大小關(guān)系,將非均勻體劃分為在空間上具有有限延續(xù)度的產(chǎn)生散射波的散射體或空間上具有一定延續(xù)度的產(chǎn)生反射波的反射體,相應(yīng)的描述擾動(dòng)波場(chǎng)的波動(dòng)方程轉(zhuǎn)化為散射波動(dòng)方程和反射波動(dòng)方程。利用小擾動(dòng)假定、Born近似(或一次波近似)和高頻近似(地震波長(zhǎng)相對(duì)于非均勻體大小或其物理性質(zhì)空間變化特征尺度為小量),得到基于散射體物性參數(shù)相對(duì)擾動(dòng)的散射數(shù)據(jù)線性正演表達(dá)和基于反射體波阻抗相對(duì)擾動(dòng)的反射數(shù)據(jù)線性正演表達(dá)。為了描述反射體邊界對(duì)反射波的作用,在高頻近似條件下通過(guò)定義反射體波阻抗相對(duì)擾動(dòng)沿入射波傳播方向的方向?qū)?shù)為反射體邊界的局部反射系數(shù),推出基于局部反射系數(shù)的反射數(shù)據(jù)線性正演表達(dá)。最后給出標(biāo)量波、聲波和彈性波方程下散射數(shù)據(jù)和反射數(shù)據(jù)的具體線性正演表達(dá)式。

2 地震波的散射和反射及其通用方程

地震勘探中的一次地震波在地下傳播的物理過(guò)程直觀上可表述為:地表震源激發(fā)的入射波場(chǎng)向下傳播,入射波場(chǎng)與非均勻體作用形成二次源,并產(chǎn)生次生波場(chǎng)(散射波場(chǎng)、反射波場(chǎng)等),次生波場(chǎng)向上傳播到地表被檢波器接收。在上述的直觀表述中,地震波的傳播和散射(反射)是相互獨(dú)立的,但在描述地震波場(chǎng)運(yùn)動(dòng)狀態(tài)的通用波動(dòng)方程(4)中,地震波的傳播與散射(反射)是交織在一起的。這是因?yàn)榉匠?4)中的模型m不僅控制了地震波的傳播,也控制了地震波的散射(反射)。如果方程(4)中的模型m為光滑模型,則波動(dòng)方程(4)就不能描述地震波的散射、反射和波型轉(zhuǎn)換,即波動(dòng)方程(5)中的波場(chǎng)ui不包含散射、反射和轉(zhuǎn)換波。為了研究光滑模型下地震波的波場(chǎng)模擬、反演和偏移成像方法,如反射波形反演、反射波阻抗反演、反射數(shù)據(jù)的偏移成像(包括最小二乘偏移成像)方法等,我們首先需要建立光滑模型下描述地震波的散射和反射的波動(dòng)方程,并以此得到適合地震數(shù)據(jù)偏移成像的地震數(shù)據(jù)嚴(yán)謹(jǐn)正演表達(dá)公式。

對(duì)通用波動(dòng)方程(4)應(yīng)用擾動(dòng)法[63],即分別引入模型m和波場(chǎng)u的擾動(dòng)表達(dá):m=mb+δm,u=ui+δu,其中,mb為光滑的背景模型(簡(jiǎn)稱為光滑模型),δm為模型擾動(dòng)(非均勻體的表示),ui為光滑模型mb下的入射波場(chǎng)(也稱為背景波場(chǎng)),δu為模型擾動(dòng)δm產(chǎn)生的擾動(dòng)波場(chǎng)。將模型與波場(chǎng)的擾動(dòng)表達(dá)代入方程(4),得到與其對(duì)應(yīng)的擾動(dòng)方程,即非均勻體產(chǎn)生的擾動(dòng)波場(chǎng)所滿足的波動(dòng)方程。

L(mb)δu=-δL(mb,δm)u

(7)

式中:L(mb)為光滑模型下的波動(dòng)算子;δL(mb,δm)為擾動(dòng)算子,有δL(mb,δm)=L(m)-L(mb)=L(mb+δm)-L(mb)。

地震數(shù)據(jù)是一種具有特定頻帶范圍的波形數(shù)據(jù),因此利用地震數(shù)據(jù)研究地下非均勻體結(jié)構(gòu)必須考慮地震波長(zhǎng)與非均勻體δm的空間尺寸或者非均勻體δm物性空間變化的特征尺度之間的相對(duì)大小關(guān)系。如果非均勻體δm的尺寸a或者其物性空間變化特征尺度a相對(duì)于地震波場(chǎng)u的波長(zhǎng)λ比較小,通常假定a/λ≤1,則非均勻體δm相對(duì)于波長(zhǎng)在空間上只有有限的延續(xù)度,可視為產(chǎn)生散射波的局部散射體,相應(yīng)的波動(dòng)方程(7)化為散射波所滿足的波動(dòng)方程,即:

L(mb)us=-δLs(mb,δm)u

(8)

式中:us表示散射體產(chǎn)生的散射波場(chǎng);δLs(mb,δm)表示與擾動(dòng)δm或相對(duì)擾動(dòng)δm/mb有關(guān)的散射算子,它與波場(chǎng)u的相互作用是激發(fā)散射波的虛源。方程(8)也稱為描述散射波的通用方程。由于δm相對(duì)波長(zhǎng)比較小,因此方程(8)的右端源可視為點(diǎn)源,產(chǎn)生的散射波沒(méi)有特定的傳播方向。利用背景模型下波動(dòng)方程(5)所對(duì)應(yīng)的Green函數(shù),將方程(8)表示為時(shí)空域積分形式,有:

(9)

式中:Ω為散射波源的分布空間;g為光滑模型mb下波動(dòng)方程(5)的Green函數(shù);“*t”表示時(shí)間褶積。

如果非均勻體δm的尺寸a或者其物性空間變化特征尺度a相對(duì)于地震波場(chǎng)u的波長(zhǎng)λ比較大,通常假定a/λ>1,則非均勻體δm相對(duì)于波長(zhǎng)在空間上具有一定的延續(xù)度,可視為產(chǎn)生反射波的反射體,相應(yīng)的波動(dòng)方程(7)化為反射波所滿足的波動(dòng)方程,即:

L(mb)ur=-δLr(mb,Ir,σ)u

(10)

式中:ur表示反射體產(chǎn)生的反射波場(chǎng);δLr(mb,Ir,σ)表示與波阻抗相對(duì)擾動(dòng)Ir和角度σ有關(guān)的反射算子,它與波場(chǎng)u的相互作用是激發(fā)反射波的虛源;σ表示入射方向與反射方向之間的開(kāi)角,它是由波動(dòng)算子中物理參數(shù)的空間導(dǎo)數(shù)和波場(chǎng)的空間導(dǎo)數(shù)在高頻近似條件下的表達(dá)式引入的;波阻抗相對(duì)擾動(dòng)Ir不同于散射體的物性參數(shù)相對(duì)擾動(dòng)δm/mb,它不僅與mb,δm有關(guān),還與反射開(kāi)角σ有關(guān)。方程(10)也稱為描述反射波的通用方程。由于δm相對(duì)于波長(zhǎng)在空間上具有一定延續(xù)度,因此方程(10)的右端源可視為由點(diǎn)源集成的局部平面波源,產(chǎn)生的反射波可視為具有特定傳播方向的局部平面波。利用背景模型下波動(dòng)方程(5)所對(duì)應(yīng)的Green函數(shù),同樣可將方程(10)表示為積分形式,有:

(11)

式中:Ω為局部平面波源的分布空間。

我們根據(jù)非均勻體δm的尺寸a或者其物性空間變化特征尺度a與地震波長(zhǎng)λ之間的相對(duì)大小關(guān)系,將δm分別劃分為產(chǎn)生散射波的散射體和產(chǎn)生反射波的反射體,并分別得到了描述散射波和反射波的通用方程(公式(8)和公式(10)),以及它們的積分形式(公式(9)和公式(11))。

利用對(duì)δm的小擾動(dòng)假定和一階Born近似(忽略多次散射波),方程(8)退化為描述一次散射波us的非齊次波動(dòng)方程,即:

L(mb)us=-δLs(mb,δm)ui

(12)

相應(yīng)的方程(9)退化為:

(13)

利用對(duì)Ir的小值假定和一次反射波近似(忽略多次反射波),方程(10)退化為描述一次反射波ur的非齊次波動(dòng)方程,即:

L(mb)ur=-δLr(mb,Ir,σ)ui

(14)

相應(yīng)的方程(11)退化為:

(15)

如果反射體δm內(nèi)部的物性變化很小或不變化,則反射主要是由反射體邊界產(chǎn)生[64]。為了描述反射體邊界對(duì)反射的作用和滿足對(duì)反射體邊界成像的需要,本文將產(chǎn)生反射的波阻抗相對(duì)擾動(dòng)Ir沿入射波傳播方向I的方向?qū)?shù)定義為反射體邊界局部切平面的局部反射系數(shù),有:

(16)

式中:r表示局部反射系數(shù),是局部平面波入射到反射體邊界位置x處局部切平面上的局部反射系數(shù),它不僅可以表征反射面的空間位置,還能反映反射面上的物性變化。Ir在數(shù)學(xué)形態(tài)上是一個(gè)階躍函數(shù),而r在數(shù)學(xué)形態(tài)上是一個(gè)沿反射體邊界分布的δ函數(shù)。圖1為分界面上點(diǎn)x處局部切平面的反射示意圖。圖中Lp為點(diǎn)x處的局部切平面;r為點(diǎn)x處的局部反射系數(shù);ui和ur分別為具有局部平面波性質(zhì)的入射波和反射波;σ為點(diǎn)x處的反射開(kāi)角。

圖1 局部切平面的反射示意

假定地震波長(zhǎng)相對(duì)于光滑模型mb空間變化特征尺度滿足高頻近似條件(即光滑模型mb的空間變化相對(duì)于地震波長(zhǎng)尺度可視為緩慢變化甚至不變化),則公式(16)中的入射波傳播方向?qū)?shù)在波數(shù)域?qū)?yīng)為iki,ki為入射波傳播方向的波數(shù),有ki=ω/[vb(x)],也稱為入射波傳播方向在x處的局部波數(shù),ω為地震波的圓頻率,vb(x)為點(diǎn)x處光滑模型的速度。因此公式(16)在高頻近似條件下的波數(shù)域表達(dá)式為:

(17)

公式(16)和公式(17)定義的局部反射系數(shù)r是地下物性參數(shù)相對(duì)擾動(dòng)空間變化的反映。不同于公式(6)中與頻率無(wú)關(guān)的無(wú)量綱鏡面反射系數(shù)R,本文定義的局部反射系數(shù)r與傳播方向有關(guān)且具有長(zhǎng)度倒數(shù)的量綱,它是高頻近似條件下,局部平面入射波作用于局部反射平面的反射系數(shù),所以與入射波傳播方向上的局部波數(shù)有關(guān)(在形式上與頻率有關(guān))。

利用公式(16)和公式(17)定義的局部反射系數(shù)r,可將方程(14)和方程(15)分別改寫(xiě)為:

L(mb)ur=-δLr(mb,r,σ)ui

(18)

(19)

式中:δLr(mb,r,σ)表示與局部反射系數(shù)r有關(guān)的反射算子。方程(18)的右端源也是一個(gè)局部平面波源。

3 地震數(shù)據(jù)線性正演表達(dá)的具體形式

利用前文推導(dǎo)出的基于散射體物性參數(shù)擾動(dòng)的一次散射波通用方程(12)和方程(13)以及基于反射體波阻抗相對(duì)擾動(dòng)的一次反射波通用方程(14)和方程(15)或基于反射體邊界局部反射系數(shù)的一次反射波通用方程(18)和方程(19),結(jié)合具體的標(biāo)量波、聲波和各向同性彈性波方程,我們給出散射數(shù)據(jù)和反射數(shù)據(jù)線性正演表達(dá)的具體形式。

3.1 標(biāo)量波地震數(shù)據(jù)的線性正演表達(dá)

如果通用波動(dòng)方程(4)中模型物性參數(shù)m僅為速度v(x),則方程(4)為非齊次標(biāo)量波動(dòng)方程,即:

(20)

根據(jù)擾動(dòng)方法,有v(x)=vb(x)+δv(x),u(x,t;xs)=ub(x,t;xs)+δu(x,t;xs),其中,vb(x)為光滑速度模型,δv(x)為速度模型擾動(dòng),ub(x,t;xs)為光滑速度模型下的背景波場(chǎng),在下文也稱為入射波場(chǎng)ui(x,t;xs),δu(x,t;xs)為速度模型擾動(dòng)產(chǎn)生的擾動(dòng)波場(chǎng)。

假定速度擾動(dòng)體δv(x)的大小或其速度空間變化特征尺度相對(duì)于地震波長(zhǎng)比較小,即滿足上文定義的散射體條件,速度擾動(dòng)體δv(x)可視為散射體。根據(jù)通用的一次散射波方程(12)和方程(13),可得到下述標(biāo)量散射波方程:

(21)

(22)

(23)

方程(22)也稱為基于速度相對(duì)擾動(dòng)的標(biāo)量波地震散射數(shù)據(jù)的線性正演表達(dá)。

假定速度擾動(dòng)體δv(x)的大小或其速度空間變化特征尺度相對(duì)于地震波長(zhǎng)比較大,即滿足前文定義的反射體條件,速度擾動(dòng)體δv(x)可視為反射體。根據(jù)通用的一次反射波方程(14),可得到下述標(biāo)量反射波方程:

(24)

式中:ur(x,t;xs)表示由速度相對(duì)擾動(dòng)av(x)產(chǎn)生的一次反射波。

根據(jù)方程(15),由方程(24)可得到下述基于速度相對(duì)擾動(dòng)的積分形式標(biāo)量波地震反射波數(shù)據(jù)線性正演表達(dá),即:

(25)

由于標(biāo)量波動(dòng)算子中不含有對(duì)速度參數(shù)的空間導(dǎo)數(shù),其對(duì)應(yīng)的反射算子也就不含有角度σ,因此反射數(shù)據(jù)的表達(dá)式(25)與散射數(shù)據(jù)的表達(dá)式(22)有相同的形式,但其中av(x)的分布范圍Ω大小不同。

如果反射體δv(x)內(nèi)部的速度變化很小或不變化,則反射可視為由反射體邊界產(chǎn)生。根據(jù)反射體邊界局部切平面的局部反射系數(shù)定義式(17),對(duì)于速度相對(duì)擾動(dòng),有下述的局部反射系數(shù)表示式:

(26)

利用公式(26)和基于局部反射系數(shù)的通用一次反射波方程(18),可得到基于局部反射系數(shù)的標(biāo)量反射波方程:

(27)

根據(jù)方程(19),由方程(27)可以得到基于局部反射系數(shù)r的積分形式標(biāo)量反射波數(shù)據(jù)線性正演表達(dá)式:

(28)

基于Kirchhoff近似,BLEISTEIN等[36]推出了與公式(28)相類似的基于公式(6)中的鏡面反射系數(shù)R的反射數(shù)據(jù)表達(dá)式,由于公式(6)中鏡面反射系數(shù)R的嚴(yán)格條件,我們認(rèn)為BLEISTEIN等推導(dǎo)出的反射數(shù)據(jù)表達(dá)式對(duì)于非平面入射波和反射面為非無(wú)窮大平邊界存在局限性,不適合用于構(gòu)建實(shí)際地震數(shù)據(jù)的偏移成像。

3.2 聲波地震數(shù)據(jù)的線性正演表達(dá)

考慮到地下密度變化對(duì)地震波傳播的影響,假定通用波動(dòng)方程(4)中模型m為包含速度v(x)和密度ρ(x)的聲學(xué)介質(zhì),則方程(4)為非齊次聲波方程,即:

=f(t)δ(x-xs)

(29)

假定δv(x)和δρ(x)的大小或其空間變化特征尺度相對(duì)于地震波長(zhǎng)滿足散射體條件,即δv(x)和δρ(x)對(duì)應(yīng)的非均勻體可視為產(chǎn)生散射波的散射體。根據(jù)通用的一次散射波方程(12)和方程(13),可得到下述聲散射波方程:

(30)

(31)

ui(x,t;xs)=f(t)δ(x-xs)

(32)

對(duì)方程(31)應(yīng)用分步積分法,有:

(33)

方程(33)也稱為基于速度、密度相對(duì)擾動(dòng)的聲散射波數(shù)據(jù)的線性正演表達(dá)。

假定δv(x)和δρ(x)的大小或其空間變化特征尺度相對(duì)于地震波長(zhǎng)滿足反射體條件,即δv(x)和δρ(x)對(duì)應(yīng)的非均勻體可視為產(chǎn)生反射波的反射體。根據(jù)方程(15)和方程(33),可得到下述基于速度和密度相對(duì)擾動(dòng)的聲反射波方程,即:

(34)

式中:ur(xg,t;xs)為由av和aρ產(chǎn)生的一次反射波數(shù)據(jù)。

將方程(34)變換到頻率域,有:

(35)

(36)

將方程(36)變換到時(shí)間域,有:

(37)

方程(37)對(duì)應(yīng)的微分形式方程為:

(38)

公式(38)中的av+aρ(1+cosσ)是一個(gè)與開(kāi)角σ有關(guān)的聲波阻抗相對(duì)擾動(dòng)。令I(lǐng)r(x,σ)=av+aρ(1+cosσ),如果取σ=0,即沿反射面法向入射,則有Ir(x,σ=0)=av+2aρ=2δv(x)/vb(x)+2δρ(x)/ρb(x)≈2δI(x)/Ib(x),Ib(x)=vb(x)ρb(x),δI(x)=δ[v(x)ρ(x)]=ρ(x)δv(x)+v(x)δρ(x)。由此可知,小角度反射主要受阻抗變化的控制。如果取σ=π,則有Ir(x,σ=π)=av=2δv(x)/vb(x)。由此可知,大角度反射主要受速度變化的控制。由方程(37),我們可以得到基于反射體聲波阻抗相對(duì)擾動(dòng)的聲波反射數(shù)據(jù)線性正演表達(dá)公式,即

(39)

為了得到基于反射體邊界局部反射系數(shù)的反射波動(dòng)方程和反射數(shù)據(jù)線性正演表達(dá)公式,我們定義聲反射的局部反射系數(shù)為聲波阻抗相對(duì)擾動(dòng)Ir(x,σ)沿入射波傳播方向I的方向?qū)?shù)。根據(jù)反射體邊界局部切平面的局部反射系數(shù)定義式(17),對(duì)于聲波阻抗相對(duì)擾動(dòng),有下述的局部反射系數(shù)表達(dá)式:

cosσ)]

(40)

利用(40)式定義的局部反射系數(shù),由方程(38)在高頻近似條件下可得到基于反射體邊界局部反射系數(shù)的聲反射波方程:

(41)

根據(jù)方程(19),由方程(41)可以得到基于反射體邊界局部反射系數(shù)的積分形式聲波反射數(shù)據(jù)線性正演表達(dá)公式,即

(42)

如果在聲波方程中不考慮密度變化,則聲波方程(29)退化為標(biāo)量波方程(20),相應(yīng)的方程(33)、方程(39)和方程(42)分別退化為方程(22)、方程(25)和方程(28)。

3.3 彈性波地震數(shù)據(jù)的線性正演表達(dá)

為了更真實(shí)地描述地下地震波的運(yùn)動(dòng),假定方程(4)中地下模型為各向同性彈性介質(zhì),其模型m的物性參數(shù)為L(zhǎng)amé系數(shù)λ(x)、μ(x)和密度ρ(x),則方程(4)為非齊次彈性波動(dòng)方程,即:

Lu(x,t;xs)=f(t,xs)

(43)

式中:u(x,t;xs)為位移矢量波場(chǎng),在笛卡爾坐標(biāo)系下,有u(x,t;xs)=[ux(x,t;xs),uy(x,t;xs),uz(x,t;xs)],其中,ux(x,t;xs),uy(x,t;xs),uz(x,t;xs)分別代表u(x,t;xs)的x,y,z方向分量;f(t,xs)為矢量震源函數(shù);L是一個(gè)3×3偏微分算子:

(44)

假定擾動(dòng)量δα(x)、δβ(x)、δρ(x)的大小或其空間變化特征尺度相對(duì)于地震波長(zhǎng)滿足前文定義的散射體條件,即δα(x)、δβ(x)、δρ(x)對(duì)應(yīng)的非均勻體可視為產(chǎn)生散射波的散射體。根據(jù)通用散射波方程(12),可得到下述通用的彈性散射波方程:

Lbus(x,t;xs)=ΔLsui(x,t;xs)

(45)

式中:us(x,t;xs)為彈性一次散射波的位移矢量波場(chǎng);Lb為光滑模型下的彈性波傳播算子;ΔLs為彈性波散射算子,有:

(46)

Lbui(x,t;xs)=f(t,xs)

(47)

利用光滑模型彈性波方程(47)所對(duì)應(yīng)的矢量Green函數(shù),可將方程(45)中的彈性散射波數(shù)據(jù)寫(xiě)成與方程(13)對(duì)應(yīng)的積分形式[65]。方程(45)也稱為基于縱、橫波速度、密度相對(duì)擾動(dòng)的彈性散射波數(shù)據(jù)的線性正演方程。

假定δα(x),δβ(x),δρ(x)的大小或其空間變化特征尺度相對(duì)于地震波長(zhǎng)滿足前文定義的反射體條件,即δα(x),δβ(x),δρ(x)對(duì)應(yīng)的非均勻體可視為產(chǎn)生反射波的反射體。根據(jù)通用反射波方程(14),利用與前兩節(jié)類似的推導(dǎo)可得到彈性反射波方程,即:

Lbur(x,t;xs)=ΔLrui(x,t;xs)

(48)

式中:ur(x,t;xs)為彈性一次反射波的位移矢量波場(chǎng);ΔLr為彈性波反射算子。

(49)

(50)

(51)

(52)

(53)

(54)

(55)

在法向入射時(shí),σ=0,有:

(56)

(57)

(58)

此時(shí)不發(fā)生波型轉(zhuǎn)換。

(59)

(60)

(61)

(62)

(63)

(64)

將(59)式代入方程(49),有:

(65)

(66)

式中:g(xg,x,t)為光滑模型下彈性波方程(47)所對(duì)應(yīng)的矢量Green函數(shù)。

(67)

對(duì)于SH波入射SH波反射的局部反射系數(shù),有:

(68)

對(duì)于SV波入射SV波反射的局部反射系數(shù),有:

(69)

對(duì)于P波入射SV波反射的局部反射系數(shù),有:

(70)

對(duì)于SV波入射P波反射的局部反射系數(shù),有:

(71)

(72)

(73)

(74)

(75)

由上述推導(dǎo)出的標(biāo)量波、聲波和彈性波的反射數(shù)據(jù)線性正演表達(dá)公式可知,基于波阻抗相對(duì)擾動(dòng)的反射數(shù)據(jù)線性正演表達(dá)公式的右端包含入射波場(chǎng)的時(shí)間二階導(dǎo)數(shù)與波阻抗相對(duì)擾動(dòng)的乘積,而基于局部反射系數(shù)的反射數(shù)據(jù)線性正演表達(dá)公式的右端包含入射波場(chǎng)的時(shí)間一階導(dǎo)數(shù)與局部反射系數(shù)的乘積。不同于平面波無(wú)窮大平邊界的無(wú)量綱鏡面反射系數(shù),本文基于高頻近似和局部平面波定義的局部切平面反射系數(shù)是一個(gè)具有長(zhǎng)度倒數(shù)量綱的變量,是真實(shí)模型與光滑模型之間模型參數(shù)相對(duì)擾動(dòng)沿入射波傳播方向的方向?qū)?shù),它不僅依賴于光滑模型,也與入射波傳播方向有關(guān)。Zoeppritz方程不僅描述了反射波與入射波之間的關(guān)系,還定義了平面波無(wú)窮大平邊界鏡面反射系數(shù)與入射方向和地下模型物理參數(shù)之間的定量關(guān)系。本文利用反射體近似和高頻近似推導(dǎo)的基于反射體彈性波阻抗相對(duì)擾動(dòng)的P/S波型彈性反射波方程與Zoeppritz方程類似,不僅建立了反射波與入射波之間的數(shù)學(xué)物理關(guān)系,也定義了反射體彈性波阻抗相對(duì)擾動(dòng)與入射方向和地下模型物理參數(shù)之間的定量關(guān)系,可為進(jìn)一步的基于偏移成像得到的角度域共成像點(diǎn)道集的巖性分析與反演提供理論基礎(chǔ)。由于本文定義的局部反射系數(shù)與入射波傳播方向的局部波數(shù)有關(guān),因此在目前我們還不能給出傾斜入射情況下局部反射系數(shù)與地下模型物理參數(shù)和入射波傳播方向之間的具體函數(shù)關(guān)系。

4 結(jié)論

當(dāng)前地震數(shù)據(jù)波動(dòng)方程偏移成像方法不區(qū)分散射波和反射波,也缺乏與之相應(yīng)的嚴(yán)謹(jǐn)正演方程,致使偏移成像不能準(zhǔn)確地利用地震數(shù)據(jù)的波形信息,得到的偏移成像結(jié)果存在相位誤差。為了使波動(dòng)方程偏移成像方法能更好地適應(yīng)地下實(shí)際情況和準(zhǔn)確地利用波形信息,彌補(bǔ)當(dāng)前偏移成像方法理論的不足,修正地震數(shù)據(jù)波動(dòng)方程偏移成像結(jié)果的相位誤差,并為發(fā)展地層物性參數(shù)和儲(chǔ)層參數(shù)成像方法打下基礎(chǔ),本文在給定滿足地震波運(yùn)動(dòng)學(xué)的準(zhǔn)確的光滑模型基礎(chǔ)上,利用擾動(dòng)方法將波動(dòng)方程轉(zhuǎn)化為描述非均勻體產(chǎn)生的擾動(dòng)波場(chǎng)的波動(dòng)方程,根據(jù)地下非均勻體的大小或非均勻體物理性質(zhì)空間變化特征尺度與地震波長(zhǎng)之間的相對(duì)大小關(guān)系,提出將非均勻體劃分為相對(duì)于地震波長(zhǎng)在空間上具有有限延續(xù)度的產(chǎn)生散射波的散射體或相對(duì)于地震波長(zhǎng)在空間上具有一定延續(xù)度的產(chǎn)生反射波的反射體,推導(dǎo)出反射波波動(dòng)方程。為區(qū)別于散射波,本文認(rèn)為反射波是由點(diǎn)源集成的局部平面波源產(chǎn)生的具有方向性的局部平面波。對(duì)于散射波,利用小擾動(dòng)假定和一階Born近似(忽略多次散射波),可得到基于散射體物性參數(shù)相對(duì)擾動(dòng)的地震一次散射數(shù)據(jù)的線性正演表達(dá)。對(duì)于反射波,利用波阻抗相對(duì)擾動(dòng)的小值假定、一次反射波近似(忽略多次反射波)和高頻近似,可得到基于反射體波阻抗相對(duì)擾動(dòng)的地震一次反射數(shù)據(jù)的線性正演表達(dá)。反射體的波阻抗相對(duì)擾動(dòng)不同于散射體的物性參數(shù)相對(duì)擾動(dòng),它不僅與反射體的物性參數(shù)相對(duì)擾動(dòng)有關(guān),還與反射開(kāi)角有關(guān)?；诜瓷潴w彈性波阻抗相對(duì)擾動(dòng)的P/S波型彈性反射波方程,不僅建立了反射波與入射波之間的數(shù)學(xué)物理關(guān)系,也定義了反射體彈性波阻抗相對(duì)擾動(dòng)與入射方向和地下模型物理參數(shù)相對(duì)擾動(dòng)之間的定量關(guān)系。對(duì)于高頻近似條件下定義的反射體與反射波,本文將波阻抗相對(duì)擾動(dòng)沿入射波傳播方向的方向?qū)?shù)定義為反射體邊界局部切平面對(duì)于局部平面波的局部反射系數(shù),得到基于局部反射系數(shù)的地震一次反射數(shù)據(jù)的線性正演表達(dá)。不同于平面波入射到無(wú)窮大平界面的與頻率無(wú)關(guān)的無(wú)量綱鏡面反射系數(shù),本文定義的局部反射系數(shù)具有長(zhǎng)度倒數(shù)的量綱且與入射波傳播方向有關(guān)(即入射波的波數(shù)),推導(dǎo)出的反射數(shù)據(jù)表達(dá)式具有更廣的適應(yīng)性。本文推導(dǎo)出的地震數(shù)據(jù)線性正演表達(dá)(特別是基于反射體波阻抗相對(duì)擾動(dòng)的地震反射數(shù)據(jù)線性正演表達(dá)),是開(kāi)展地震數(shù)據(jù)波形成像和其它波形線性反演方法研究的基礎(chǔ)方程。