陳生昌,魯方正,劉亞楠,周華敏
(1.浙江大學(xué)地球科學(xué)學(xué)院,浙江杭州310027;2.長江水利委員會長江科學(xué)院,湖北武漢430010)
一次波(一次反射波、一次散射波)地震數(shù)據(jù)的波動方程偏移成像起源于20世紀(jì)70年代初期[1-3],是當(dāng)前地震數(shù)據(jù)處理中最重要的方法技術(shù)之一,也是獲得地下三維構(gòu)造圖像的主要技術(shù)手段[4-7]。然而,隨著油氣勘探開發(fā)目標(biāo)的日趨復(fù)雜化和精細(xì)化(如復(fù)雜構(gòu)造巖性油氣藏和非常規(guī)油氣藏的勘探開發(fā)),需要我們發(fā)展高分辨、高保真的反演成像方法。寬方位、寬頻帶、高密度數(shù)據(jù)采集方法技術(shù)的應(yīng)用,為我們利用波形(即振幅、相位)信息進行高分辨、高保真的反演成像方法研究提供了數(shù)據(jù)基礎(chǔ)。地震數(shù)據(jù)全波形反演方法理論[8-10]的發(fā)展也為我們開展高分辨、高保真的反演成像方法研究提供了借鑒。此外,地震數(shù)據(jù)的波動方程偏移成像方法理論也需要從當(dāng)前的構(gòu)造成像走向以地層物性參數(shù)和儲層參數(shù)為目標(biāo)的巖性成像。
地震數(shù)據(jù)是地下介質(zhì)變化情況的客觀反映,不僅包含地震波的走時信息,還含有地震波的波形信息。地震勘探中地震波的頻帶具有特定的范圍,相應(yīng)的地震波長也有一定的范圍,因此在利用地震數(shù)據(jù)研究地下介質(zhì)變化時,須考慮地震波長的長短。地震數(shù)據(jù)的正演表達是構(gòu)建地震波反演成像方法的基礎(chǔ),不同的正演表達會導(dǎo)致不同的反演目標(biāo)和反演方法。從地震數(shù)據(jù)的全波形反演可知,要實現(xiàn)地震數(shù)據(jù)的巖性成像不僅需要利用好地震數(shù)據(jù)的走時信息,更要準(zhǔn)確利用地震數(shù)據(jù)的波形信息。為了在地震數(shù)據(jù)的偏移成像中準(zhǔn)確地利用波形信息,使當(dāng)前的構(gòu)造成像走向以地層物性參數(shù)和儲層參數(shù)為目標(biāo)的巖性成像,我們必須從波動理論出發(fā)為地震數(shù)據(jù)構(gòu)建一個嚴(yán)謹(jǐn)?shù)木€性正演表達。
基于上述認(rèn)識,我們從地震波的控制方程(波動方程)出發(fā),利用擾動理論,并根據(jù)地下非均勻體的尺寸或其物理性質(zhì)空間變化特征尺度與地震波長之間的相對大小關(guān)系,將非均勻體劃分為產(chǎn)生散射波的散射體和產(chǎn)生反射波的反射體,推導(dǎo)出控制散射波的散射波方程和控制反射波的反射波方程,然后在一定的近似條件下分別得到散射地震數(shù)據(jù)和反射地震數(shù)據(jù)的線性正演表達式,再利用線性反演理論構(gòu)建基于散射地震數(shù)據(jù)和反射地震數(shù)據(jù)的線性正演表達式的地震數(shù)據(jù)波形成像方法理論。上述工作也是對我們近年相關(guān)研究工作[11-19]的補充、完善和系統(tǒng)化。
根據(jù)CLAERBOUT[1]提出的地震數(shù)據(jù)波動方程偏移成像理論和BERKHOUT[20-21]提出的有關(guān)地震波傳播的WRW概念模型,在給定滿足地震波運動學(xué)的準(zhǔn)確的光滑模型(也稱為偏移的宏觀模型)下,當(dāng)前的地震數(shù)據(jù)波動方程逆時偏移主要由兩步組成:第一步是波場外推,即利用宏觀模型下的波動方程分別對震源波場和記錄波場進行順時外推和逆時外推;第二步是波場成像,即利用CLAERBOUT提出的時間一致性成像原理[22-26]建立的成像條件(公式)對外推得到的記錄波場和震源波場進行成像,得到由地下反射系數(shù)估計的偏移成像結(jié)果,以此實現(xiàn)對反射面或散射(繞射)體空間位置的構(gòu)造成像。經(jīng)過近50年的發(fā)展,當(dāng)前的地震數(shù)據(jù)偏移成像方法主要包括:①利用描述波場傳播的Kirchhoff積分公式的射線Kirchhoff偏移方法和射線Beam偏移方法[27-32];②在Born近似和WKBJ近似建立的散射數(shù)據(jù)表達式或Kirchhoff近似建立的基于反射系數(shù)的反射數(shù)據(jù)表達式的基礎(chǔ)上,利用漸進反演理論構(gòu)建的積分形式的射線Kirchhoff真振幅偏移方法[33-36];③利用單程波波動方程和Claerbout成像條件的單程波偏移方法[37-40];④利用雙程波波動方程和Claerbout成像條件的逆時偏移方法[41-43];⑤借助射線Kirchhoff真振幅偏移方法或利用采集照明補償?shù)牟▌臃匠陶嬲穹品椒╗44-53];⑥利用Born近似的地震數(shù)據(jù)正演方程建立的波動方程最小二乘偏移方法[54-55];⑦利用反射波方程構(gòu)建的波動方程最小二乘偏移方法[13,17]。
由于利用了準(zhǔn)確的波動方程進行波場外推,逆時偏移被認(rèn)為是當(dāng)前理論上最為完善的偏移方法。對于標(biāo)量波方程的地震數(shù)據(jù)常規(guī)逆時偏移方法有以下計算公式。
1)光滑偏移速度模型下的地下入射波場計算。
(1)
2)光滑偏移速度模型下利用逆時外推重建地下反射波場。
(2)
3)利用成像條件(反射系數(shù)成像條件)進行波場成像。
(3)
當(dāng)前地震數(shù)據(jù)波動方程偏移成像方法存在的不足主要有以下3方面。①不區(qū)分散射數(shù)據(jù)和反射數(shù)據(jù),將散射和反射視為同義詞,用同一種公式對散射數(shù)據(jù)和反射數(shù)據(jù)進行偏移成像[7,11]。我們知道,在波動力學(xué)中,散射波與反射波具有不同的波場特征,散射波是由局部非均勻體形成的二次源產(chǎn)生的,而反射波是由空間上具有一定延續(xù)度的非均勻體形成的二次源產(chǎn)生的,可視為散射波的疊加,所以不能將基于散射概念建立的散射數(shù)據(jù)表達應(yīng)用于反射數(shù)據(jù)的偏移成像,同樣也不能將基于平面波平界面的反射概念的波場關(guān)系應(yīng)用于散射數(shù)據(jù)的偏移成像。②成像公式?jīng)]有正確考慮散射波、反射波的不同產(chǎn)生機制,公式(3)所表示的成像公式僅適合平面波小角度入射到無窮大平反射面的特殊情況。③缺乏與波動方程偏移成像相適應(yīng)的地震數(shù)據(jù)嚴(yán)謹(jǐn)正演表達式,使得偏移成像不能準(zhǔn)確地利用地震數(shù)據(jù)的波形信息,得到的偏移成像結(jié)果會出現(xiàn)相位誤差。在當(dāng)前的偏移成像中,對于地震數(shù)據(jù)的正演表達有基于介質(zhì)模型參數(shù)擾動的Born近似表達[33,56-57]和基于反射面反射率的Kirchhoff近似表達[58-60]。前者適合相對地震波長為小尺度的局部非均勻體產(chǎn)生的散射數(shù)據(jù),而不適用于相對地震波長為大尺度的非均勻體產(chǎn)生的反射數(shù)據(jù),后者主要用來表達反射面產(chǎn)生的反射數(shù)據(jù),但也存在問題(見下文)?;诘卣饠?shù)據(jù)Born近似表達而建立的偏移成像方法得到的成像結(jié)果是有關(guān)模型參數(shù)擾動的近似估計,如果非均勻體為小尺度的局部散射體,則該近似估計可有效反映散射體的位置;如果非均勻體為大尺度的反射體,則該近似估計難以準(zhǔn)確反映反射體的邊界,因此基于模型參數(shù)擾動Born近似的地震數(shù)據(jù)表達是不適合反射數(shù)據(jù)的構(gòu)造(反射面位置)成像。根據(jù)上述認(rèn)識,我們認(rèn)為在偏移成像中對于地震散射數(shù)據(jù)和反射數(shù)據(jù)應(yīng)該有不同的正演表達公式,并由此建立其不同的偏移成像計算公式。
在應(yīng)力與應(yīng)變間滿足線性近似的假定下,地震震源激發(fā)的地震波在地下的運動可以用下述二階微分算子方程描述:
L(m)u=f
(4)
式中:L為與地下模型物性參數(shù)m(也稱為彈性參數(shù))有關(guān)的二階偏微分算子,也稱為波動算子;u表示地震波在地下運動狀態(tài)的狀態(tài)變量,也稱為地震波場;f為激發(fā)地震波的震源函數(shù)。方程(4)可描述震源和地下模型m變化產(chǎn)生的各種波,如直達波、散射波(繞射波)、反射波、透射波和轉(zhuǎn)換波(對于彈性波動方程)等等,也是標(biāo)量波方程、聲波方程和彈性波方程等的抽象形式,因此被稱為一般(或通用)波動方程,是地震波場模擬、地震數(shù)據(jù)全波形反演的基礎(chǔ)方程。
給定滿足地震波運動學(xué)的準(zhǔn)確的光滑偏移模型mb(也稱為用于偏移的宏觀模型),逆時偏移中用于波場傳播的通用波動方程為:
L(mb)ui=f
(5)
由于mb的光滑性,方程(5)和其對應(yīng)的齊次方程只能描述地震波的傳播,而不能描述地震波的散射、反射和波型轉(zhuǎn)換,因此入射波場ui中不包含散射波場、反射波場和轉(zhuǎn)換波場。
逆時偏移中Claerbout成像公式(條件)(3)的反射波場除以入射波場的成像公式是平面波小角度入射到無窮大平反射面的反射系數(shù)計算公式。對于無窮大平邊界和平面波,在入射角小于臨界角的條件下,反射系數(shù)為與頻率無關(guān)的實數(shù)[61-62],有:
ur=Rui
(6)
式中:R為與入射角有關(guān)的鏡面反射系數(shù)。如果入射波不是平面波或反射面不是無窮大平邊界或入射角大于臨界角,則公式(6)不成立。因為公式(6)中的反射系數(shù)R為實數(shù),則要求式中的反射波與入射波應(yīng)有相同的波形(即有相同的相位)。我們將(6)式所表示的反射波與入射波之間的關(guān)系稱為鏡面反射方程,它也被用于反射數(shù)據(jù)的Kirchhoff近似表達和基于Kirchhoff近似表達的偏移成像中[58-59]。因此,我們認(rèn)為當(dāng)前廣泛應(yīng)用的Claerbout成像條件和Kirchhoff近似對于非無窮大平邊界和非平面波地震數(shù)據(jù)的偏移成像存在不足,它無法考慮入射波與反射波、散射波之間的相位差異。
在當(dāng)前的偏移成像方法中,將方程(5)和方程(6)作為描述地震波傳播和反射(散射)的通用方程,并將它們聯(lián)合作為偏移成像方法中地震數(shù)據(jù)的正演表達。由上述的分析可知,方程(5)和方程(6)是不同條件下相互獨立的波動方程,它們的聯(lián)合不能作為非無窮大平邊界和非平面波情況下地震數(shù)據(jù)的正演表達。因此,我們認(rèn)為當(dāng)前的偏移成像缺乏與之適應(yīng)的地震數(shù)據(jù)嚴(yán)謹(jǐn)正演表達公式,致使其不能準(zhǔn)確地利用地震數(shù)據(jù)的波形信息。也由于方程(6)的嚴(yán)格假定條件,導(dǎo)致當(dāng)前的偏移成像方法對于散射數(shù)據(jù)和反射數(shù)據(jù)的偏移成像通常都會出現(xiàn)相位誤差。
為了使偏移成像方法更適合地下實際情況和能準(zhǔn)確地利用地震數(shù)據(jù)波形信息,彌補當(dāng)前地震數(shù)據(jù)偏移成像方法理論中存在的不足和修正當(dāng)前波動方程偏移成像方法對散射數(shù)據(jù)和反射數(shù)據(jù)偏移成像存在的相位誤差,有必要為偏移成像方法提供與之相適應(yīng)的地震數(shù)據(jù)嚴(yán)謹(jǐn)正演表達,以滿足地下復(fù)雜構(gòu)造油氣藏、巖性油氣藏和非常規(guī)油氣藏勘探開發(fā)對高分辨、高保真偏移成像的要求。本文在滿足地震波運動學(xué)的準(zhǔn)確的光滑模型下,首先將擾動法應(yīng)用于通用波動方程,得到描述非均勻體產(chǎn)生的擾動波場的波動方程。依據(jù)地震波長與非均勻體大小或非均勻體物理性質(zhì)空間變化特征尺度之間的相對大小關(guān)系,將非均勻體劃分為在空間上具有有限延續(xù)度的產(chǎn)生散射波的散射體或空間上具有一定延續(xù)度的產(chǎn)生反射波的反射體,相應(yīng)的描述擾動波場的波動方程轉(zhuǎn)化為散射波動方程和反射波動方程。利用小擾動假定、Born近似(或一次波近似)和高頻近似(地震波長相對于非均勻體大小或其物理性質(zhì)空間變化特征尺度為小量),得到基于散射體物性參數(shù)相對擾動的散射數(shù)據(jù)線性正演表達和基于反射體波阻抗相對擾動的反射數(shù)據(jù)線性正演表達。為了描述反射體邊界對反射波的作用,在高頻近似條件下通過定義反射體波阻抗相對擾動沿入射波傳播方向的方向?qū)?shù)為反射體邊界的局部反射系數(shù),推出基于局部反射系數(shù)的反射數(shù)據(jù)線性正演表達。最后給出標(biāo)量波、聲波和彈性波方程下散射數(shù)據(jù)和反射數(shù)據(jù)的具體線性正演表達式。
地震勘探中的一次地震波在地下傳播的物理過程直觀上可表述為:地表震源激發(fā)的入射波場向下傳播,入射波場與非均勻體作用形成二次源,并產(chǎn)生次生波場(散射波場、反射波場等),次生波場向上傳播到地表被檢波器接收。在上述的直觀表述中,地震波的傳播和散射(反射)是相互獨立的,但在描述地震波場運動狀態(tài)的通用波動方程(4)中,地震波的傳播與散射(反射)是交織在一起的。這是因為方程(4)中的模型m不僅控制了地震波的傳播,也控制了地震波的散射(反射)。如果方程(4)中的模型m為光滑模型,則波動方程(4)就不能描述地震波的散射、反射和波型轉(zhuǎn)換,即波動方程(5)中的波場ui不包含散射、反射和轉(zhuǎn)換波。為了研究光滑模型下地震波的波場模擬、反演和偏移成像方法,如反射波形反演、反射波阻抗反演、反射數(shù)據(jù)的偏移成像(包括最小二乘偏移成像)方法等,我們首先需要建立光滑模型下描述地震波的散射和反射的波動方程,并以此得到適合地震數(shù)據(jù)偏移成像的地震數(shù)據(jù)嚴(yán)謹(jǐn)正演表達公式。
對通用波動方程(4)應(yīng)用擾動法[63],即分別引入模型m和波場u的擾動表達:m=mb+δm,u=ui+δu,其中,mb為光滑的背景模型(簡稱為光滑模型),δm為模型擾動(非均勻體的表示),ui為光滑模型mb下的入射波場(也稱為背景波場),δu為模型擾動δm產(chǎn)生的擾動波場。將模型與波場的擾動表達代入方程(4),得到與其對應(yīng)的擾動方程,即非均勻體產(chǎn)生的擾動波場所滿足的波動方程。
L(mb)δu=-δL(mb,δm)u
(7)
式中:L(mb)為光滑模型下的波動算子;δL(mb,δm)為擾動算子,有δL(mb,δm)=L(m)-L(mb)=L(mb+δm)-L(mb)。
地震數(shù)據(jù)是一種具有特定頻帶范圍的波形數(shù)據(jù),因此利用地震數(shù)據(jù)研究地下非均勻體結(jié)構(gòu)必須考慮地震波長與非均勻體δm的空間尺寸或者非均勻體δm物性空間變化的特征尺度之間的相對大小關(guān)系。如果非均勻體δm的尺寸a或者其物性空間變化特征尺度a相對于地震波場u的波長λ比較小,通常假定a/λ≤1,則非均勻體δm相對于波長在空間上只有有限的延續(xù)度,可視為產(chǎn)生散射波的局部散射體,相應(yīng)的波動方程(7)化為散射波所滿足的波動方程,即:
L(mb)us=-δLs(mb,δm)u
(8)
式中:us表示散射體產(chǎn)生的散射波場;δLs(mb,δm)表示與擾動δm或相對擾動δm/mb有關(guān)的散射算子,它與波場u的相互作用是激發(fā)散射波的虛源。方程(8)也稱為描述散射波的通用方程。由于δm相對波長比較小,因此方程(8)的右端源可視為點源,產(chǎn)生的散射波沒有特定的傳播方向。利用背景模型下波動方程(5)所對應(yīng)的Green函數(shù),將方程(8)表示為時空域積分形式,有:
(9)
式中:Ω為散射波源的分布空間;g為光滑模型mb下波動方程(5)的Green函數(shù);“*t”表示時間褶積。
如果非均勻體δm的尺寸a或者其物性空間變化特征尺度a相對于地震波場u的波長λ比較大,通常假定a/λ>1,則非均勻體δm相對于波長在空間上具有一定的延續(xù)度,可視為產(chǎn)生反射波的反射體,相應(yīng)的波動方程(7)化為反射波所滿足的波動方程,即:
L(mb)ur=-δLr(mb,Ir,σ)u
(10)
式中:ur表示反射體產(chǎn)生的反射波場;δLr(mb,Ir,σ)表示與波阻抗相對擾動Ir和角度σ有關(guān)的反射算子,它與波場u的相互作用是激發(fā)反射波的虛源;σ表示入射方向與反射方向之間的開角,它是由波動算子中物理參數(shù)的空間導(dǎo)數(shù)和波場的空間導(dǎo)數(shù)在高頻近似條件下的表達式引入的;波阻抗相對擾動Ir不同于散射體的物性參數(shù)相對擾動δm/mb,它不僅與mb,δm有關(guān),還與反射開角σ有關(guān)。方程(10)也稱為描述反射波的通用方程。由于δm相對于波長在空間上具有一定延續(xù)度,因此方程(10)的右端源可視為由點源集成的局部平面波源,產(chǎn)生的反射波可視為具有特定傳播方向的局部平面波。利用背景模型下波動方程(5)所對應(yīng)的Green函數(shù),同樣可將方程(10)表示為積分形式,有:
(11)
式中:Ω為局部平面波源的分布空間。
我們根據(jù)非均勻體δm的尺寸a或者其物性空間變化特征尺度a與地震波長λ之間的相對大小關(guān)系,將δm分別劃分為產(chǎn)生散射波的散射體和產(chǎn)生反射波的反射體,并分別得到了描述散射波和反射波的通用方程(公式(8)和公式(10)),以及它們的積分形式(公式(9)和公式(11))。
利用對δm的小擾動假定和一階Born近似(忽略多次散射波),方程(8)退化為描述一次散射波us的非齊次波動方程,即:
L(mb)us=-δLs(mb,δm)ui
(12)
相應(yīng)的方程(9)退化為:
(13)
利用對Ir的小值假定和一次反射波近似(忽略多次反射波),方程(10)退化為描述一次反射波ur的非齊次波動方程,即:
L(mb)ur=-δLr(mb,Ir,σ)ui
(14)
相應(yīng)的方程(11)退化為:
(15)
如果反射體δm內(nèi)部的物性變化很小或不變化,則反射主要是由反射體邊界產(chǎn)生[64]。為了描述反射體邊界對反射的作用和滿足對反射體邊界成像的需要,本文將產(chǎn)生反射的波阻抗相對擾動Ir沿入射波傳播方向I的方向?qū)?shù)定義為反射體邊界局部切平面的局部反射系數(shù),有:
(16)
式中:r表示局部反射系數(shù),是局部平面波入射到反射體邊界位置x處局部切平面上的局部反射系數(shù),它不僅可以表征反射面的空間位置,還能反映反射面上的物性變化。Ir在數(shù)學(xué)形態(tài)上是一個階躍函數(shù),而r在數(shù)學(xué)形態(tài)上是一個沿反射體邊界分布的δ函數(shù)。圖1為分界面上點x處局部切平面的反射示意圖。圖中Lp為點x處的局部切平面;r為點x處的局部反射系數(shù);ui和ur分別為具有局部平面波性質(zhì)的入射波和反射波;σ為點x處的反射開角。
圖1 局部切平面的反射示意
假定地震波長相對于光滑模型mb空間變化特征尺度滿足高頻近似條件(即光滑模型mb的空間變化相對于地震波長尺度可視為緩慢變化甚至不變化),則公式(16)中的入射波傳播方向?qū)?shù)在波數(shù)域?qū)?yīng)為iki,ki為入射波傳播方向的波數(shù),有ki=ω/[vb(x)],也稱為入射波傳播方向在x處的局部波數(shù),ω為地震波的圓頻率,vb(x)為點x處光滑模型的速度。因此公式(16)在高頻近似條件下的波數(shù)域表達式為:
(17)
公式(16)和公式(17)定義的局部反射系數(shù)r是地下物性參數(shù)相對擾動空間變化的反映。不同于公式(6)中與頻率無關(guān)的無量綱鏡面反射系數(shù)R,本文定義的局部反射系數(shù)r與傳播方向有關(guān)且具有長度倒數(shù)的量綱,它是高頻近似條件下,局部平面入射波作用于局部反射平面的反射系數(shù),所以與入射波傳播方向上的局部波數(shù)有關(guān)(在形式上與頻率有關(guān))。
利用公式(16)和公式(17)定義的局部反射系數(shù)r,可將方程(14)和方程(15)分別改寫為:
L(mb)ur=-δLr(mb,r,σ)ui
(18)
(19)
式中:δLr(mb,r,σ)表示與局部反射系數(shù)r有關(guān)的反射算子。方程(18)的右端源也是一個局部平面波源。
利用前文推導(dǎo)出的基于散射體物性參數(shù)擾動的一次散射波通用方程(12)和方程(13)以及基于反射體波阻抗相對擾動的一次反射波通用方程(14)和方程(15)或基于反射體邊界局部反射系數(shù)的一次反射波通用方程(18)和方程(19),結(jié)合具體的標(biāo)量波、聲波和各向同性彈性波方程,我們給出散射數(shù)據(jù)和反射數(shù)據(jù)線性正演表達的具體形式。
如果通用波動方程(4)中模型物性參數(shù)m僅為速度v(x),則方程(4)為非齊次標(biāo)量波動方程,即:
(20)
根據(jù)擾動方法,有v(x)=vb(x)+δv(x),u(x,t;xs)=ub(x,t;xs)+δu(x,t;xs),其中,vb(x)為光滑速度模型,δv(x)為速度模型擾動,ub(x,t;xs)為光滑速度模型下的背景波場,在下文也稱為入射波場ui(x,t;xs),δu(x,t;xs)為速度模型擾動產(chǎn)生的擾動波場。
假定速度擾動體δv(x)的大小或其速度空間變化特征尺度相對于地震波長比較小,即滿足上文定義的散射體條件,速度擾動體δv(x)可視為散射體。根據(jù)通用的一次散射波方程(12)和方程(13),可得到下述標(biāo)量散射波方程:
(21)
(22)
(23)
方程(22)也稱為基于速度相對擾動的標(biāo)量波地震散射數(shù)據(jù)的線性正演表達。
假定速度擾動體δv(x)的大小或其速度空間變化特征尺度相對于地震波長比較大,即滿足前文定義的反射體條件,速度擾動體δv(x)可視為反射體。根據(jù)通用的一次反射波方程(14),可得到下述標(biāo)量反射波方程:
(24)
式中:ur(x,t;xs)表示由速度相對擾動av(x)產(chǎn)生的一次反射波。
根據(jù)方程(15),由方程(24)可得到下述基于速度相對擾動的積分形式標(biāo)量波地震反射波數(shù)據(jù)線性正演表達,即:
(25)
由于標(biāo)量波動算子中不含有對速度參數(shù)的空間導(dǎo)數(shù),其對應(yīng)的反射算子也就不含有角度σ,因此反射數(shù)據(jù)的表達式(25)與散射數(shù)據(jù)的表達式(22)有相同的形式,但其中av(x)的分布范圍Ω大小不同。
如果反射體δv(x)內(nèi)部的速度變化很小或不變化,則反射可視為由反射體邊界產(chǎn)生。根據(jù)反射體邊界局部切平面的局部反射系數(shù)定義式(17),對于速度相對擾動,有下述的局部反射系數(shù)表示式:
(26)
利用公式(26)和基于局部反射系數(shù)的通用一次反射波方程(18),可得到基于局部反射系數(shù)的標(biāo)量反射波方程:
(27)
根據(jù)方程(19),由方程(27)可以得到基于局部反射系數(shù)r的積分形式標(biāo)量反射波數(shù)據(jù)線性正演表達式:
(28)
基于Kirchhoff近似,BLEISTEIN等[36]推出了與公式(28)相類似的基于公式(6)中的鏡面反射系數(shù)R的反射數(shù)據(jù)表達式,由于公式(6)中鏡面反射系數(shù)R的嚴(yán)格條件,我們認(rèn)為BLEISTEIN等推導(dǎo)出的反射數(shù)據(jù)表達式對于非平面入射波和反射面為非無窮大平邊界存在局限性,不適合用于構(gòu)建實際地震數(shù)據(jù)的偏移成像。
考慮到地下密度變化對地震波傳播的影響,假定通用波動方程(4)中模型m為包含速度v(x)和密度ρ(x)的聲學(xué)介質(zhì),則方程(4)為非齊次聲波方程,即:
=f(t)δ(x-xs)
(29)
假定δv(x)和δρ(x)的大小或其空間變化特征尺度相對于地震波長滿足散射體條件,即δv(x)和δρ(x)對應(yīng)的非均勻體可視為產(chǎn)生散射波的散射體。根據(jù)通用的一次散射波方程(12)和方程(13),可得到下述聲散射波方程:
(30)
(31)
ui(x,t;xs)=f(t)δ(x-xs)
(32)
對方程(31)應(yīng)用分步積分法,有:
(33)
方程(33)也稱為基于速度、密度相對擾動的聲散射波數(shù)據(jù)的線性正演表達。
假定δv(x)和δρ(x)的大小或其空間變化特征尺度相對于地震波長滿足反射體條件,即δv(x)和δρ(x)對應(yīng)的非均勻體可視為產(chǎn)生反射波的反射體。根據(jù)方程(15)和方程(33),可得到下述基于速度和密度相對擾動的聲反射波方程,即:
(34)
式中:ur(xg,t;xs)為由av和aρ產(chǎn)生的一次反射波數(shù)據(jù)。
將方程(34)變換到頻率域,有:
(35)
(36)
將方程(36)變換到時間域,有:
(37)
方程(37)對應(yīng)的微分形式方程為:
(38)
公式(38)中的av+aρ(1+cosσ)是一個與開角σ有關(guān)的聲波阻抗相對擾動。令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),我們可以得到基于反射體聲波阻抗相對擾動的聲波反射數(shù)據(jù)線性正演表達公式,即
(39)
為了得到基于反射體邊界局部反射系數(shù)的反射波動方程和反射數(shù)據(jù)線性正演表達公式,我們定義聲反射的局部反射系數(shù)為聲波阻抗相對擾動Ir(x,σ)沿入射波傳播方向I的方向?qū)?shù)。根據(jù)反射體邊界局部切平面的局部反射系數(shù)定義式(17),對于聲波阻抗相對擾動,有下述的局部反射系數(shù)表達式:
cosσ)]
(40)
利用(40)式定義的局部反射系數(shù),由方程(38)在高頻近似條件下可得到基于反射體邊界局部反射系數(shù)的聲反射波方程:
(41)
根據(jù)方程(19),由方程(41)可以得到基于反射體邊界局部反射系數(shù)的積分形式聲波反射數(shù)據(jù)線性正演表達公式,即
(42)
如果在聲波方程中不考慮密度變化,則聲波方程(29)退化為標(biāo)量波方程(20),相應(yīng)的方程(33)、方程(39)和方程(42)分別退化為方程(22)、方程(25)和方程(28)。
為了更真實地描述地下地震波的運動,假定方程(4)中地下模型為各向同性彈性介質(zhì),其模型m的物性參數(shù)為Lamé系數(shù)λ(x)、μ(x)和密度ρ(x),則方程(4)為非齊次彈性波動方程,即:
Lu(x,t;xs)=f(t,xs)
(43)
式中:u(x,t;xs)為位移矢量波場,在笛卡爾坐標(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是一個3×3偏微分算子:
(44)
假定擾動量δα(x)、δβ(x)、δρ(x)的大小或其空間變化特征尺度相對于地震波長滿足前文定義的散射體條件,即δα(x)、δβ(x)、δρ(x)對應(yīng)的非均勻體可視為產(chǎn)生散射波的散射體。根據(jù)通用散射波方程(12),可得到下述通用的彈性散射波方程:
Lbus(x,t;xs)=ΔLsui(x,t;xs)
(45)
式中:us(x,t;xs)為彈性一次散射波的位移矢量波場;Lb為光滑模型下的彈性波傳播算子;ΔLs為彈性波散射算子,有:
(46)
Lbui(x,t;xs)=f(t,xs)
(47)
利用光滑模型彈性波方程(47)所對應(yīng)的矢量Green函數(shù),可將方程(45)中的彈性散射波數(shù)據(jù)寫成與方程(13)對應(yīng)的積分形式[65]。方程(45)也稱為基于縱、橫波速度、密度相對擾動的彈性散射波數(shù)據(jù)的線性正演方程。
假定δα(x),δβ(x),δρ(x)的大小或其空間變化特征尺度相對于地震波長滿足前文定義的反射體條件,即δα(x),δβ(x),δρ(x)對應(yīng)的非均勻體可視為產(chǎn)生反射波的反射體。根據(jù)通用反射波方程(14),利用與前兩節(jié)類似的推導(dǎo)可得到彈性反射波方程,即:
Lbur(x,t;xs)=ΔLrui(x,t;xs)
(48)
式中:ur(x,t;xs)為彈性一次反射波的位移矢量波場;ΔLr為彈性波反射算子。
(49)
(50)
(51)
(52)
(53)
(54)
(55)
在法向入射時,σ=0,有:
(56)
(57)
(58)
此時不發(fā)生波型轉(zhuǎn)換。
(59)
(60)
(61)
(62)
(63)
(64)
將(59)式代入方程(49),有:
(65)
(66)
式中:g(xg,x,t)為光滑模型下彈性波方程(47)所對應(yīng)的矢量Green函數(shù)。
(67)
對于SH波入射SH波反射的局部反射系數(shù),有:
(68)
對于SV波入射SV波反射的局部反射系數(shù),有:
(69)
對于P波入射SV波反射的局部反射系數(shù),有:
(70)
對于SV波入射P波反射的局部反射系數(shù),有:
(71)
(72)
(73)
(74)
(75)
由上述推導(dǎo)出的標(biāo)量波、聲波和彈性波的反射數(shù)據(jù)線性正演表達公式可知,基于波阻抗相對擾動的反射數(shù)據(jù)線性正演表達公式的右端包含入射波場的時間二階導(dǎo)數(shù)與波阻抗相對擾動的乘積,而基于局部反射系數(shù)的反射數(shù)據(jù)線性正演表達公式的右端包含入射波場的時間一階導(dǎo)數(shù)與局部反射系數(shù)的乘積。不同于平面波無窮大平邊界的無量綱鏡面反射系數(shù),本文基于高頻近似和局部平面波定義的局部切平面反射系數(shù)是一個具有長度倒數(shù)量綱的變量,是真實模型與光滑模型之間模型參數(shù)相對擾動沿入射波傳播方向的方向?qū)?shù),它不僅依賴于光滑模型,也與入射波傳播方向有關(guān)。Zoeppritz方程不僅描述了反射波與入射波之間的關(guān)系,還定義了平面波無窮大平邊界鏡面反射系數(shù)與入射方向和地下模型物理參數(shù)之間的定量關(guān)系。本文利用反射體近似和高頻近似推導(dǎo)的基于反射體彈性波阻抗相對擾動的P/S波型彈性反射波方程與Zoeppritz方程類似,不僅建立了反射波與入射波之間的數(shù)學(xué)物理關(guān)系,也定義了反射體彈性波阻抗相對擾動與入射方向和地下模型物理參數(shù)之間的定量關(guān)系,可為進一步的基于偏移成像得到的角度域共成像點道集的巖性分析與反演提供理論基礎(chǔ)。由于本文定義的局部反射系數(shù)與入射波傳播方向的局部波數(shù)有關(guān),因此在目前我們還不能給出傾斜入射情況下局部反射系數(shù)與地下模型物理參數(shù)和入射波傳播方向之間的具體函數(shù)關(guān)系。
當(dāng)前地震數(shù)據(jù)波動方程偏移成像方法不區(qū)分散射波和反射波,也缺乏與之相應(yīng)的嚴(yán)謹(jǐn)正演方程,致使偏移成像不能準(zhǔn)確地利用地震數(shù)據(jù)的波形信息,得到的偏移成像結(jié)果存在相位誤差。為了使波動方程偏移成像方法能更好地適應(yīng)地下實際情況和準(zhǔn)確地利用波形信息,彌補當(dāng)前偏移成像方法理論的不足,修正地震數(shù)據(jù)波動方程偏移成像結(jié)果的相位誤差,并為發(fā)展地層物性參數(shù)和儲層參數(shù)成像方法打下基礎(chǔ),本文在給定滿足地震波運動學(xué)的準(zhǔn)確的光滑模型基礎(chǔ)上,利用擾動方法將波動方程轉(zhuǎn)化為描述非均勻體產(chǎn)生的擾動波場的波動方程,根據(jù)地下非均勻體的大小或非均勻體物理性質(zhì)空間變化特征尺度與地震波長之間的相對大小關(guān)系,提出將非均勻體劃分為相對于地震波長在空間上具有有限延續(xù)度的產(chǎn)生散射波的散射體或相對于地震波長在空間上具有一定延續(xù)度的產(chǎn)生反射波的反射體,推導(dǎo)出反射波波動方程。為區(qū)別于散射波,本文認(rèn)為反射波是由點源集成的局部平面波源產(chǎn)生的具有方向性的局部平面波。對于散射波,利用小擾動假定和一階Born近似(忽略多次散射波),可得到基于散射體物性參數(shù)相對擾動的地震一次散射數(shù)據(jù)的線性正演表達。對于反射波,利用波阻抗相對擾動的小值假定、一次反射波近似(忽略多次反射波)和高頻近似,可得到基于反射體波阻抗相對擾動的地震一次反射數(shù)據(jù)的線性正演表達。反射體的波阻抗相對擾動不同于散射體的物性參數(shù)相對擾動,它不僅與反射體的物性參數(shù)相對擾動有關(guān),還與反射開角有關(guān)?;诜瓷潴w彈性波阻抗相對擾動的P/S波型彈性反射波方程,不僅建立了反射波與入射波之間的數(shù)學(xué)物理關(guān)系,也定義了反射體彈性波阻抗相對擾動與入射方向和地下模型物理參數(shù)相對擾動之間的定量關(guān)系。對于高頻近似條件下定義的反射體與反射波,本文將波阻抗相對擾動沿入射波傳播方向的方向?qū)?shù)定義為反射體邊界局部切平面對于局部平面波的局部反射系數(shù),得到基于局部反射系數(shù)的地震一次反射數(shù)據(jù)的線性正演表達。不同于平面波入射到無窮大平界面的與頻率無關(guān)的無量綱鏡面反射系數(shù),本文定義的局部反射系數(shù)具有長度倒數(shù)的量綱且與入射波傳播方向有關(guān)(即入射波的波數(shù)),推導(dǎo)出的反射數(shù)據(jù)表達式具有更廣的適應(yīng)性。本文推導(dǎo)出的地震數(shù)據(jù)線性正演表達(特別是基于反射體波阻抗相對擾動的地震反射數(shù)據(jù)線性正演表達),是開展地震數(shù)據(jù)波形成像和其它波形線性反演方法研究的基礎(chǔ)方程。