復(fù)雜地震波場的自適應(yīng)流預(yù)測插值方法

劉 一，劉 財，劉 洋，勾福巖，李炳秀

吉林大學(xué)地球探測科學(xué)與技術(shù)學(xué)院,長春 130026

0 引言

地震數(shù)據(jù)在采集過程中總是表現(xiàn)出空間不連續(xù)性：一方面是由于實(shí)際勘探環(huán)境非常復(fù)雜，如村莊、河流等障礙物的存在造成觀測系統(tǒng)的改變；另一方面是由于經(jīng)濟(jì)因素的制約，造成檢波器和炮點(diǎn)不可能連續(xù)布置，使地震數(shù)據(jù)缺失。但是完整的地震數(shù)據(jù)是許多重要處理方法的前提，如表層相關(guān)多次波消除、波動方程偏移和時移地震等。常規(guī)的地震數(shù)據(jù)插值方法往往對地震數(shù)據(jù)的平穩(wěn)性要求較高，如地震同相軸為平穩(wěn)平面波假設(shè)，但地震數(shù)據(jù)本質(zhì)上是非平穩(wěn)的，開發(fā)能夠有效解決非平穩(wěn)復(fù)雜波場的缺失數(shù)據(jù)重建問題具有重要意義。

預(yù)測濾波器在地震數(shù)據(jù)處理和分析中具有重要的作用，例如在地震數(shù)據(jù)反褶積[1]和壓制隨機(jī)噪聲[2-3]等方面。預(yù)測濾波器可以有效地表征地震數(shù)據(jù)能量譜[4-5]，為解決地球物理反演問題提供了一種有效的近似估計(jì)方法，也為恢復(fù)缺失的地震信息指明了方向。地球物理學(xué)者對預(yù)測濾波器的研究取得了一定的進(jìn)展。預(yù)測濾波器的本質(zhì)是信號自回歸，在時間-空間域和頻率-空間域都具有適用性，而時間-空間域的濾波能夠避免出現(xiàn)假頻現(xiàn)象[6]。可以利用預(yù)測誤差濾波器(由預(yù)測濾波系數(shù)構(gòu)成)的尺度縮放不變性計(jì)算反空間假頻地震傾角模式，對缺失的地震數(shù)據(jù)進(jìn)行反假頻插值[7]。

隨著勘探區(qū)域地下介質(zhì)復(fù)雜性的增加，地震信號的非平穩(wěn)性往往也愈加明顯，如何將傳統(tǒng)的平穩(wěn)預(yù)測濾波處理方法應(yīng)用到非平穩(wěn)數(shù)據(jù)是一個重要的研究方向。常規(guī)的解決方式將數(shù)據(jù)進(jìn)行重疊時窗處理，或者假設(shè)數(shù)據(jù)是局部平穩(wěn)[8-9]。Fomel[10]提出了基于整形正則化算子的非平穩(wěn)自回歸方法，并用于解決地球物理反演問題，其擴(kuò)展可應(yīng)用于時頻分析和非平穩(wěn)多項(xiàng)式擬合。Liu 等[11]提出應(yīng)用非因果正則化非平穩(wěn)自回歸(NRNA)進(jìn)行頻率-空間域的隨機(jī)噪聲壓制，在二維地震數(shù)據(jù)中取得了理想的效果，繼而將濾波系數(shù)在空間方向進(jìn)行擴(kuò)展，開發(fā)了三維非平穩(wěn)自回歸對地震數(shù)據(jù)的隨機(jī)噪聲壓制方法。Liu 等[12]利用整形正則化條件約束預(yù)測誤差濾波系數(shù)，對帶有空間假頻的地震數(shù)據(jù)進(jìn)行了有效的重建。但是這些方法在解決濾波器系數(shù)的時空變化屬性表征問題時都增加了額外的計(jì)算成本。流預(yù)測濾波器[13]可以使濾波系數(shù)隨著數(shù)據(jù)的變化同時更新，運(yùn)算形式僅僅表現(xiàn)為信號的褶積，并且不需多次迭代，在提高運(yùn)算效率的同時降低了內(nèi)存成本。流預(yù)測濾波器對非平穩(wěn)數(shù)據(jù)具有預(yù)測逆運(yùn)算的特點(diǎn)，可以快速地完成缺失數(shù)據(jù)重建。但是，缺失的地震道對濾波器系數(shù)的準(zhǔn)確估計(jì)有較大影響，很難保證準(zhǔn)確的插值效果。

近年來，很多地球物理學(xué)家從不同的角度對多次波的問題進(jìn)行探討和嘗試，并不局限于視其為噪聲從地震數(shù)據(jù)中去除，而是如何更好地利用多次波中所包含的地下反射信息。多次波在偏移成像[14]和聯(lián)合反演[15]中都得到了一定的應(yīng)用。構(gòu)建準(zhǔn)一次波的插值方法[16-17]以及直接利用含有多次波的地震數(shù)據(jù)構(gòu)建擬地震數(shù)據(jù)的插值方法[18]，主要是利用地震數(shù)據(jù)的動力學(xué)可預(yù)測特性，將多次波中蘊(yùn)含的有效信息提取出來。利用多次波的信息形成虛擬一次波，能夠?yàn)榛诹黝A(yù)測濾波的插值方法提供濾波器系數(shù)估計(jì)的數(shù)據(jù)來源。

本文針對流預(yù)測濾波器具有較準(zhǔn)確非平穩(wěn)波場能量譜估計(jì)和較高運(yùn)算效率的特征，通過虛擬一次波為缺失數(shù)據(jù)提供合理的濾波估計(jì)，利用自回歸預(yù)測理論對缺失數(shù)據(jù)進(jìn)行重建，以解決復(fù)雜疊前地震波場的數(shù)據(jù)缺失問題。

1 理論基礎(chǔ)

1.1 預(yù)測濾波插值方法

預(yù)測濾波器基于自回歸理論，利用預(yù)測濾波器進(jìn)行地震數(shù)據(jù)的重建可以表述為兩步求解反問題的過程。第一步，利用最小二乘方法求解反問題，估計(jì)濾波系數(shù)。例如，給定一組數(shù)據(jù)d(其中可能包含缺失數(shù)據(jù))，利用預(yù)測濾波器估計(jì)濾波系數(shù)，當(dāng)dn+1已知時，基于預(yù)測誤差近似為零的自回歸表達(dá)公式為

(1)

式中：ak為預(yù)測濾波系數(shù)；n、k為自然數(shù)。其矩陣表達(dá)形式為

dn+1-dTa=r≈0。

(2)

式中：r為殘差；d為數(shù)據(jù)的矩陣形式；a為預(yù)測濾波系數(shù)的矩陣形式。

利用最小二乘反演估計(jì)濾波系數(shù)，表達(dá)式為

a=(ddT)-1dn+1d。

(3)

第二步，利用求得的濾波系數(shù)對缺失的地震數(shù)據(jù)進(jìn)行恢復(fù)。當(dāng)dn+1為缺失數(shù)據(jù)時，通過自回歸方程再次求解公式(2)，并利用第一步估計(jì)出的預(yù)測濾波系數(shù)可以完成對缺失數(shù)據(jù)的重建，其計(jì)算表達(dá)式為

dn+1=r+dTa。

(4)

由于傳統(tǒng)的預(yù)測濾波方法很難表征復(fù)雜的非平穩(wěn)地震數(shù)據(jù)，直接利用傳統(tǒng)的預(yù)測濾波器進(jìn)行數(shù)據(jù)插值會引起較大的誤差；而自適應(yīng)預(yù)測濾波器[13]的結(jié)構(gòu)造成計(jì)算資源占用率過高，很難滿足實(shí)際生產(chǎn)的需求。因此，利用改進(jìn)的流預(yù)測濾波方法進(jìn)行濾波器系數(shù)的非迭代計(jì)算。

1.2 流預(yù)測濾波器估計(jì)

(5)

powerGUN鋁點(diǎn)焊自動焊鉗采用一個或者并聯(lián)兩個變壓器，用于高電流和高節(jié)拍焊接；增大了次級回路元件截面積；優(yōu)化了冷卻水系統(tǒng)，提高了焊鉗的冷卻散熱能力；高強(qiáng)度鋁合金結(jié)構(gòu)增加了焊鉗的剛性，獨(dú)立平衡補(bǔ)償、伺服電動機(jī)驅(qū)動和壓力傳感器保證焊接質(zhì)量，再配合HWH/NIMAK組合控制柜，讓其更完美地達(dá)到了焊接過程的控制。

在求解方程組的過程中，采用Sherman-Morrison公式[19]進(jìn)行解析表達(dá)式推導(dǎo)，濾波系數(shù)估計(jì)的表達(dá)式為

(6)

新的濾波器系數(shù)通過加上一個按一定尺度縮放的數(shù)據(jù)而不斷更新，實(shí)現(xiàn)自適應(yīng)變化。這個尺度與前一個濾波器系數(shù)和殘差存在比例關(guān)系(殘差由前一個濾波器系數(shù)和數(shù)據(jù)點(diǎn)積(dTa)求得)，公式(6)得到的更新濾波器系數(shù)只進(jìn)行代數(shù)運(yùn)算(向量的點(diǎn)積)而沒有迭代的需求，此時預(yù)測濾波器的非平穩(wěn)數(shù)據(jù)表征問題得到解決。但是在數(shù)據(jù)重建過程中，公式(2)中d所包含的缺失數(shù)據(jù)會造成預(yù)測濾波器系數(shù)的較大偏差，進(jìn)而導(dǎo)致插值精度的降低，直接利用缺失數(shù)據(jù)的濾波系數(shù)估計(jì)可以轉(zhuǎn)換為利用近似數(shù)據(jù)的濾波估計(jì)，只要近似數(shù)據(jù)具有與原始數(shù)據(jù)相似的預(yù)測特征即可。

1.3 基于虛擬一次波的預(yù)測濾波估計(jì)

多次波可以通過地震數(shù)據(jù)本身的褶積預(yù)測出來，將一次波“升階”為多次波。反之，可以利用互相關(guān)代替褶積，將多次波“降階”，利用地震數(shù)據(jù)本身的多次波提取出有效的一次波信息，稱為虛擬一次波。構(gòu)建虛擬一次波的具體實(shí)現(xiàn)算法是，對同一個炮集的原始記錄和多次波記錄進(jìn)行互相關(guān)，一次波的接收點(diǎn)則成為虛擬一次波的震源點(diǎn)，對同一震源進(jìn)行疊加，可以提高虛擬一次波的信噪比，繼而得到完整的道集記錄。圖1給出了構(gòu)建“虛擬一次波”的示意圖。虛擬一次波具有與原始一次波相似的動力學(xué)特征，由于其利用觀測系統(tǒng)的坐標(biāo)關(guān)系進(jìn)行疊加計(jì)算，因此不受缺失數(shù)據(jù)的影響，可以代替缺失數(shù)據(jù)進(jìn)行預(yù)測濾波器的系數(shù)估算。

圖1 虛擬一次波射線軌跡示意圖Fig.1 Raypaths of the pseudoprimaries

在二維觀測系統(tǒng)中(圖1)，利用SR1的一次波數(shù)據(jù)與SR2的多次波數(shù)據(jù)進(jìn)行互相關(guān)運(yùn)算，可以得到R1R2的虛擬一次數(shù)據(jù)，其與原始數(shù)據(jù)中以R1為震源點(diǎn)在R2接收的一次波有著相似的動力學(xué)特征，表達(dá)式為

(7)

虛擬一次波雖然不能完全與原始數(shù)據(jù)相匹配，在振幅、相位、子波等方面存在誤差，但不影響預(yù)測濾波器的估計(jì)。將虛擬一次波作為初始數(shù)據(jù)，利用公式(6)估計(jì)濾波系數(shù)，并將其代入公式(4)求解插值表達(dá)式

(8)

(9)

2 理論模型測試

首先選取Abma等[20]提出的平面波模型數(shù)據(jù)進(jìn)行測試，如圖2所示。該模型提供了兩個數(shù)據(jù)體，分別為標(biāo)準(zhǔn)數(shù)據(jù)和偏差數(shù)據(jù)(即振幅、相位與標(biāo)準(zhǔn)數(shù)據(jù)不匹配)。這里利用該模型驗(yàn)證本文方法的有效性。圖2a為標(biāo)準(zhǔn)數(shù)據(jù)，時間長度為1 s，共計(jì)120道。圖2b為缺失數(shù)據(jù)，在第20道的位置出現(xiàn)較嚴(yán)重的缺失(15道)，在第60道位置出現(xiàn)較小的缺失(5道)。圖2c為直接使用圖2b進(jìn)行流預(yù)測濾波器系數(shù)估計(jì)、再對缺失位置進(jìn)行數(shù)據(jù)重建的結(jié)果，可見在缺失較嚴(yán)重的位置并不能得到理想的重建結(jié)果。圖2d為利用偏差數(shù)據(jù)進(jìn)行濾波系數(shù)估計(jì)、再對缺失數(shù)據(jù)進(jìn)行重建的結(jié)果，圖中同相軸的連續(xù)性和能量得到了比較準(zhǔn)確的恢復(fù)。

a.標(biāo)準(zhǔn)數(shù)據(jù)；b.缺失數(shù)據(jù)；c.使用缺失數(shù)據(jù)濾波系數(shù)的重建結(jié)果；d.使用偏差數(shù)據(jù)濾波系數(shù)的重建結(jié)果。圖2 Abma模型測試Fig.2 Abma model test

為了驗(yàn)證本文方法對復(fù)雜地震數(shù)據(jù)的有效性，接下來測試Sigsbee 2B模型，利用有限差分正演疊前地震數(shù)據(jù)。模型中包含一個不規(guī)則的鹽丘體反射層，使得正演數(shù)據(jù)體現(xiàn)出復(fù)雜的地震波場特征。圖3為Sigsbee 2B的速度模型，正演數(shù)據(jù)時間采樣間隔為8 ms，共計(jì)1 000個采樣點(diǎn)；由于空氣、水和鹽丘體存在較大的波阻抗差，所以多次波發(fā)育明顯，為虛擬一次波的構(gòu)建提供了數(shù)據(jù)保證。圖4a為去除20%地震道的正演記錄，圖4b為利用圖4a構(gòu)建的虛擬一次波數(shù)據(jù)，其與原始數(shù)據(jù)(圖4a)具有相似的動力學(xué)特征。

圖3 Sigsbee 2B速度模型Fig.3 Velocity of Sigsbee2B model

選取測線坐標(biāo)為14 km的完整單炮數(shù)據(jù)(圖5a)測試插值效果，缺失數(shù)據(jù)中有3個較大的數(shù)據(jù)缺口。利用流預(yù)測濾波器對虛擬一次波數(shù)據(jù)進(jìn)行濾波估計(jì)，該方法只需要兩個向量的褶積而不需迭代，大大降低了運(yùn)算成本，根據(jù)公式(9)可以高效地完成缺失數(shù)據(jù)的重建。濾波器的尺寸為20(時間)×5(空間)，近似系數(shù)為1.22。圖5b為插值結(jié)果，結(jié)果表明本文的方法在處理復(fù)雜地震數(shù)據(jù)時，可以快速和準(zhǔn)確地對缺失地震道完成重建工作。

3 實(shí)際數(shù)據(jù)測試

選取墨西哥灣某地區(qū)的三維地震數(shù)據(jù)進(jìn)一步測試本文方法的穩(wěn)定性。圖6a為隨機(jī)缺失40%的炮集數(shù)據(jù)體，上部為3.8 s的時間切片，左下部為距離位置在1.7 km的單炮記錄，右下部為0.5 km的共炮檢距剖面。炮間距為0.05 km，共計(jì)80炮。實(shí)際地震資料中多次波和散射波都有明顯的體現(xiàn)，在數(shù)據(jù)缺失較嚴(yán)重的位置很難判斷同相軸的走向。

在實(shí)際資料的處理中，利用Radon變換[21]或者在逆數(shù)據(jù)域下[22]有效地分離一次波和多次波，再通過多次波與一次波的互相關(guān)(公式(7))構(gòu)建虛擬一次波(圖6b)。結(jié)果雖然出現(xiàn)了一定的隨機(jī)噪聲， 但主要的同相軸信息得到了恢復(fù)，可以為本文的插值方法提供合理的濾波估計(jì)。

利用本文方法對三維地震缺失道進(jìn)行數(shù)據(jù)重建。應(yīng)用流預(yù)測濾波器對構(gòu)建的虛擬一次波進(jìn)行濾波估計(jì)，濾波器尺寸參數(shù)為25(時間)×5(空間)×10(空間)，通過計(jì)算公式(9)對缺失數(shù)據(jù)進(jìn)行重建，結(jié)果如圖7所示。淺層的海底反射以及非平穩(wěn)的同相軸信息都得到了合理恢復(fù)，隨機(jī)噪聲并沒有對重建結(jié)果產(chǎn)生明顯的影響，在圖中難以區(qū)分插值數(shù)據(jù)與原始數(shù)據(jù)的位置，說明插值結(jié)果具有較好的合理性。

圖4 原始缺失數(shù)據(jù)(a)和虛擬一次波數(shù)據(jù)(b)Fig.4 The missing data (a) and the pseudoprimaries data of Sigsbee2B model(b)

a.缺失數(shù)據(jù)單炮記錄；b.本文方法數(shù)據(jù)重建結(jié)果。圖5 Sigsbee2B模型數(shù)據(jù)重建結(jié)果Fig.5 Result of the reconstruction of Sigsbee2B model

圖6 隨機(jī)缺失40%的地震數(shù)據(jù)(a)和構(gòu)建的虛擬一次波(b)Fig.6 Data with 40%randomly missing traces (a) and the pseudoprimaries data (b)

圖7 本文方法數(shù)據(jù)重建結(jié)果Fig.7 Result of the data reconstruction

4 結(jié)論

1)本文提出基于流預(yù)測濾波器結(jié)合虛擬一次波的技術(shù)方案對缺失地震數(shù)據(jù)進(jìn)行重建。在應(yīng)用流預(yù)測濾波器估計(jì)濾波系數(shù)的過程中，將相鄰濾波器系數(shù)具有比例相似性作為約束條件可快速有效地完成濾波估計(jì)，并利用數(shù)據(jù)本身多次波所蘊(yùn)含的動力學(xué)信息預(yù)測缺失的地震數(shù)據(jù)，在數(shù)據(jù)插值過程中更合理地為濾波估計(jì)提供數(shù)據(jù)分析基礎(chǔ)。

2)理論模型測試結(jié)果表明，在缺失較嚴(yán)重的數(shù)據(jù)中，改進(jìn)的流預(yù)測濾波插值方法可以獲得更加理想的數(shù)據(jù)重建效果；實(shí)際地震數(shù)據(jù)處理結(jié)果表明，該方法可以有效地針對復(fù)雜非平穩(wěn)地震波場的缺失問題進(jìn)行求解。