摘 要： 本文實(shí)現(xiàn)了時(shí)間—空間域和頻率—波數(shù)域的地震波場(chǎng)外推，在時(shí)間—空間域?qū)崿F(xiàn)了65度波動(dòng)方程有限差分偏移，在頻率—波數(shù)域?qū)崿F(xiàn)了只適合地震波速度，僅隨深度變化的相移法偏移和適應(yīng)地震波速度橫向變化的分裂步有限差分偏移。
　　關(guān)鍵詞： 65度方程 有限差分 相移法偏移 分裂步傅立葉偏移（SSF）
　　
　　一、引言
　　1972年克萊鮑特（J.Clearbout）首先提出波動(dòng)方程有限差分偏移法，給出了至今仍然廣泛應(yīng)用的15度方程。雖然這種偏移技術(shù)只能對(duì)傾角很小地質(zhì)構(gòu)造進(jìn)行成像，但是我們今天所用的各種偏移技術(shù)都借鑒了克萊鮑特對(duì)于波動(dòng)方程的處理方法。1978年蓋茲達(dá)格（Gazdag）提出相移法偏移，同樣是在頻率—波數(shù)域偏移，相移法能夠適應(yīng)地震波速度隨深度變化的情況，也是到目前為止速度僅隨深度變化地質(zhì)模型的最佳偏移方法，其計(jì)算量遠(yuǎn)小于波動(dòng)方程有限差分偏移技術(shù)，沒有傾角限制，沒有假頻的危險(xiǎn)，可是這種方法還不能適應(yīng)地震波速度的橫向變化。1990年Stoffa提出分裂步傅立葉變換法（SSF）偏移，也叫相移加校正偏移，對(duì)相移法進(jìn)行改進(jìn)適應(yīng)地震波速度橫向變化，在頻率域進(jìn)行相位移動(dòng)以后，在頻率—空間域再做一次相位校正。這種方法沒有傾角限制，計(jì)算量也較小，但是它要求速度函數(shù)是平滑的（一階可導(dǎo)），不能適應(yīng)地震波速度橫向變化劇烈的地區(qū)。
　　我從地震波動(dòng)理論出發(fā)，根據(jù)前人所做的工作，分別在時(shí)間—空間域和頻率—波數(shù)域，實(shí)現(xiàn)了單程波動(dòng)方程地震波場(chǎng)延拓。運(yùn)用爆炸反射模型的概念，在時(shí)間—空間域和頻率—空間域?qū)崿F(xiàn)了零偏移距（疊后）65度波動(dòng)方程有限差分偏移。在頻率—波數(shù)域?qū)崿F(xiàn)了相移法偏移，進(jìn)一步的，根據(jù)Stoffa的推導(dǎo)實(shí)現(xiàn)了分裂步傅立葉（SSF）偏移以適應(yīng)地震波速度的橫向變化。在此基礎(chǔ)上，設(shè)計(jì)了向斜構(gòu)造地層模型來(lái)驗(yàn)證相移法偏移。
　　二、理論基礎(chǔ)與程序?qū)崿F(xiàn)方法
　　1.爆炸反射面成像原理
　　圖1顯示了兩種波動(dòng)傳播情形，第一種是假象的地下反射層突然爆炸的實(shí)驗(yàn)，地表布置一排假想的檢波器，記錄從假想的爆炸反射面向上傳播到地表的波動(dòng)。
　　圖1在地表所有的點(diǎn)布置炮—檢對(duì)觀測(cè)回聲信號(hào)（左邊）與“爆炸反射面”概念模型（右邊）
　　注意圖中野外記錄情況下的射線路徑與爆炸反射面情形的射線路徑看起來(lái)是一樣的。想象這兩種波場(chǎng)（觀測(cè)波場(chǎng)和假想波場(chǎng)）是一樣的，這將給我們帶來(lái)概念上的方便。如果它們是相同的，我們就可以忽略數(shù)以千計(jì)的業(yè)已做過的試驗(yàn)，而把注意力集中到一個(gè)假想的試驗(yàn)上。兩種情況的一個(gè)明顯的區(qū)別是野外觀測(cè)系統(tǒng)中的波必須先向下傳播然后又沿相同的路經(jīng)向上傳播，而假想實(shí)驗(yàn)只需要向上傳播。野外實(shí)驗(yàn)的旅行時(shí)除以2，實(shí)際上，通常在分析野外實(shí)驗(yàn)數(shù)據(jù)（雙程走時(shí)）時(shí)，都假設(shè)其聲波速度為實(shí)際速度的一半。類比零偏移數(shù)據(jù)和爆炸反射模型，我們用爆炸反射模型零時(shí)刻的波場(chǎng)來(lái)表示地下地層的構(gòu)造圖像。
　　2.波動(dòng)方程有限差分偏移
　　根據(jù)Clerbout的15度近似單程波方程
　　+=0 （1）
　　利用地面記錄的零偏移距地震波場(chǎng)外推得到地下任意點(diǎn)的地震波場(chǎng)，根據(jù)爆炸反射面概念，認(rèn)為零時(shí)刻地下任意點(diǎn)的波場(chǎng)振幅即為該點(diǎn)的反射強(qiáng)度，即為地下反射界面圖像。
　　將方程（1）進(jìn)行差分離散，并根據(jù)初始條件和邊界條件：
　　p=φ，z=0的邊界條件 （2）
　　并且，p=p=0側(cè)邊界條件（可以用吸收邊界等替代） （3）
　　還有初始條件p=0 （4）
　　我們從j=1，n=N開始，考慮所有的1≤i≤I，以及側(cè)邊界條件（3），就可以形成一個(gè)三對(duì)角方程組：
　　1+2a-a0… 0 0-a1+2a-a … 000 -a1+2a …0 0M MMO M M000… 1+2a -a000… -a 1+2a pppMpp
　　=cccMcc（5）
　　方程（5）可以通過追趕法求解，得到p，1，2，…，I。對(duì)所有j和n循環(huán)可求出所有p，i=1，2，…，I；j=1，2，…，J；n=1，2，…，N，應(yīng)用成像條件t′=τ，p即為最終輸出圖像。
　　有限差分45度和65度方程偏移法：
　　考慮（10）式時(shí)我們不是直接把?鄣p/?鄣z(mì)′忽略，而是進(jìn)一步對(duì)z′求導(dǎo)，利用（10）式消去?鄣p/?鄣z(mì)′，然后略去?鄣p/?鄣z(mì)′得到45度方程:
　　-+=0 （6）
　　用系數(shù)優(yōu)化的方法把上式中的1/4和1/2分別由0.3767和0.4761代替可以得到65度方程，它們的差分求解方法類似于15度方程。
　　3.速度隨深度變化的相移法偏移
　　相移法用exp（ikz）直接進(jìn)行向下外推，然后估算t=0時(shí)（反射面在t=0時(shí)激發(fā)）的波場(chǎng)。先對(duì)時(shí)間剖面進(jìn)行二維傅立葉變換，然后把所有在（k，ω）平面內(nèi)的變換后數(shù)據(jù)值乘以下式進(jìn)行相位移：
　　е=exp{-i［1-（）］△z} （7）
　　總結(jié)起來(lái)，相移法的計(jì)算步驟為：
　?。?）對(duì)輸入疊加剖面進(jìn)行二維傅立葉變化，把p（x，z，t）變換為p（k，0，ω）；
　　（2）對(duì)變換域波場(chǎng)進(jìn)行相位移，即乘以因子C，即：
　　p（k，j△z，ω）=C×p（k，（j-1）△z，ω），j=1，2，…，J；
　?。?）對(duì)所有頻率ω求和，得到p（k，j△z，t=0）=p（k，j△z，ω），j=1，2，…，J；
　?。?）沿k進(jìn)行反傅立葉變換，得到p（x，z）即為地層映像。
　　4.分裂步傅立葉偏移（SSF）
　　Stoffa把慢速度場(chǎng)（速度的倒數(shù)）分為兩部分，一部分稱為參考慢度u，另一部分稱為慢度擾動(dòng)△u=u（，z）。與速度隨深度變化的相移法一樣，先用參考慢速做一次相移，再在頻率空間域在每一空間位置用慢速擾動(dòng)做一次相位校正，從而適應(yīng)地震波速的橫向變化。具體偏移步驟為：
　?。?）首先把地面波場(chǎng)變換到頻率-波數(shù)域p（，z=0，ω）；
　?。?）用參考慢度場(chǎng)把地面波場(chǎng)向下延拓，即p（，z+△z，ω）=p（，z，ω）е其中k=ω；
　　（3）把延拓后的波場(chǎng)變換到頻率空間域p（，z+△z，ω）；
　?。?）根據(jù)慢度擾動(dòng)△u再次相移，即p（，z+△z，ω）=p（，z+△z，ω）е；
　?。?）對(duì)所有的頻率求和得到地下任意位置的地層圖像p（，z）。