第1章地震偏移成像基礎(chǔ)

第一章 地震偏移成像基礎(chǔ)地震偏移技術(shù)是現(xiàn)代地震勘探數(shù)據(jù)處理的三大基本技術(shù)之一。它是在過去的古典技術(shù)上發(fā)展起來的，其它兩大技術(shù)都是從其它相關(guān)學(xué)科引進(jìn)到地震中來的。所以，偏移技術(shù)具有地震勘探本身的特征。但是，地震偏移方法本身由于使用計算機(jī)而引起了許多革命性的變化。這就是把它從研究簡單的探測目標(biāo)的幾何圖形進(jìn)而發(fā)展成研究反射界面空間的波場特征、振幅變化和反射率等。本章主要介紹地震偏移成像技術(shù)的基礎(chǔ)知識。首先給出偏移成像的概念；第二節(jié)介紹有限差分法的基礎(chǔ)知識；第三節(jié)敘述基于波動方程的波場外推與地震成像原理；第四節(jié)討論波場外推的Kirchhoff積分法；第五節(jié)簡單分析Born近似和Rytov近似；最后闡述基于DeWolf近似、薄板近似、屏近似和相屏傳播算子計算反向散射波場的方法?！?.1 偏移成像的概念反射地震方法是根據(jù)在地面上以一定方式進(jìn)行彈性波激發(fā)，并在地面的一定范圍（孔徑）內(nèi)記錄來自地下彈性分界面的反射波來研究地下地質(zhì)巖層結(jié)構(gòu)及其物性特征的一種方法。因此，也可以把它看做是一種反散射問題。就反射地震觀測方式的特點，它的成像問題要分做兩步，第一步是按照一定的方式記錄到達(dá)地面的反射波，第二步用計算機(jī)按一定的計算方法對觀測數(shù)據(jù)進(jìn)行處理，使之成為反映地下地質(zhì)分層面位置及反射系數(shù)值的反射界面的像。而地震偏移技術(shù)就是在第二步過程使反射界面最佳地成像的一種技術(shù)。地震偏移可在疊前做也可在疊后做。疊前偏移是把共炮點道集記錄或共偏移距道集記錄中的反射波歸位到產(chǎn)生它們的反射界面上并使繞射波收斂到產(chǎn)生它的繞射點上。在把反射波回投到反射界面上和繞射波收斂到繞射點上時要去掉傳播過程的效應(yīng)， 如擴(kuò)散與衰減等。最后得到能夠反映界面反射系數(shù)特點的并正確歸位了的地震波形剖面，即偏移剖面。疊后偏移是在水平疊加剖面的基礎(chǔ)上進(jìn)行的，針對水平疊加剖面上存在的傾斜反射層不能正確地歸位和繞射波不能完全收斂的問題， 采用了爆炸反射面的概念來實現(xiàn)傾斜反射層的正確歸位和繞射波的完全收斂。地震偏移的效果見圖 1-1和圖1-2。地震偏移的類型見表 1-1。地震偏移技術(shù)在二十世紀(jì)六十年代以前是用手工操作的一種制圖技術(shù)， 只是用來求得反射點的空間位置，而不考慮反射波的特點。它是一種古典的偏移方法。早期的計算機(jī)偏移方法是在古典的偏移方法的基礎(chǔ)上提出來的。其中有的成功了，有的失敗了。成功的是那些符合波的傳播特征的方法。盡管這些方法使用了波前、繞射等地震波傳播的惠更斯原理，但只是定性的、概念性的。偏移剖面的質(zhì)量雖然能夠滿足最基本的要求，但歸位的精度和成像時的波形特征都不是很準(zhǔn)確的。因此，研究更有效的地震偏移方法是很迫切的。二十世紀(jì)七十年代初J.Claerbout教授首先提出了用有限差分法解單程波動方程的近似式，用地面觀測的地震數(shù)據(jù)重建地震波在地下傳播過程中的波場，從這些傳播過程的波場中提取使地震界面成像的那些數(shù)據(jù)，組成地震偏移剖面。由于這種偏移方法在計算過程中要解波動方程或其近似式，所以被稱為波動方程法偏移技術(shù)。以后，F(xiàn)rench和Schneider等在繞射偏移法的基礎(chǔ)上使用了波動方程解的Kirchhoff積分公式，發(fā)展為地震偏移的波動方程積分法。使繞射偏移建立在可靠的波的基本原理上。因而改善了偏移剖面，取得了良好的效果。圖1-1(a)共中心點疊加剖面， (b)偏移剖面，(c)明顯的繞射波 D、偏移前的傾斜同相軸 B和偏移后的傾斜同相軸A。偏移使傾斜同相軸 B歸位到它的真實地下界面 A，并使繞射波 D收斂到其頂點 P。點畫線指出了鹽丘的邊界。圖1-2偏移前(a)及偏移后(b)彎曲反射界面 (向斜和背斜 )的形狀。詳細(xì)情況見正文 (模型據(jù) Union OilCompany)在二十世紀(jì)七十年代后期， Stolt和Gazdag等又先后提出了在頻率 -波數(shù)域解波動方程，外推地震波場的方法。這種方法被稱為F-K域偏移方法。由于該方法計算簡單，效率高，因而很快得到了推廣。上述三種波動方程偏移方法是同時并存的，因為它們各有自己的特點，因而不能用一個方法來取代其它方法。使用時視具體條件和要求決定采用何種方法。波動方程偏移方法在最近20年間迅速發(fā)展并不斷完善，許多人對此做出了有益的貢獻(xiàn)。其中，Loewenthal等人的爆炸反射面的概念對于理解疊加剖面的偏移成像具有很大價值，Hubral，Larner等人提出的深度偏移的概念具有很大意義，Berkhout提出的偏移過程是一個空間褶積的概念對于偏移的橫向分辨力的理解很有益處，馬在田院士提出的高階方程的分裂算法對提高有限差分法偏移的精度有很大貢獻(xiàn)，Yilmaz等提出的雙平方根法為解決疊前偏移奠定了基礎(chǔ)?，F(xiàn)在仍有許多學(xué)者還在探索波動方程偏移技術(shù)，以期更加完善該方法。表1-1偏移方法分類類型論述疊加解釋人員總是想要的剖面。法向射線深度轉(zhuǎn)換嚴(yán)格適用于沒有構(gòu)造傾角且速度只隨深度變化的情況。時間偏移適用于疊加剖面上有繞射波或構(gòu)造傾角以及速度有垂向變化的情況；速度的橫向變化不大時也能用。深度偏移用于疊加剖面上有構(gòu)造傾角和強(qiáng)橫向變速的情況。疊前部分偏移疊后偏移適用于疊加剖面與零炮檢距剖面等價的情況，但不適合具（PSPM）有不同疊加速度的地層傾角不一致或強(qiáng)橫向變速的地區(qū)，疊前部分偏移(傾角時差(DMO)校正)能夠為疊后偏移提供更好的疊加剖面，但疊前部分偏移只解決具有不同疊加速度的地層傾角不一致的問題。疊前全時間偏移輸出偏移剖面，不產(chǎn)生未經(jīng)偏移的中間疊加剖面，所以不太受歡迎，因為解釋人員普遍喜歡既有疊加剖面又有偏移剖面的解釋方式。但不論如何這是解決傾角不一致地層問題的最精確方法。疊前部分偏移是這種處理方法的一種簡化。疊前深度偏移用于嚴(yán)重橫向變速的情況，這時已無法作合適的疊加處理。三維疊后時間偏移用于疊加剖面上出現(xiàn)來自射線平面以外的傾斜同相軸(即垂直測線方向)的情況，這是疊后最常用的一種三維偏移方法。三維疊后深度偏移用來解決與三維地下復(fù)雜構(gòu)造有關(guān)的強(qiáng)橫向變速問題。三維疊前時間偏移用于疊前部分偏移不適用且疊加剖面上有橫向傾斜層反射的情況。三維疊前深度偏移只要計算機(jī)機(jī)時允許，并且又能精確知道三維速度模型，這是人人樂于接受的處理方法?！?.2 有限差分法的基礎(chǔ)知識在計算機(jī)上進(jìn)行數(shù)值運(yùn)算，使用的是離散的和有限的數(shù)值，而不是連續(xù)的和無限的函數(shù)，為此要為離散數(shù)值的計算建立基本方法。最基本和廣泛使用的方法就是有限差分法，借助這一方法可以研究連續(xù)物理系統(tǒng)的性質(zhì)，近似地、但相當(dāng)精確地解出各種數(shù)理方程問題。本節(jié)將概要地介紹有限差分法的基礎(chǔ)知識，供以后各章解偏移方程使用。一．差分方程的建立與求解1．有限差分的概念我們感興趣的是用有限差分法解各類微分方程，因此，要把導(dǎo)數(shù)用有限差分來近似，所以我們這里只研究用有限差分近似導(dǎo)數(shù)的方法。（1）一元函數(shù)差分法當(dāng)一個函數(shù)u和它的各階導(dǎo)數(shù)是變量x的單值的、有限的和連續(xù)的函數(shù)時，可以用泰勒原理展開為：u(xx)u(x)xdu1x2d2u1x3d3u（1-1）dx2dx26dx3和u(xx)u(x)xdu1x2d2u1x3d3u（1-2）dx2dx26dx3以上二式相加，得：u(xx)u(xx)2u(x)x2d2u(x4)（1-3）dx2式中O(x4)表示包含有x的四階和高于四階以上的項。從（1-1）~（1-3）式我們可以導(dǎo)出（圖1-3）下列表示式：一階向前差分：xu1[u(xx)u(x)]xdu(x)（1-4）dx一階向后差分：xu1[u(x)u(xx)]xdu(x)（1-5）dx一階中心差分：xu1[u(xx)u(xx)]2xdu(x2)dx二階差分：21xxuxux2[u(xx)d2u2（1-7）dx2(x)（2）多元函數(shù)差分法圖1-3一元函數(shù)差分法

（1-6）2u(x) u(x x)]圖1-4二元函數(shù)差分網(wǎng)格設(shè)有一個二元函數(shù)u(x,t)，我們用網(wǎng)格把它們離散（圖1-4），令x=ix，t=jti和j為整數(shù)。用ui,j=u(ix,jt)表示各網(wǎng)格點上的函數(shù)值?，F(xiàn)在我們用泰勒級數(shù)展開下面各點之值為：ui1,ju[(i1)x,jt]u(xix,tj)ui,jx(u)i,j1x2(2u2)i,j1x3(3u3)i,j（1-8）x2x6xui1,ju[(i1)x,jt]u(xix,tj)ui,ju122u13(3u（1-9）x()i,jx(2)i,jxx3)i,jx2x6由此可以求出對 x的一階和二階差分為：xu1(ui1,jui,j)xu（1-10）(x)x1xu(ui,jui1,j)xu（1-11）(x)xxu1ui1,j)(ui1,j2 x(x2)（1-12）xxxux2u1(ui1,j2ui,jui1.j)2x2u(x2)（1-13）x2同理可求出對t的一階和二階差分：tu1(ui,j1ui,j)tu(t)（1-14）ttu1(ui,jui,j1)tu(t)（1-15）ttu1(ui,j1ui,j1)2tu(t2)（1-16）tttut2u12(ui,j12ui,jui,j1)t2u(t2)（1-17）t22．建立差分方程建立差分方程的方法很多，有積分法，物理量守恒法，變分法和最小平方法等。對于常系數(shù)的微分方程來說， 積分法是最簡便的。用積分法構(gòu)造差分方程的過程見參考文獻(xiàn) [1]。3．差分方程的格式差分方程的格式基本可分為兩大類：即顯式格式和隱式格式。在實際工作中又可以衍生出許多格式，甚至可以用顯式與隱式聯(lián)合形式的差分格式。在這里我們僅以拋物型偏微分方程為例說明常用的顯式格式和 Crank-Nicolson格式。1）顯式格式設(shè)有拋物型偏微分方程：u2u0（1-18）tx2經(jīng)推導(dǎo)可具體寫出差分方程為：ui,j1ui,jx2t(ui1,j2ui,jui1,j)（1-19）t所用的差分網(wǎng)格如圖1-5所示。令，則（1-19）式可寫為：x2ui,j1ui,j(ui1,j2ui,jui1,j)（1-20）其中ui,j，ui1,j和ui1,j為已知值。從第j時間層上的已知值， 可用（1-20）式直接計算出第 j 1時間層的值。以此類推，可解出全時間上的物理量值。圖1-5顯式差分格式 圖1-6隱式差分格式2）Crank-Nicolson隱式格式雖然顯式法計算上很簡單，但它有一個嚴(yán)重的缺點，即時間步長 t一定要很小，必須滿足0t1，才能保持計算上穩(wěn)定并達(dá)到必要的精度。Crank-Nicolson在1947年x22提出了一個使所有的有限r(nóng)值都滿足計算要求（收斂性和穩(wěn)定性）的方法。他們把2u/x2項用第j和第j+1時間行上的平均差分來逼近。求出下列的差分方程（圖1-6）：ui,j1ui,jui1,j12ui,j1ui1,j1ui1,j2ui,jui1,j（1-21）t2x2x2由此得出下列等式：ui1,j1 2(1 )ui,j1 ui1,j1 ui1,j 2(1 )ui,j ui1,j （1-22）其中，t2。（1-22）式的左邊包含有3個未知數(shù)值，右邊3個值是已知的。x如果在每個時間層上有I個格點，則方程（1-22）根據(jù)零時間層(j=0)計算第一時間層(j=1)的網(wǎng)格點上的函數(shù)值ui1時可列出I個聯(lián)立的方程組。由于j=1時間層上的初值和邊值是給出的，因此方程組表現(xiàn)為三對角線陣。對于這樣的方程組，可用代數(shù)的追趕法求解，也可以用矩陣方程來分析。4．差分方程的解法1）追趕法1-22）方程組可用（1-23）式表示：iui 1 iui iui 1 i其中，0 i I。補(bǔ)充的邊界條件為：0u1 0u0 0IuI 1 IuI I

1-23）1-24）我們要求出未知值ui。用遞推法求解，即已知i點上的ui值可推出i+1點上的ui1值。這就要先求出中間助算未知數(shù)xi和yi，其關(guān)系式為：ui1xiuiyi（1-25）為此要求出xi,yi。（1-25）式對任何i點都是成立的。我們把它代入（1-23）式，得：i(xiui yi) iui iui 1 i或者寫為：uiiui1iiyi（1-26）ixiixiii這個公式與（1-25）式相同，它把i點和i1點用（1-27）式：uixi1ui1yi1（1-27）連結(jié)起來。由此求出下列等式：xi1i，yi1iiyi（1-28）xixiiiii用這個關(guān)系式可以求出所有的xi和yi的值。這樣，（1-25）和（1-28）就給出解三對角線方程組的兩個遞推過程。利用（1-28）式，我們在網(wǎng)格上從iI到i0求出所有的xi和yi值。在求出所有的xi和yi之后，從i0向iI用（1-25）式求出所有的ui1值。在iI點確定初始值xI1和yI1。在i0點從邊界條件確定u01值。把第I點的邊界條件（1-24）與遞推關(guān)系式（1-25）比較，則得：xI1IIyI1I（1-29）I在用遞推關(guān)系求出所有的xi和yi，一直到x0和y0之后，則用（1-24）和（1-25）式可求出未知函數(shù)值u01為：u01x000y0（1-30）0x00從u01出發(fā)用（1-25）式可求出所有的ui1值。根據(jù)上述遞推原理可求出所有時間層上的ui,j值。（2）矩陣法現(xiàn)在從矩陣觀點來分析三對角線方程組的求解問題。令w為已知量，需要求出的向量為u,它們滿足矩陣方程：Auw（1-31）其中A是三對角矩陣。我們引入平移的矩陣算子L和L，它們對向量u的分量有下列關(guān)系：L{ui}{ui1}（1-32）L{ui}{ui1}（1-33）由于A是三對角矩陣，故（1-31）可寫為：(PLQRL)uw（1-34）其中P，Q和R是對角矩陣。（1-34）式給出雙向遞推過程，可分為兩步。首先有下列關(guān)系式：LuXuy（1-35）求對角矩陣X和向量y。把（1-35）式與對角陣P相乘，寫為下列形式：(PLPXO)uPy（1-36）式中，O為零陣。從（1-34）式減去（1-36）式，得(QPX)uRLuwPy（1-37）或者寫為：u(QPX)1RLu(QPX)1(wPy)（1-38）1-38）式與（1-35）式在形式上是一致的。為了使它與所討論的問題相聯(lián)系，必須使它們等價。因此得出下列等式：LX(QPX)1RLy(QPX)1(wPy)（1-39）方程（1-35）和（1-39）組成兩步遞推過程。這些結(jié)果與標(biāo)量方法所得結(jié)果一致。因為P、Q、R和X是對角矩陣，所以它們每列的對角線上的元素相應(yīng)地是 i、i、i和xi。顯然，矩陣結(jié)果與代數(shù)方程是一致的。二．收斂性、相容性和穩(wěn)定性為了使差分方程的解能夠相當(dāng)精確地逼近于相應(yīng)的微分方程的解，我們必須討論它所需要的條件。這是不同的，但又互相關(guān)聯(lián)的兩個條件。第一個是近似的差分方程的精確解向微分方程的解收斂的條件， 第二個是在解差分方程過程中任何誤差不可無界地增長的條件。1．收斂性令u是變量x和t的偏微分方程的精確解，u?是用來近似偏微分方程的差分方程的精確解。當(dāng)t和xuu，就說有限差分方程是收斂的。趨近于零時，?趨近于某個固定點上的e u u?之差稱為離散誤差。討論收斂性問題的簡單例子請看參考文獻(xiàn)[1]。一般來說收斂性問題常常是難以研究的。因為離散誤差的最終表達(dá)式經(jīng)常含有未知導(dǎo)數(shù)項，無法估計其邊界值。但是，對于線性偏微分方程中的拋物型與雙曲型方程來說可以用穩(wěn)定性和相容性來研究收斂性。對此存在Lax等價原理。Lax等價原理該定理說，當(dāng)給出了一個適定的初值問題和它的有限差分方程時，如果該差分方程滿足相容性條件，則穩(wěn)定性是收斂性的充分而且必要的條件。關(guān)于這個定理的證明請看參考文獻(xiàn) [8]。下面我們給出相容性和穩(wěn)定性的定義， 研究相容性和穩(wěn)定性問題的實例和方法請看參考文獻(xiàn) [1]。2．相容性設(shè)F(?)0代表在(i,j)點上的差分方程。如果把?用u來代替，則()就稱為i,jFi,juu在網(wǎng)格點(i,j)上的局部截斷誤差。如果在網(wǎng)格步長趨于零時這個誤差亦趨于零，則稱差分方程與偏微分方程相容。3．穩(wěn)定性在解差分方程過程中，除了離散誤差之外，還有舍入誤差。這就是差分方程的理論精確解ui,j與計算機(jī)輸出的實際值Ni,j之差。全部誤差有離散誤差和總舍入誤差組成。??Ni,j)（1-40）i,jui,jNi,j(ui,jui,j)(ui,j離散誤差在一定的差分格式和差分網(wǎng)格步長下是可以控制的。而總的舍入誤差不僅與差分方程有關(guān)，而且與計算機(jī)舍入誤差有關(guān)，也與算法有關(guān)。差分解法是以步進(jìn)方式工作的。在逐步推進(jìn)過程中這種誤差也會逐漸積累。 這種誤差積累是保持有界還是惡性增長以致把解淹沒在積累的誤差中，這就是數(shù)值解的穩(wěn)定性問題。數(shù)值穩(wěn)定性是差分格式的必須條件。也是差分方程收斂性的必要條件。因此，在確定一個差分格式作為計算手段之前，一定要研究該差分格式的穩(wěn)定性。常用的有兩種研究穩(wěn)定性的方法。一個是 Fourier分析法，這是一種比較簡單的，但不考慮邊界條件的方法。另一種方法是矩陣法。它把方程組表示為矩陣，并考查與矩陣有關(guān)的本征值。此外還有能量法。§1.3 基于波動方程的波場外推與地震成像原理使用波動方程進(jìn)行偏移，首先就是要重建反射波的原來波場。反射界面上剛剛產(chǎn)生的反射波，就認(rèn)為是該反射面的像。為了使用波動方程公式進(jìn)行波場外推，一般是把波動方程分解為上行波方程和下行波方程。這是由于波動方程一般式既滿足上行波也滿足下行波。我們在成像中使用的是上行波方程。在這一節(jié)里我們將討論上行波與下行波方程，波場外推方法和地震成像的原理。一．上行波和下行波波動方程有兩個解，一般表示為 {exp[ i (t r/v)]}/r，其中 r是從源到觀測點的距離，是正實數(shù)。在地震勘探中一般取深度方向向下為正 z的方向。向正 z方向傳播的地震波稱為下行波，即用{exp[i(tr/v)]}/r所表示的波。向負(fù)z方向傳播的波為上行波，即用{exp[i(tr/v)]}/r代表的波。下行波即入射波，上行波為反射波。為了利用波動方程外推波場要求把它分解為上行波方程和下行波方程，然后才能利用它們進(jìn)行波場外推。否則將難以進(jìn)行外推的計算。從Claerbout導(dǎo)出的結(jié)果可以看出，既使只存在z方向介質(zhì)不均勻的情況下，上行波和下行波也是相互耦合的。只有在均勻各向同性完全彈性介質(zhì)的情況下上行波和下行波才是分離的。由于現(xiàn)代地震研究的主要是縱波，因此，我們可以用聲波方程來代替彈性波方程。設(shè)u和w表示x和z方向上的質(zhì)點位移分量，v表示縱波的傳播速度?？紤]到在地震勘探中我們研究的是無源問題，因此在均勻各向同性完全彈性介質(zhì)中有如下的二維波動方程：2u2u12u（1-41）x2z2v2t2對（1-41）式相對x和t做二維Fourier正變換，并進(jìn)行算子分解得到：2~22~2~2~dd~du(duikz)(dz22kx)udz2kzu(dzikz)u0（1-42）vdz其中利用了波散關(guān)系：kx2kz22v2（1-43）式中，表示圓頻率，kx和kz分別表示水平方向和深度方向上的波數(shù)。由（1-42）式得出：~~22~du（1-44dzikzui2kxu）v其中，正號代表上行波方程，負(fù)號代表下行波方程，~(kx,z,)表示波場ux,z,t的二維uFourier正變換。二．波場外推無論上行波還是下行波都可以用數(shù)學(xué)手段進(jìn)行外推，外推方向可以是正向的，也可以是反向的。所謂正向外推就是根據(jù)波在當(dāng)前位置上的振動情況向波的自然傳播方向用計算手段預(yù)測出波場。反向外推是向波的自然傳播方向的反方向上重建原來的波場。這樣，對每一種波來說，無論是上行波，還是下行波，都可以進(jìn)行正向外推和反向外推。對一個波場應(yīng)是進(jìn)行正向外推還是反向外推均有物理問題決定。1．上行波的外推把（1-44）式中的上行波方程改寫為：~22dui（1-45）~v2kxdzu~1-45）式兩邊取積分，其積分限為z和zz，有其中u表示上行波波場。對（zz~22zzduz~iv2kxzdzu積分結(jié)果為：~2kx2z)izu(zev2（1-46）~u(z)由此得出上行波的正、反向外推式。1）上行波正向外推公式上行波的正向外推式就是向負(fù)z方向的外推公式。從（1-46）式可求出為：2kx2~~iv2zz)e（1-47）u(z)u(z根據(jù)這個公式可以計算模擬反射波的地震記錄（地震圖）。2）上行波反向外推公式上行波的反向外推式就是向正z方向的外推公式。從（1-46）式可得出為：2kx2z~~iv2（1-48）u(zz)u(z)e根據(jù)這個公式可以進(jìn)行地震記錄的向下半空間延拓，求出地下任何一點的波場，實現(xiàn)地震波偏移的目的。2．下行波的外推把（1-44）式中的下行波方程改寫為：~2ddi2dz（1-49）~2kxdv~1-49）式進(jìn)行積分，仍取積分限為z和zz，得其中d表示下行波波場。對（~2z)ikx2zd(zev2（1-50）~d(z)據(jù)此可以得出下行波的正、反向外推公式。1）下行波正向外推公式下行波的正向外推式是指沿正z方向的外推。其外推式為：2~~ikx2zv2（1-51）d(zz)d(z)e這個方程可用來模擬下行波的地震記錄。2）下行波反向外推公式下行波的反向外推是指沿負(fù)z方向的外推。其外推式為：2~~z)eiv2kx2z（1-52）d(z)d(z1-52）式可用來從下行波場進(jìn)行反向求源的計算工作。從上面的外推公式中可以看出，上行波的正向外推公式與下行波的正向外推公式在形式上是一致的。上行波的反向外推公式與下行波的反向外推公式在形式上也是一致的。盡管它們在形式上一致，但外推方向不同，是正相反的（從坐標(biāo)系來看），同時各自的初始條件也不相同。這是應(yīng)當(dāng)注意的?，F(xiàn)在來分析波場本身的條件對外推結(jié)果的影響。這里所說的波場本身的條件就是外推方程中指數(shù)項的根號項：kzk2kx2（1-53）的情況，式中k222。當(dāng)kkx時，kz為正或負(fù)的實數(shù)，這時所有外推公式中存在虛v指數(shù)。說明在外推過程中波場發(fā)生相位變化。一般都能得出正確的結(jié)果。當(dāng)kxk時，kz值為虛數(shù)：kzikx2k2（1-54）當(dāng)把它代入各個外推公式中時，將得到如下形式的外推方程：~~kx2k2z（1-55）u(zz)u(z)e式中 z為正值，即 z zn1 zn。從（1-55）式可以看出，波場外推時只有振幅變化，而無相位變化。當(dāng)指數(shù)項取負(fù)號時，外推的波場迅速衰減，稱這種波為倏逝波。當(dāng)指數(shù)項取正號時，外推波場迅速增大，這是一種實際不存在的波，只是進(jìn)行波場計算時發(fā)生，我們稱它為耗損波。在計算中要避免發(fā)生這種情況?，F(xiàn)在我們來看一看上行波和下行波外推時產(chǎn)生波場損耗增大的情況。當(dāng)kx k時，上行波的外推式可寫為：~z)2k2zu(zkx~e（1-56）u(z)從（1-56）式可以看出，此時反向外推遇到倏逝波，正向外推發(fā)生耗損波。分別表示為：~z)~u(zu(z)e~~z)eu(z)u(z

kx2 k2 zk2x k2 z

1-57）1-58）由此可以得出結(jié)論，用上行波方程進(jìn)行向下波場外推永遠(yuǎn)是計算穩(wěn)定的。而用上行波方程進(jìn)行正向外推就可能遇到耗損波，因此有可能是不穩(wěn)定的。除非在計算中不斷地把kxk的波場濾除掉。同理可求出kxk時下行波的外推式為：~d(zz)ekx2k2z（1-59）~d(z)此時也是反向外推遇到倏逝波，正向外推遇到耗損波。上面是在波數(shù)-頻率域進(jìn)行的波場外推的討論。很容易把討論的結(jié)論轉(zhuǎn)用于時間-空間域，這只要把（1.3.4）式進(jìn)行某種展開，代入外推式中，并經(jīng)過Fourier反變換回到時間-空間域中即可得到相應(yīng)的外推方程。除了用波動方程分解后的上行波和下行波方程進(jìn)行波場的正、反向外推之外，還可以用Kirchhoff積分公式進(jìn)行波場的正反向外推。三．地震反射波場成像在這里我們將從波動場的觀點來敘述反射波成像的一般原理。關(guān)于反射波成像的概念和幾何做圖的各種方法請見參考文獻(xiàn) [1]。在介質(zhì)的速度參數(shù)為已知的條件下，確定反射圖象的任務(wù)就是求反射點的空間位置及其反射系數(shù)。由于我們現(xiàn)在還無法求出確切的反射系數(shù)值，成像的反射系數(shù)這一要求實際上是用能反映該反射點反射系數(shù)相對值的反射波振幅來表示的。因此，在目前階段，反射成像實際上就是把地面上觀測到的反射波歸位到產(chǎn)生它的反射點上去。能做到這一點就算實現(xiàn)了成像。這實際就是地震偏移問題。因此地震偏移與地震成像在現(xiàn)階段可以視為同一概念。為了實現(xiàn)地震偏移成像， 首先要進(jìn)行上行波場的反向外推。 外推后求出的各點波場值，有的是來自本點的反射波，有的是該點下方許多點上的反射波。因此，要在外推波場中提取成像值。Claerbout提出下述反射波成像原則：反射面位于這些點上，其入射波的初至與反射波的產(chǎn)生時間相同。如圖1-7所示，入射波從O點發(fā)射，反射波可設(shè)想從虛源O*點出射。在p1點上入射波先到達(dá)，反射波后到達(dá)。在p3點上設(shè)想的反射波先到達(dá)，入射波后到達(dá)。在這兩點上反射波與入射波不同時，因此它們都不是反射圖象的位置。只有在p2點上入射波和反射波才同時到達(dá)，而該點正好位于反射界面上。在這點圖1-7入射波與反射波上的反射系數(shù)為反射波振幅u(x,z,td)除以入射波的振幅d(x,z,td)。因此，反射波成像的基本公式可寫為：Map(x,z) u(x,z,td)d(x,z,td) （1-60）1-60）式?jīng)]有考慮反射系數(shù)隨著入射角變化的情況，它實質(zhì)上是相位信息的公式?；蛘哒f，它對接近法線入射的情況時基本是正確的，能夠反映反射系數(shù)在各點上的變化情況。應(yīng)用（1-60）式涉及到要選擇下行波的初始時間。這是一個困難問題。我們通過假設(shè)下行波是最小相位而避開這個問題。我們把td作為初始時間，則上行波和下行波對這個時間寫成z變換形式為：U(Z)u0u1Zu2Z2（1-61a）D(Z)d0d1Zd2Z2（1-61b）如果D（Z）是最小相位序列，則圖象函數(shù)可以表示為：Map(x,z)U(z)D(z)d（1-62）把1/D(z)展開，我們得到：Map(x,z)u01f1(ui,di)zf2(ui,di)z2dd0式中f1,f2為系數(shù)。由于z的各冪次項在單位圓上擺動，取0~2的積分限時它們等于零。在不考慮一個常系數(shù)的情況下上式等于：Map(x,z)u0d0（1-63）由于（1-60）式本質(zhì)上也是等于u0d0。因此，（1-62）式與（1-60）式等價。（1-62）式比（1-60）式有優(yōu)點。前者不需要考慮拾取發(fā)生的時間，而是使用u和d的全頻譜。但是要求下行波是最小相位的。我們把（1-62）式的分子和分母都乘上下行波的復(fù)共軛函數(shù)D*，則有：Map(x,z)U(x,z,)D*(x,z,)d（1-64）D(x,z,)D*(x,z,)考慮到分母中的譜密度DD*不含有相位信息，因此反射圖象公式可寫為：Map(x,z)U(x,z,)D*(x,z,)d（1-65）根據(jù)Parseval公式，（1-65）右端項有：1 U(x,z, )D*(x,z, )d u(x,z,t)d*(x,z,t)dt2考慮到下行波是實數(shù)序列，并且不考慮積分號前的常系數(shù)，則反射圖象可寫為：Map(x,z) u(x,z,t)d(x,z,t)dt （1-66）當(dāng)下行波是脈沖波時， （1-66）式是很精確的。但是，如果 d(x,z,t)是一個短延續(xù)長度的子波時，它只是一個很好的近似成像公式。為了驗證上面反射成像公式的正確性，我們著重從相位角度來討論它。為此取任意反射界面 f(x,z)，或可寫為：Zf(x)（1-67）取一單震源在O點激發(fā)（圖1-8），設(shè)激發(fā)波為脈沖，則從O點向界面發(fā)射的入射波顯然可寫為（設(shè)介質(zhì)是均勻的，速度為v）：d(x,z,t)A0(tx2z2v)（1-68）這個任意界面的反射波方程是很復(fù)雜的。為了求出反射波的方程，我們用射線法的虛震源來導(dǎo)出。在彎曲界面的情況下虛震源是變動的。在二維情況下它是一條曲線。取Z f(x)界面上任意點 R的導(dǎo)數(shù)為：dz（1-69）f'(xR)dxR求其切線方程為：zzRf'(x)(xxR)（1-70）通過O點做垂直切線的直線方程為：xz （1-71）f'(x)聯(lián)立解（1-70）和（1-71）式，則得二直線的交點為：

圖1-8任意反射界面的入射波和反射波xcf'2(xR)xRf(xR)f'(xR)1f'2(xR)（1-72a）f(xR)f'(xR)xR（1-72b）zc2(xR)1f'由此求出虛震源 O*(x0,z0)的坐標(biāo)方程為：x02xcu(xR)（1-73a）z02zcw(xR)（1-73b）現(xiàn)在我們來研究Zf(x)界面的反射波的傳播方程。為此，首先要求出反射波的波前面的表達(dá)式，即要求出下列方程的包絡(luò)：xu(xR)2zw(xR)2v2t2（1-74）把（1-74）式對xR微分，得xu(xR)u'(xR)zw(xR)w'(xR)0（1-75）聯(lián)立解（1-74）和（1-75）方程可求出包絡(luò)，其參數(shù)式為：xu(xR)vtw'(xR)（1-76a）u'2(xR)w'2(xR)zw(xR)vtu'(xR)（1-76b）u'2(xR)w'2(xR)由此求出反射波波前面的方程為（式中一撇表示導(dǎo)數(shù)）：tpxu(xR)r(xR)zw(xR)r(xR)vw'(xR)（1-77）vu'(xR)式中：r(xR)u'2'2(xR)w(xR)反射成像式為：Map(x,z)u(x,z,ttp)(tx2z2v)dt（1-78）從而求得成像位置的公式為：tpx2z2v（1-79）把（1-77）式代入（1-79）式，經(jīng)過整理后，得v2(x)x2z22r2(x)xz2u(x)1v2(x)w(x)r2(x)xv2(x)u'(x)w'(x)u'(x)w'(x)2w(x)11u(x)r2(x)zu2(x)1v2(x)（1-80）v2(x)u'(x)w'(x)w2(x)1v102(x)式中：v(x)u'(x)，r(x)u'2'2w'(x)(x)w(x)為了直觀驗證（ 1-80）式，我們設(shè)任意界面為：z f(x) h則：f'(x) 0將它代入（1-80）式，求出界面成像位置方程為：4hz 4h2 0z h與原設(shè)界面方程一致。這說明反射成像的普遍公式（ 1-66）是正確的。以上是反射波成像一般原理的討論。在地震勘探中常常要對水平疊加的地震記錄進(jìn)行偏移成像處理。水平疊加剖面可以看做是沿測線單點自激自收的觀測剖面。在水平層均勻介質(zhì)覆蓋的情況下水平疊加地震剖面和沿測線進(jìn)行的單點自激自收的剖面完全等價。 在非均勻介質(zhì)和非水平層界面的情況下水平疊加剖面與自激自收剖面只是近似相等。 由于經(jīng)濟(jì)和效率的原因，地震勘探從來沒有使用過自激自收觀測方式。所以常常就把水平疊加地震剖面看做是自激自收的地震剖面。自激自收剖面的特點就是各個點上的入射波和反射波沿著同一條法線射線傳播的。因此，在地面上觀測到的反射波可以認(rèn)為是在反射界面上同時激發(fā)地震波在地面上觀測的結(jié)果。也就是自激自收地震剖面又與在反射界面上同時爆炸產(chǎn)生地震波，并以半速度向外傳播，在地面上觀測到的上行波是等價的。這就是 Loewenthal等人首先提出來的爆炸反射面的概念（圖 1-9）。這個概念對于理解水平疊加剖面的偏移成像是很重要的。因為它比較直觀地說明了這種剖面的成像原理。比前面所述的反射波成像的一般原理要容易理解得多。在這里我們引入了兩種等價。一個是水平疊加剖面等價于自激自收剖面。另一個是自激自 圖1-9爆炸反射面的概念收剖面等價于界面上同時激發(fā)在地面記錄的上行波剖面。這種等價只是概念上的，實際上只有一種水平疊加剖面，并沒有三種剖面?！?.4波場外推的Kirchhoff 積分法Kirchhoff 積分法并不直接解波動方程， 而是用數(shù)學(xué)方法來描述關(guān)于波的傳播的惠更斯原理，從而求出空間上任一點波場值的。而這個波場值正好滿足波動方程。Kirchhoff積分最初是為了描述波場從一個波前向傳播方向上任何一點傳播結(jié)果而導(dǎo)出的，它描述的正是一個實際物理過程，這也正是我們所說的正向外推，這樣的正向外推的Kirchhoff公式對任何方向傳播的波都是適用的。Kirchhoff早在1883年就證明了，從擾動區(qū)向外某點M(x1,y1,z1)傳播的波的

t時刻的波場

u(x1,y1

,z1,t)

，可以從擾動區(qū)封閉表面上的波場ru(x,y,z,t )以及該波場對時間和表面法線方向的導(dǎo)數(shù)通過積分式求出來。因此要假定cu(x,y,z,t)在封閉面上和封閉面內(nèi)有直至二階導(dǎo)數(shù)的連續(xù)性。Kirchhoff 利用了格林定理：(u2vv2u)dVuvvudS（1-81）VSnv式中u，v為具有直至二階導(dǎo)數(shù)的連續(xù)性的函數(shù)。在我們所討論的情況下，u取為波場函數(shù)，v1。其中R2(xx1)2(yy1)2(zz1)2。R1-10（a）,（b）那樣的封當(dāng)把觀測點用包含有波前面在內(nèi)的封閉曲面包圍起來，如圖閉時，這樣的封閉面S和它所包圍的體積V作為（1-81）式的積分限，經(jīng)過一定的推導(dǎo)后得出M(x,y,z)點的正向外推波場為：u(x1,y1,z1,t)11Ruu11u4SvRntnRRdS（1-82）n這里n的方向取封閉表面的外法線方向。如果把觀測點M移至封閉面外，則有：u(x1,y1,z1,t)0（1-83）（1-82）式中u(,,,R)，v為傳v播速度。在實用中所取的計算封閉面一定要把需計算波場的那一點或所有計算點包圍在圖1-10計算正向外推波場的封閉面內(nèi)。同時，計算點M一定要滿足處在波前推進(jìn)方向的前方這一物理要求。否則，既使能夠計算出某個數(shù)值，也不代表那一點的波場。1-82）式就是著名的Kirchhoff積分。它描述了物理波場傳播的過程，也滿足齊次波動方程，是它的積分形式解。對我們來說，也可以稱它為正向外推公式，它既可以用于上行波的正向外推，也可以用于下行波的正向外推。只是外推方向由物理條件而定。這一點與前一節(jié)把波動方程分解為上行波方程和下行波方程后它們的外推公式在形式上是相同的情況完全符合。還要指出，（1-82）式的v是介質(zhì)的傳播速度，在積分式的推導(dǎo)中令其為常數(shù)。因此Kirchhoff積分也只滿足均勻介質(zhì)的情況。下面討論用 Kirchhoff 積分進(jìn)行波場反向外推問題。對偏移來說，需要利用 Kirchhoff積分進(jìn)行反向外推。 這時，我們所取的封閉體積 V應(yīng)當(dāng)在波前傳播方向的反方向， 計算點(x1,y2,z3)就在這個封閉體內(nèi)。根據(jù)格林定理同樣可求出形式上相同的反向外推的Kirchhoff積分式：11Ruu11uu(x1,y1,z1,t)vRntnRRdS'（1-84）4S'n式中的[[u]]不再是推遲場，而是超前場(,,,R)。uxyztv（1-84）式為用于波場反向外推的Kirchhoff積分式。它可用于上行波的反向外推，也可用于下行波的反向外推。當(dāng)然，這種外推與正向外推不同，它不代表一個物理過程，而只是一種重建波場的計算過程。1.5Born近似和Rytov近似本節(jié)首先討論 Born近似，然后討論 Rytov基本特點。為以下討論的方便，給出如下的目標(biāo)函數(shù)

近似，最后討論

Born

近似和

Rytov

近似的O(r)v01（1-85）v(r)則有uk02u2k02O(r)u(r)k02O2(r)u(r)（1-86）在波場弱散射的假設(shè)條件下討論Born近似和Rytov近似表達(dá)式。一．Born近似設(shè)波函數(shù)具有形式解uum（1-87）m0其中，um關(guān)于O(r)是m階的，把式（1-87）代入式（1-86）比較O(r)的階數(shù)，得um滿足方程組u0k02u00（1-88）uk2u12k02O(r)u0（1-89）10222（1-90）um1k0um12k0O(r)umk0O2(r)um1m=1,2,3,。從上面的方程組可得各階Born近似u1(r)2k02drG(rr)O(r)u0(r)（1-91）式中，u0為入射波場，um1(r)2k02drG(rr)O(r)um(r)k02drG(rr)O2(r)um1(r)（1-92）特別是u2(r)2k02drG(rr)O(r)u1(r)k02drG(rr)O(r)u0(r)（1-93）二．Rytov近似設(shè)波函數(shù)ue，則方程（1-86）可寫為(222（1-94）)2k02k0O(r)k0O2(r)設(shè)相位有級數(shù)解，即m（1-95）m0其中，m關(guān)于O(r)是m階的。將式（1-95）代入式（1-94）。比較O(r)的階數(shù)，得到012

(0)2k0200（1-96）22k020（1-97）01O(r)202(1)2k02O2(r)0（1-98）這里，方程（ 1-96）對應(yīng)方程（ 1-94）在O（r）=0時的情形，即無擾動情形。令u0 e0，則入射波u0滿足u0 k02u0 0 （1-99）設(shè)入射波u0為u0e0，0ikr（1-100）利用方程（1-97）、（1-99）和自由空間的 Green函數(shù)G，有2k02r)O(r)u0(r)（1-101）1drG(ru0(r)同理，利用式（1-98）和（1-100），得到21drG(rr)[k02O2(r)(1)2]u0(r)（1-102）u0(r)一般地，對m3，有1m1mdrG(rr)[pmp]u0(r)（1-103）u0(r)p1其中，m稱為第m階的Rytov近似。在上面表達(dá)式的推導(dǎo)過程中，利用式（1-100），則關(guān)于的一階Rytov近似的方程式（1-97）就可以成為(1u0)22（1-104）k01u02k0O(r)u0比較式（1-104）和式（1-89），除相差一個因子u0外，方程形式完全相同，但這并不意味著Born近似和Rytov近似是等同的。應(yīng)該指出，Born近似利用的是波場的振幅，而Rytov近似則利用的是波場的相位。利用波場相位的Rytov近似在某些假定之下可與旅行時反演相比較。三．Born近似和Rytov近似的基本特點在研究Born近似時，假設(shè)總波場是由入射波場u0與散射波場us疊加而成，即uu0us（1-105）這時方程（1-85）變?yōu)閡sk02us2k02O(r)u(r)k02O2(r)u(r)（1-106）利用自由空間的Green函數(shù)G，得到散射波場us為us(r)2k02drG(rr)O(r)(u0us)（1-107）Born近似是在假定usu0之下，在目標(biāo)函數(shù)O(r)內(nèi)測不到散射場，用u0代替u u0 us，即us(r)2drG(rr)O(r)u0(r)（1-108）2k0但是當(dāng)目標(biāo)O(r)的尺度較大時，目標(biāo)O(r)內(nèi)的波場uu0us是不能簡單地用入射波u0來代替的。如假定目標(biāo)O(r)是半徑為a的圓柱體，在n(v0v1)/v1條件下，沿k方向傳播的平面波u0eikr在圓柱體內(nèi)部并不等于入射波u0，特別在通過圓柱體對稱軸的、柱體內(nèi)部的路徑上的波場為u0ei(1n)kr，尤其是通過圓柱體的波場的相位變化是2an4na（1-109）v0其中圓頻率2/，為波長。Born近似有效的一個必要條件是入射波通過目標(biāo)時的波長變化小于，即4na/或na1（1-110）4式（1-110）表明，Born近似的有效區(qū)域依賴于目標(biāo)尺度的大小。在目標(biāo)尺度較小時，Born近似有效。關(guān)于Rytov近似，令相位0s（1-111）則式（1-97）變成s20s(s)22k02O(r)k02O2(r)0（1-112）求解方程（1-112），得s1drG(rr)[2k02O(r)k02O2(r)()2]u0(r)（1-113）u0Rytov近似是在假定(s)222(r)22（1-114）koO2koO(r)2koO(r)下進(jìn)行的，這就要求( s)2 ko2n即n1（1-115）(s)242式（1-115）表明，Rytov近似對目標(biāo)尺度的依賴性不像 Born近似那樣緊密； Rytov近似的有效性依賴于波長下復(fù)相位的變化情況。換言之， s越小，Rytov近似效果越好，這個特點是Born近似所沒有的?！?.6DeWolf近似、薄板近似和屏近似基于DeWolf近似和相屏單程傳播，我們可以在地面反射地震觀測中計算一次反射波或反向散射波的合成地震記錄。具體實現(xiàn)時一般有兩種不同的算法：一種是直接利用多次正向散射-單次反向散射（ MFSB）近似，其中計算反向散射波的步長等于網(wǎng)格間隔，但計算正向傳播的步長使用遠(yuǎn)遠(yuǎn)大于網(wǎng)格間隔的屏間隔；另一種是對正向和反向散射波的屏近似，這樣兩種步長都能采用屏間隔，從而大大提高了計算效率。在復(fù)雜不均勻介質(zhì)（尤其是3-D介質(zhì)）中，快速模擬方法和算法是發(fā)展復(fù)雜構(gòu)造地震成像、反演和解釋方法的關(guān)鍵。有限差分和有限元法從原理上能模擬在任意不均勻介質(zhì)中波的傳播，但計算工作量大、效率低。在大型3-D彈性波問題的情況下更是不可想象。下面給出了一種基于多次正向散射-單次反向散射（Multiple-forescatteringsingle-backscatteringMFSB）近似，即DeWolf近似，在地面反射地震觀測中計算反向散射場的方法。為快速實現(xiàn)該方法，我們給出了雙域算法。當(dāng)不均勻體尺度大于地震波主波長時，這一理論可進(jìn)一步通過薄板近似和屏散射近似得到合理近似和簡化。薄板近似和屏近似能減少計算工作量、提高計算效率。當(dāng)介質(zhì)體內(nèi)的不連續(xù)性不是很強(qiáng)或不均勻體的參數(shù)擾動不太大時，在不均勻體和諧振散射之間的混響通?？梢院雎圆挥?。但正向散射的累加效應(yīng)一般不能忽略掉。實際上，對于大型不均勻體介質(zhì)或長傳播距離，多次散射對正演模擬和反問題來講是很重要的。在這種情況下，Born近似不滿足，但可應(yīng)用 DeWolf近似。一．Lipmann-Schwinger 方程從標(biāo)量波動方程出發(fā)2( 2 )u(x) 0 （1-116）v2(x)定義F(x)v021s2