基于壓縮感知的地震觀測(cè)系統(tǒng)設(shè)計(jì)及數(shù)據(jù)重建

崔慶輝， 芮擁軍， 秦寧， 劉鵬翔， 隆文韜

(中國(guó)石油化工股份有限公司勝利油田分公司物探研究院， 東營(yíng) 257022)

隨著社會(huì)經(jīng)濟(jì)的快速發(fā)展，當(dāng)前在東部老油區(qū)城鎮(zhèn)、工廠、村莊、道路等障礙物密集分布，地震勘探的地表?xiàng)l件發(fā)生了巨大變化，造成基于規(guī)則采樣理論即Shannon采樣定理[1]設(shè)計(jì)的觀測(cè)系統(tǒng)受密集障礙物影響，在實(shí)際布設(shè)時(shí)點(diǎn)位發(fā)生嚴(yán)重偏離，特別是炮點(diǎn)位置變化更大，呈現(xiàn)出隨機(jī)分布特點(diǎn)，導(dǎo)致成像結(jié)果偏離預(yù)期。另一方面，一些老油田勘探開發(fā)程度較高，需要解決的地質(zhì)目標(biāo)表現(xiàn)出“深、小、薄、碎”等特征，在常規(guī)規(guī)則采樣理論下只能不斷提高采樣密度即炮道密度，這無疑會(huì)造成采集成本的急劇上升。

壓縮感知理論是一種新的采樣理論，基于信號(hào)的稀疏特性，在遠(yuǎn)小于Nyquist采樣率的條件下，通過隨機(jī)采樣獲取信號(hào)的離散樣點(diǎn)值，然后進(jìn)行信號(hào)的重建[2-3]。壓縮感知理論一經(jīng)提出，就引起學(xué)術(shù)界和工業(yè)界的廣泛關(guān)注，在各個(gè)領(lǐng)域得到迅速發(fā)展。

國(guó)外較早將壓縮感知理論應(yīng)用于海上地震勘探，取得了較好的應(yīng)用效果[4-5]。馬堅(jiān)偉[6-7]較早研究了壓縮感知理論在地震勘探領(lǐng)域應(yīng)用的可行性，周松等[8]和呂公河等[9]在西部沙漠地區(qū)開展了基于壓縮感知的地震采集試驗(yàn)，為中國(guó)陸上基于壓縮感知的地震勘探實(shí)際應(yīng)用奠定了基礎(chǔ)。黃小剛[10-11]提出了一套海上壓縮感知采集設(shè)計(jì)及數(shù)據(jù)重構(gòu)方法，并用模型數(shù)據(jù)展示了海上壓縮感知采集的效果。總之，在當(dāng)前低油價(jià)背景下，壓縮感知理論在地震勘探領(lǐng)域的研究得到了廣泛的重視[12]。

利用壓縮感知理論實(shí)現(xiàn)隨機(jī)地震采樣并重構(gòu)規(guī)則數(shù)據(jù)，首先要分析地震數(shù)據(jù)的稀疏性，地震信號(hào)(或波場(chǎng))的時(shí)空復(fù)雜度(與地下及地表速度結(jié)構(gòu)、噪音特征等相關(guān))會(huì)顯著影響地震信號(hào)的稀疏性，進(jìn)而導(dǎo)致很多情形下隨機(jī)采集的實(shí)際數(shù)據(jù)的重建效果無法達(dá)到規(guī)則采樣數(shù)據(jù)的質(zhì)量[13]。中國(guó)陸上油田地表噪聲發(fā)育，地下波場(chǎng)復(fù)雜，基于壓縮感知理論進(jìn)行隨機(jī)地震數(shù)據(jù)采集及規(guī)則化重建具有更大難度。

首先研究基于壓縮感知理論的隨機(jī)地震觀測(cè)系統(tǒng)設(shè)計(jì)方法及基于地震數(shù)據(jù)高維空間低秩約束的隨機(jī)地震數(shù)據(jù)規(guī)則化重建方法，并通過三維模型數(shù)值模擬證明了以上方法在同等采樣密度下比常規(guī)規(guī)則采樣具有更好的成像效果。

1 隨機(jī)地震觀測(cè)系統(tǒng)設(shè)計(jì)基本理論

隨機(jī)地震采集需要滿足三個(gè)條件：地震數(shù)據(jù)具有稀疏性、合理的隨機(jī)采樣方法和有效的數(shù)據(jù)規(guī)則化重建算法。關(guān)于地震數(shù)據(jù)的稀疏性，很多學(xué)者的研究表明實(shí)際地震數(shù)據(jù)在時(shí)間域或其他域具有稀疏特征，即地震數(shù)據(jù)可稀疏表達(dá)[13-15]?，F(xiàn)重點(diǎn)研究隨機(jī)地震觀測(cè)系統(tǒng)設(shè)計(jì)及數(shù)據(jù)規(guī)則化方法。

1.1 隨機(jī)地震采樣矩陣

野外地震數(shù)據(jù)觀測(cè)實(shí)際上是對(duì)信號(hào)的一種采樣，地震數(shù)據(jù)采集過程可以表示為

d=Df

(1)

式(1)中：f∈Rm為原始地震數(shù)據(jù)；d∈Rn為經(jīng)過采樣得到的地震數(shù)據(jù)；D∈Rn×m為采樣矩陣；通常有n

若f本身不稀疏，則可以使用f在變換域中的稀疏性來得到一種等價(jià)的表示。這里用C表示f的稀疏變換，則式(1)就可以轉(zhuǎn)換為

d=Ak

(2)

A=DCH

(3)

式中：H代表反變換矩陣；k∈Rn為f在稀疏域中的系數(shù)；n為變換域的維度(若C為正交變換，則n=m；若C為冗余變換，n>m)，A∈Rn×n。

要滿足A為隨機(jī)矩陣的要求，在變換矩陣C確定的情況下，采樣矩陣R也必須是隨機(jī)分布矩陣，那么就可以利用采樣數(shù)據(jù)高概率地恢復(fù)完整地震數(shù)據(jù)。

1.2 隨機(jī)地震觀測(cè)系統(tǒng)設(shè)計(jì)方法

通過高斯隨機(jī)函數(shù)進(jìn)行完全的隨機(jī)采樣，被基于壓縮感知的技術(shù)廣泛采用。但是，這種方法容易造成采樣點(diǎn)高度聚集或者高度分散的情況，導(dǎo)致對(duì)信號(hào)的某些部分采樣冗余或者采樣不足，不僅不符合高精度地震數(shù)據(jù)采集的要求，而且也會(huì)給后續(xù)數(shù)據(jù)恢復(fù)處理帶來困難，難以達(dá)到理想的重建效果。

計(jì)算機(jī)圖像學(xué)領(lǐng)域應(yīng)用的泊松碟隨機(jī)采樣方法[16]、抖動(dòng)采樣方法(Jittered)[17]以及分段采樣方法不僅能保證測(cè)點(diǎn)的隨機(jī)性，還能保證測(cè)點(diǎn)的空間均勻分布，并具有明顯的“藍(lán)色噪聲”頻譜特征[14]?！八{(lán)色噪聲”頻譜特征是描述隨機(jī)類型特征的一個(gè)統(tǒng)計(jì)模型，其特點(diǎn)是在主要頻率周圍分布的其他頻率值都非常小，所以真值或者主頻率就很容易被檢測(cè)到。

結(jié)合實(shí)際地震數(shù)據(jù)特點(diǎn)，基于Jittered采樣原理提出了以下隨機(jī)地震觀測(cè)系統(tǒng)設(shè)計(jì)技術(shù)流程。

步驟1按照常規(guī)的三維地震觀測(cè)系統(tǒng)設(shè)計(jì)方法確定理想的觀測(cè)參數(shù)，在正交觀測(cè)系統(tǒng)下可以用炮點(diǎn)和檢波點(diǎn)的點(diǎn)距、線距、點(diǎn)數(shù)和線數(shù)表示，即{x0,y0,Δx,Δy,nx,ny}，其中(x0,y0)為采樣點(diǎn)的起始坐標(biāo)，Δx和Δy分別為采樣點(diǎn)之間的橫向和縱向間隔，nx和ny分別為橫向和縱向采樣點(diǎn)數(shù)。基于以上參數(shù)可得到規(guī)則的采樣點(diǎn)坐標(biāo){x,y}。

步驟2根據(jù)工區(qū)地表及地下復(fù)雜程度設(shè)置縱橫向壓縮采樣比例{rx,ry}，基于該壓縮采樣比例計(jì)算隨機(jī)采樣點(diǎn)坐標(biāo)的期望值{xe,ye}。

步驟3設(shè)置隨機(jī)采樣點(diǎn)坐標(biāo)相對(duì)于期望值的擾動(dòng)量，即高斯隨機(jī)分布的標(biāo)準(zhǔn)差(σx,σy)，一般取(σx,σy)=(3Δx,3Δy)。

最后，根據(jù)步驟2中期望值，以及步驟3中標(biāo)準(zhǔn)差利用隨機(jī)函數(shù)生成隨機(jī)觀測(cè)網(wǎng)格{xrand,yrand}。

2 高維空間低秩約束的隨機(jī)地震數(shù)據(jù)重建

地震數(shù)據(jù)的偏移對(duì)數(shù)據(jù)規(guī)則性要求很高，實(shí)際采集的地震數(shù)據(jù)受觀測(cè)系統(tǒng)、地表變觀等影響需要進(jìn)行規(guī)則化處理，因此產(chǎn)生了很多的地震數(shù)據(jù)插值算法[18-19]，主要包括基于函數(shù)變換的重建、基于波場(chǎng)延拓的重建、低秩約束的地震數(shù)據(jù)重建。低秩約束的地震數(shù)據(jù)重建算法應(yīng)用最為廣泛，其基本原理為：完整地震數(shù)據(jù)在頻率域可以被一個(gè)低秩矩陣或者數(shù)組表示，不規(guī)則缺道數(shù)據(jù)被看作隨機(jī)噪聲，而隨機(jī)噪聲會(huì)增加矩陣或者高維數(shù)組的秩。地震數(shù)據(jù)重建就可以歸結(jié)為一個(gè)高秩矩陣或者數(shù)組的降秩問題，通過降秩實(shí)現(xiàn)重建。

地震數(shù)據(jù)中的反射波同相軸在局部表現(xiàn)為線性特征，利用數(shù)據(jù)的局部相似性，抽取局部五維數(shù)據(jù)，在高維空間進(jìn)行“稀疏-低秩聯(lián)合反演重建”，從而在壓縮感知理論下保證了重建數(shù)據(jù)的精度，算法主要流程如圖1所示。

圖1 低秩約束的高維空間地震數(shù)據(jù)重建算法流程Fig.1 Low-rank constrained high-dimensional spatial seismic data reconstruction algorithm flow

步驟1定義一個(gè)理想的規(guī)則工區(qū)網(wǎng)格G，按照規(guī)則觀測(cè)系統(tǒng)設(shè)計(jì)方法[20]在網(wǎng)格G上定義規(guī)則觀測(cè)系統(tǒng)，記為Ωreg。

步驟2獲取原始輸入數(shù)據(jù)集Dirr的非規(guī)則地震數(shù)據(jù)的觀測(cè)系統(tǒng)。以炮集為單位獲取炮檢點(diǎn)坐標(biāo)，取整到臨近網(wǎng)格，得到Dirr對(duì)應(yīng)的非規(guī)則觀測(cè)系統(tǒng)Ωirr。

步驟3以窗口滑動(dòng)方式從原始數(shù)據(jù)集取得一個(gè)局部高維數(shù)據(jù)子集。通過在炮端檢端分別取局部窗口，抽取原始道集中落在窗口中的所有地震道，從而獲得一個(gè)不完備的五維局部數(shù)據(jù)體。通過滑動(dòng)中心炮檢點(diǎn)，可以循環(huán)獲取不同位置的局部高維數(shù)據(jù)。

步驟4在步驟3獲取的不完備五維局部數(shù)據(jù)體上，在頻率域進(jìn)行規(guī)則數(shù)據(jù)重建，輸出一個(gè)完備且規(guī)則的局部五維局部數(shù)據(jù)體。計(jì)算過程即求解優(yōu)化問題，即

(4)

該優(yōu)化問題實(shí)現(xiàn)了對(duì)重建結(jié)果的三種約束。

(1)要求在觀測(cè)點(diǎn)位上的重建結(jié)果與原始觀測(cè)相吻合(L2項(xiàng))。

(2)要求重建數(shù)據(jù)在四維波數(shù)具有稀疏性(L1項(xiàng))。

(3)要求重建數(shù)據(jù)在進(jìn)行Hankel排列后具有低秩性(核范數(shù)項(xiàng))。

式(4)所示優(yōu)化問題的求解，可通過以下迭代算法實(shí)現(xiàn)。

(1)關(guān)于L2項(xiàng)的梯度更新為

(5)

(6)

(7)

(2)關(guān)于L1項(xiàng)的軟閾值收縮為

(8)

(3)關(guān)于核范數(shù)項(xiàng)的奇異值收縮：

(9)

式(9)中：H(·)和H-1(·)分別是Hankel矩陣化正、反變換算子；Sμ(x)表示對(duì)x進(jìn)行奇異值軟閾值收縮，閾值取μ。

步驟5對(duì)步驟4重建的局部五維規(guī)則數(shù)據(jù)體做余弦窗函數(shù)加權(quán)，然后以累加方式更新到規(guī)則數(shù)據(jù)集Dreg中，并記錄輸出數(shù)據(jù)集中各個(gè)地震道的累加總權(quán)重Csum(xs,xr)。

更新規(guī)則數(shù)據(jù)集之前需要對(duì)重建數(shù)據(jù)的余弦加權(quán)操作，即每一道都需要乘以一個(gè)標(biāo)量系數(shù)，對(duì)于震源位于xs，檢波器位于xr的地震道d(xs,ys,xr,yr,t) ，它對(duì)應(yīng)的加權(quán)系數(shù)C計(jì)算公式為

(10)

(11)

同時(shí)，對(duì)式(10)進(jìn)行更新，每一次更新之后都需要將各個(gè)地震道的加權(quán)系數(shù)累加到總權(quán)重?cái)?shù)組Csum中，以便后續(xù)消除加權(quán)幅值不均勻?qū)е碌幕?。更新公式?/p>

(12)

步驟6重復(fù)步驟3～步驟5，直到完成遍歷整個(gè)數(shù)據(jù)集。通過滑動(dòng)局部窗遍歷數(shù)據(jù)的過程中，保持?jǐn)?shù)據(jù)窗尺寸不變，各方向的滑動(dòng)步長(zhǎng)取相應(yīng)方向的窗口尺寸的一半。

步驟7通過除以各個(gè)地震道的累加總權(quán)重Csum(xs,xr)，將計(jì)算得到的規(guī)則數(shù)據(jù)集Dreg權(quán)重歸一化。將重建的規(guī)則數(shù)據(jù)的各個(gè)地震道Dreg(xs,xr,t)除以其相應(yīng)的累加總幅值Csum(xs,xr)，以消除加權(quán)幅值不均勻?qū)е碌幕儭?/p>

3 模型數(shù)值模擬

選取東部某油田一個(gè)三維模型(圖2)，通過炮集數(shù)據(jù)的數(shù)值模擬及疊前深度偏移成像，驗(yàn)證文中隨機(jī)地震觀測(cè)系統(tǒng)設(shè)計(jì)及數(shù)據(jù)重建算法的有效性。

圖2 東部油田某工區(qū)三維地質(zhì)模型Fig.2 Three-dimensional geological model of a work area in the eastern oilfield

首先設(shè)計(jì)如圖3(a)所示的三維規(guī)則觀測(cè)系統(tǒng)，炮點(diǎn)橫向間隔50 m?？v向間隔400 m，檢波點(diǎn)橫向間隔400 m，縱向間隔50 m。考慮到目前實(shí)際生產(chǎn)中單方向隨機(jī)采集可行性最高，因此本次隨機(jī)觀測(cè)系統(tǒng)設(shè)計(jì)對(duì)炮點(diǎn)只在橫向上按照壓縮率50%進(jìn)行隨機(jī)采集設(shè)計(jì)[圖3(b)]，對(duì)檢波點(diǎn)只在縱向上按照壓縮率50%進(jìn)行隨機(jī)采集設(shè)計(jì)[圖3(c)]。

圖3 用于數(shù)值模擬的觀測(cè)系統(tǒng)Fig.3 Geometry used in numerical simulation

利用三維波動(dòng)方程正演得到模擬的規(guī)則炮集和隨機(jī)炮集數(shù)據(jù)[圖4(a)]，分別加入高斯噪聲與椒鹽噪聲[圖5(a)和圖6(a)]，利用本文提出的高維空間低秩約束的隨機(jī)地震數(shù)據(jù)重建算法得到重建后的規(guī)則地震數(shù)據(jù)，對(duì)比重建地震數(shù)據(jù)的波數(shù)譜。在波數(shù)譜中可以看出，模型正演地震記錄的波數(shù)譜展布區(qū)域較大[圖5(a)]，表明了地震信號(hào)的復(fù)雜度較高。此外，無噪聲及加入不同噪聲情況下，隨機(jī)地震數(shù)據(jù)重建結(jié)果較好且波數(shù)譜也得到了較好的恢復(fù)，并未在波數(shù)譜中引入明顯的空間假頻(圖4～圖6)，充分展示了高維空間低秩約束的隨機(jī)地震數(shù)據(jù)重建的性能。

將規(guī)則的地震數(shù)據(jù)按照地震道均勻抽取25%，得到稀疏的規(guī)則地震數(shù)據(jù)，基于完全相同的處理流程和速度模型進(jìn)行疊前深度偏移，對(duì)比在同樣采集密度下稀疏規(guī)則采樣和隨機(jī)采樣得到的數(shù)據(jù)成像效果(圖7、圖8)。

圖4 無噪聲情況下的隨機(jī)地震數(shù)據(jù)重建Fig.4 Reconstruction of random seismic data without noise

圖5 加入高斯噪聲的隨機(jī)地震數(shù)據(jù)重建Fig.5 Reconstruction of random seismic data with gause noise

圖6 加入椒鹽噪聲的隨機(jī)地震數(shù)據(jù)重建效果Fig.6 Reconstruction of random seismic data with salt and pepper noise

稀疏規(guī)則觀測(cè)系統(tǒng)造成了成像質(zhì)量明顯降低[圖 7(a)]，而隨機(jī)采集重建數(shù)據(jù)的成像結(jié)果與理想數(shù)據(jù)的成像效果基本相當(dāng)，但成像噪聲稍強(qiáng)[圖7(b)]。

為對(duì)比目標(biāo)層成像結(jié)果(大致位于X=2.4～4.5 km,Y=1.2～2.2 km的區(qū)域)，進(jìn)一步放大目標(biāo)層成像結(jié)果，如圖8所示，可以看出，相對(duì)于稀疏規(guī)則數(shù)據(jù)成像結(jié)果，隨機(jī)采集數(shù)據(jù)經(jīng)過重建后成像結(jié)果的信噪比更高，且對(duì)薄層的分辨能力更好。

圖7 不同數(shù)據(jù)的成像對(duì)比Fig.7 Imaging comparison of different data

圖8 不同觀測(cè)系統(tǒng)數(shù)據(jù)的局部成像對(duì)比Fig.8 Zoomed in imaging comparison of different data

4 結(jié)論

(1)針對(duì)具體的地質(zhì)任務(wù)，先進(jìn)行常規(guī)的規(guī)則地震觀測(cè)系統(tǒng)設(shè)計(jì)，在此基礎(chǔ)上結(jié)合工區(qū)地表及地下復(fù)雜程度，對(duì)規(guī)則的采樣點(diǎn)數(shù)進(jìn)行一定比例的壓縮，基于Jittered采樣原理進(jìn)行隨機(jī)觀測(cè)系統(tǒng)設(shè)計(jì)，是目前較為可行的方法。

(2)在進(jìn)行隨機(jī)觀測(cè)系統(tǒng)設(shè)計(jì)時(shí)，炮點(diǎn)和檢波點(diǎn)的隨機(jī)方式可以靈活選擇，在復(fù)雜地區(qū)可對(duì)炮點(diǎn)和檢波點(diǎn)只進(jìn)行線內(nèi)隨機(jī)；對(duì)于簡(jiǎn)單的地區(qū)，可以進(jìn)行完全隨機(jī)。

(3)本文提出的基于低秩約束的高維地震數(shù)據(jù)重建算法，充分利用地震數(shù)據(jù)在局部具有的線性特征，在不同背景噪聲下均能夠完成數(shù)據(jù)的無假頻重建，可進(jìn)一步應(yīng)用于實(shí)際地震數(shù)據(jù)重建。

(4)模型數(shù)值模擬結(jié)果證明，在同樣采集密度下，經(jīng)過規(guī)則化數(shù)據(jù)重建后，隨機(jī)地震數(shù)據(jù)比稀疏的規(guī)則采樣地震數(shù)據(jù)具有更好的成像效果，為進(jìn)一步在實(shí)際地震勘探中開展隨機(jī)采集試驗(yàn)奠定了基礎(chǔ)。