劉 一，劉 財(cái)，劉 洋，勾福巖，李炳秀

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

0 引言

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

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

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

近年來(lái)，很多地球物理學(xué)家從不同的角度對(duì)多次波的問(wèn)題進(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ù)的動(dòng)力學(xué)可預(yù)測(cè)特性，將多次波中蘊(yùn)含的有效信息提取出來(lái)。利用多次波的信息形成虛擬一次波，能夠?yàn)榛诹黝A(yù)測(cè)濾波的插值方法提供濾波器系數(shù)估計(jì)的數(shù)據(jù)來(lái)源。

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

1 理論基礎(chǔ)

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

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

(1)

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

dn+1-dTa=r≈0。

(2)

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

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

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

(3)

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

dn+1=r+dTa。

(4)

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

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

(5)

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

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

(6)

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

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

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

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

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

(7)

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

(8)

(9)

2 理論模型測(cè)試

首先選取Abma等[20]提出的平面波模型數(shù)據(jù)進(jìn)行測(cè)試，如圖2所示。該模型提供了兩個(gè)數(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ù)，時(shí)間長(zhǎng)度為1 s，共計(jì)120道。圖2b為缺失數(shù)據(jù)，在第20道的位置出現(xiàn)較嚴(yán)重的缺失(15道)，在第60道位置出現(xiàn)較小的缺失(5道)。圖2c為直接使用圖2b進(jìn)行流預(yù)測(cè)濾波器系數(shù)估計(jì)、再對(duì)缺失位置進(jìn)行數(shù)據(jù)重建的結(jié)果，可見(jiàn)在缺失較嚴(yán)重的位置并不能得到理想的重建結(jié)果。圖2d為利用偏差數(shù)據(jù)進(jìn)行濾波系數(shù)估計(jì)、再對(duì)缺失數(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模型測(cè)試Fig.2 Abma model test

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

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

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

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

選取墨西哥灣某地區(qū)的三維地震數(shù)據(jù)進(jìn)一步測(cè)試本文方法的穩(wěn)定性。圖6a為隨機(jī)缺失40%的炮集數(shù)據(jù)體，上部為3.8 s的時(shí)間切片，左下部為距離位置在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ò)多次波與一次波的互相關(guān)(公式(7))構(gòu)建虛擬一次波(圖6b)。結(jié)果雖然出現(xiàn)了一定的隨機(jī)噪聲， 但主要的同相軸信息得到了恢復(fù)，可以為本文的插值方法提供合理的濾波估計(jì)。

利用本文方法對(duì)三維地震缺失道進(jìn)行數(shù)據(jù)重建。應(yīng)用流預(yù)測(cè)濾波器對(duì)構(gòu)建的虛擬一次波進(jìn)行濾波估計(jì)，濾波器尺寸參數(shù)為25(時(shí)間)×5(空間)×10(空間)，通過(guò)計(jì)算公式(9)對(duì)缺失數(shù)據(jù)進(jìn)行重建，結(jié)果如圖7所示。淺層的海底反射以及非平穩(wěn)的同相軸信息都得到了合理恢復(fù)，隨機(jī)噪聲并沒(méi)有對(duì)重建結(jié)果產(chǎn)生明顯的影響，在圖中難以區(qū)分插值數(shù)據(jù)與原始數(shù)據(jù)的位置，說(shuō)明插值結(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ù)測(cè)濾波器結(jié)合虛擬一次波的技術(shù)方案對(duì)缺失地震數(shù)據(jù)進(jìn)行重建。在應(yīng)用流預(yù)測(cè)濾波器估計(jì)濾波系數(shù)的過(guò)程中，將相鄰濾波器系數(shù)具有比例相似性作為約束條件可快速有效地完成濾波估計(jì)，并利用數(shù)據(jù)本身多次波所蘊(yùn)含的動(dòng)力學(xué)信息預(yù)測(cè)缺失的地震數(shù)據(jù)，在數(shù)據(jù)插值過(guò)程中更合理地為濾波估計(jì)提供數(shù)據(jù)分析基礎(chǔ)。

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