基于壓縮感知的地震觀測系統設計及數據重建

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

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

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

壓縮感知理論是一種新的采樣理論，基于信號的稀疏特性，在遠小于Nyquist采樣率的條件下，通過隨機采樣獲取信號的離散樣點值，然后進行信號的重建[2-3]。壓縮感知理論一經提出，就引起學術界和工業界的廣泛關注，在各個領域得到迅速發展。

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

利用壓縮感知理論實現隨機地震采樣并重構規則數據，首先要分析地震數據的稀疏性，地震信號(或波場)的時空復雜度(與地下及地表速度結構、噪音特征等相關)會顯著影響地震信號的稀疏性，進而導致很多情形下隨機采集的實際數據的重建效果無法達到規則采樣數據的質量[13]。中國陸上油田地表噪聲發育，地下波場復雜，基于壓縮感知理論進行隨機地震數據采集及規則化重建具有更大難度。

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

1 隨機地震觀測系統設計基本理論

隨機地震采集需要滿足三個條件：地震數據具有稀疏性、合理的隨機采樣方法和有效的數據規則化重建算法。關于地震數據的稀疏性，很多學者的研究表明實際地震數據在時間域或其他域具有稀疏特征，即地震數據可稀疏表達[13-15]。現重點研究隨機地震觀測系統設計及數據規則化方法。

1.1 隨機地震采樣矩陣

野外地震數據觀測實際上是對信號的一種采樣，地震數據采集過程可以表示為

d=Df

(1)

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

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

d=Ak

(2)

A=DCH

(3)

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

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

1.2 隨機地震觀測系統設計方法

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

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

結合實際地震數據特點，基于Jittered采樣原理提出了以下隨機地震觀測系統設計技術流程。

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

步驟2根據工區地表及地下復雜程度設置縱橫向壓縮采樣比例{rx,ry}，基于該壓縮采樣比例計算隨機采樣點坐標的期望值{xe,ye}。

步驟3設置隨機采樣點坐標相對于期望值的擾動量，即高斯隨機分布的標準差(σx,σy)，一般取(σx,σy)=(3Δx,3Δy)。

最后，根據步驟2中期望值，以及步驟3中標準差利用隨機函數生成隨機觀測網格{xrand,yrand}。

2 高維空間低秩約束的隨機地震數據重建

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

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

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

步驟1定義一個理想的規則工區網格G，按照規則觀測系統設計方法[20]在網格G上定義規則觀測系統，記為Ωreg。

步驟2獲取原始輸入數據集Dirr的非規則地震數據的觀測系統。以炮集為單位獲取炮檢點坐標，取整到臨近網格，得到Dirr對應的非規則觀測系統Ωirr。

步驟3以窗口滑動方式從原始數據集取得一個局部高維數據子集。通過在炮端檢端分別取局部窗口，抽取原始道集中落在窗口中的所有地震道，從而獲得一個不完備的五維局部數據體。通過滑動中心炮檢點，可以循環獲取不同位置的局部高維數據。

步驟4在步驟3獲取的不完備五維局部數據體上，在頻率域進行規則數據重建，輸出一個完備且規則的局部五維局部數據體。計算過程即求解優化問題，即

(4)

該優化問題實現了對重建結果的三種約束。

(1)要求在觀測點位上的重建結果與原始觀測相吻合(L2項)。

(2)要求重建數據在四維波數具有稀疏性(L1項)。

(3)要求重建數據在進行Hankel排列后具有低秩性(核范數項)。

式(4)所示優化問題的求解，可通過以下迭代算法實現。

(1)關于L2項的梯度更新為

(5)

(6)

(7)

(2)關于L1項的軟閾值收縮為

(8)

(3)關于核范數項的奇異值收縮：

(9)

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

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

更新規則數據集之前需要對重建數據的余弦加權操作，即每一道都需要乘以一個標量系數，對于震源位于xs，檢波器位于xr的地震道d(xs,ys,xr,yr,t) ，它對應的加權系數C計算公式為

(10)

(11)

同時，對式(10)進行更新，每一次更新之后都需要將各個地震道的加權系數累加到總權重數組Csum中，以便后續消除加權幅值不均勻導致的畸變。更新公式為

(12)

步驟6重復步驟3～步驟5，直到完成遍歷整個數據集。通過滑動局部窗遍歷數據的過程中，保持數據窗尺寸不變，各方向的滑動步長取相應方向的窗口尺寸的一半。

步驟7通過除以各個地震道的累加總權重Csum(xs,xr)，將計算得到的規則數據集Dreg權重歸一化。將重建的規則數據的各個地震道Dreg(xs,xr,t)除以其相應的累加總幅值Csum(xs,xr)，以消除加權幅值不均勻導致的畸變。

3 模型數值模擬

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

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

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

圖3 用于數值模擬的觀測系統Fig.3 Geometry used in numerical simulation

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

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

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

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

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

稀疏規則觀測系統造成了成像質量明顯降低[圖 7(a)]，而隨機采集重建數據的成像結果與理想數據的成像效果基本相當，但成像噪聲稍強[圖7(b)]。

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

圖7 不同數據的成像對比Fig.7 Imaging comparison of different data

圖8 不同觀測系統數據的局部成像對比Fig.8 Zoomed in imaging comparison of different data

4 結論

(1)針對具體的地質任務，先進行常規的規則地震觀測系統設計，在此基礎上結合工區地表及地下復雜程度，對規則的采樣點數進行一定比例的壓縮，基于Jittered采樣原理進行隨機觀測系統設計，是目前較為可行的方法。

(2)在進行隨機觀測系統設計時，炮點和檢波點的隨機方式可以靈活選擇，在復雜地區可對炮點和檢波點只進行線內隨機；對于簡單的地區，可以進行完全隨機。

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

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