段中鈺 李婷婷 肖 勇 王云雷 鄭桂娟
(①北京信息科技大學信息與通信工程學院,北京 100101;②東方地球物理公司國際勘探事業部,河北涿州 072751;③東方地球物理公司研究院,河北涿州,072751)
地震數據重建是地震信號分析領域的一個重要研究方向。地震數據采集時由于采集環境以及采集成本等因素[1],采集到的數據通常呈現不規則或稀疏分布,不完整的地震數據缺失部分地球物理信息,據其不能準確地預測地下巖層的性質。因此,為保證滿足后續地震資料處理對數據完整性的要求,必須恢復缺失的地震數據,以提高儲層預測的準確性,促進油氣的勘探與開發。顯然,重建缺失的地震數據具有非常重要的現實意義。
傳統的地震數據重建方法主要有基于濾波的重建方法[2-4]、基于波動方程的重建方法[5-7]以及基于變換函數的重建方法[8-11]等三種類型。濾波重建法是通過褶積插值濾波器實現地震數據重建,主要用于重建規則采樣地震數據,且計算量大。波動方程法是通過Kirchhoff積分算子實現地震數據重建[12],該方法需地下信息作為先驗條件,在地下信息未知或精度低時對重建結果影響很大。變換函數法是利用最小二乘法將信號變換到某一新的變換域,再通過反變換重構地震數據。常用的變換函數有傅里葉變換、曲波變換、離散余弦變換等。變換函數法結構簡單,計算復雜度低,可實現規則和非規則采樣地震數據的重建,且不需將地下地層信息作為先驗條件,但該方法對采樣率要求較高,造成采集成本高昂,且要求地震數據的采集規模較大。
壓縮感知理論[13-16](Compressive sensing,CS)的應用,克服了傳統地震數據重建方法的缺陷,解決了對采樣率要求較高的問題[17],使稀疏信號可在低于奈奎斯特采樣頻率的條件下用稀疏優化算法重建完整信號。壓縮感知地震數據重建主要分為稀疏表示階段和稀疏重建階段,重建階段相當于稀疏優化問題。白蘭淑等[18]結合稀疏促進的曲波恢復迭代算法和Bregman迭代閾值算法提出改進的聯合迭代算法,在壓縮感知理論框架下,用Curvelet變換稀疏表示地震數據,該改進算法能稀疏優化重建地震數據,具有收斂快、重建信噪比高的優點,但其計算量大、計算復雜度高。周亞同等[19]提出壓縮感知下K-SVD字典學習自適應稀疏表示地震數據、正則化正交匹配追蹤(Regularized orthogonal matching pursuit,ROMP)稀疏優化重建地震數據的方法,能加快復雜實際地震數據重建速度。但該稀疏優化算法得到的不是全局最優解,故重建精度較低。
交替乘子方向算法(Alternating direction method of multipliers,ADMM)[20-22]具有簡化復雜性、高重建速度和良好重建質量的優點,在地震數據領域引起了廣泛關注。Sternfels等[23]將ADMM算法與壓縮感知理論結合實現含噪地震數據的降噪。Zhang等[24]將ADMM算法應用于缺失地震數據的恢復,實現壓縮感知地震數據的精確重構。李慧等[25]在此基礎上提出壓縮感知框架下K-SVD字典學習稀疏表示地震數據,利用ADMM算法優化重建地震數據,達到提高地震數據重建精度的目的。
盡管ADMM算法在地震數據重建中效果顯著,但由于該算法中平衡因子無法根據地震數據自適應地進行調整,易損害重建精度。另外,在基于ADMM算法的壓縮感知地震數據重建模型中,由于輸入訓練樣本較少,算法參數較多以及模型復雜度較高等問題,容易出現過擬合現象。針對上述問題,本文對ADMM算法進行改進,提出平方正則交替乘子方向算法(Square regular-alternating direction method of multipliers,SR-ADMM)。一方面在ADMM算法迭代中加入平方正則項;另一方面自適應選取平衡因子。最終將SR-ADMM算法與壓縮感知理論相結合,實現基于壓縮感知的SR-ADMM算法地震數據重建。通過模擬地震數據、實際單炮數據以及大慶油田實際地震數據對本文所提算法進行了驗證。
壓縮感知理論具有三個重要的前提條件:地震數據的稀疏性、測量矩陣的不相干性以及合適的優化算法。壓縮感知的數據重建模型為
y=Φx
(1)
式中:x∈RN,為原始信號;Φ∈RM×N,為測量矩陣;y∈RM,為測量矩陣下原始信號的觀測值。該式表示原始信號從N維信號降維到M維獲得觀測值,信號降維過程可用圖1所示圖像表示。

圖1 原始信號降維
由于式(1)為欠定方程,有無窮多的解,很難直接通過觀測值y恢復出原始信號x。因此可用一個與測量矩陣不相關的變換矩陣對原始信號做稀疏表示,得到稀疏矩陣。其數學表達式為
x=DΘ
(2)
式中:D∈RN×N,為稀疏變換矩陣;Θ為稀疏系數,是長度為N的一維向量。式(2)的原始信號稀疏表示過程可用圖2所示圖像表示。

圖2 信號的稀疏表示
因此觀測值y可以表示為
y=Φx=ΦDΘ=AΘ
(3)
式中A=ΦD,為M×N階感知矩陣。A滿足約束等距性質,即可通過求解L0范數優化問題很好地從觀測值y中恢復出原始信號x
(4)
式中稀疏系數Θ的L0范數表示該向量中所有非零元素的個數。但L0范數的最優化問題是一個非確定性多項式問題,在稀疏變換矩陣D與測量矩陣Φ不相干時,求解Θ的L0范數問題等價于求解Θ的L1范數,故可將式(4)轉化為求解L1范數優化問題
(5)


(6)
適用的重構算法是壓縮感知精確重構的重要前提。ADMM算法的優點在于將對偶上升法的可分解性與乘子法的收斂性質結合起來,通過對偶上升法使線性稀疏項不斷接近最優解,保證在懲罰參數較小的情況下滿足精度要求;克服乘子法中變量二次項無法拆分的問題,將原始問題拆分求解,可實現并行運算且提高運算效率。使用ADMM算法求解壓縮感知地震數據重建問題的數學模型為
(7)
式中λ表示平衡因子,用于控制前、后兩項之間的權重。
引入輔助變量Z
(8)
該式的增廣拉格朗日函數為
(9)
根據上式可知,Θ的更新就是在固定參數Z和λ的情況下,求解增廣拉格朗日函數的最小值問題,即
L′Θ=AT(Aθk+1-y)+λk+ρ(θk+1-zk)
(10)
求解可知ADMM算法第k+1次的迭代表達式為
(11)
式中:u=λ/ρ,其中ρ表示懲罰參數;Sλ/ρ為軟閾值函數,定義為
(12)
以K-SVD字典學習對地震數據進行稀疏表示,并用ADMM算法解決稀疏優化問題,構成基于ADMM算法和K-SVD字典學習的壓縮感知地震數據重建模型。在該模型中,由于存在輸入訓練樣本較少、ADMM算法迭代參數較多及模型復雜度較高等問題,容易出現過擬合現象,因此增加平方鄰近項正則化式中的子問題,防止過擬合。正則化交替乘子方向算法公式如下
(13)
另外,由于ADMM算法中平衡因子的取值只能根據經驗選取,若平衡因子λ較大,則說明稀疏項權重更大,所求值更稀疏而不會過多考慮約束;反之,重建數據更接近于原始信號。為了使λ更精確,從而獲得更高重建精度,本文將λ的選取改進為可自適應方式。對平衡因子的改進方案如下
λ1,a=λ×αl,λ2,a=λ×αl-1,…
(14)
λ1,b=λ/αl,λ2,b=λ/αl-1,…
(15)
首先根據經驗設定初始平衡因子λ。兩式中α的作用是選取λ的鄰近值,可選擇略大于1或略小于1的數,如α=0.95。式(14)和式(15)分別表示對初始平衡因子做乘法和除法,獲得多個近似值。接著分別用所獲近似值以式(13)進行迭代實現重建,選取使地震數據重建誤差最小的平衡因子為最優平衡因子,通過一輪迭代就可得到平衡因子最優值;之后,繼續完成后續迭代過程。兩式中常數l表征近似解的個數,若l取值越大,就能獲得數量更多的平衡因子近似值,理論上重建效果更好。但l過大,會造成時間空間負荷,因此需根據地震數據實情選取。
本文基于壓縮感知的SR-ADMM算法地震數據重建的流程如圖3所示。

圖3 重建流程圖
在壓縮感知框架下,以K-SVD字典學習對地震數據做稀疏表示,用SR-ADMM算法稀疏優化實現缺失地震數據重建。缺失數據重建具體步驟為:
(1)輸入訓練樣本,對重建模型進行訓練,得到初始字典D0、初始稀疏系數Θ0;
(2)輸入原始缺失地震數據x,用初始字典D0對該地震數據進行稀疏表示;

(4)以K-SVD字典學習算法更新稀疏變換字典D,用更新后的字典D對步驟(3)初步重建地震數據做稀疏表示;

將信噪比作為地震數據重建性能的評判指標
(16)

驗證模型的模擬地震數據(圖4)共76道,522個采樣點,采樣間隔為1ms。圖5為該炮點集隨機采樣數據。ADMM算法中平衡因子根據多次實驗后得出,選取λ=0.15;SR-ADMM算法平衡因子則根據模擬地震數據自適應地選取。

圖4 原始模擬地震數據
OMP算法是壓縮感知框架下常用于解決優化問題的一種稀疏優化算法。分別采用OMP(圖6a)、ADMM(圖6b)和SR-ADMM(圖6c)三種算法對圖5所示隨機采樣模擬地震數據進行重建。

圖5 隨機采樣模擬地震數據
上述三種算法的重建信噪比和均方誤差如表1所示。將這三種算法重建結果(圖6,深藍色波形)與原始地震數據(黑色波形)進行對比,并將缺失較嚴重的第19~第33道數據進行局部放大。可見對于該隨機缺失模擬地震數據,基于壓縮感知理論的三種算法都能實現數據重建,OMP算法(圖6a)重建效果一般,從隨機缺失的第21道~第26道數據,易見缺失部分重建振幅與原始振幅存在誤差,且該算法對于未缺失部分(如第45道)數據產生影響,重建的均方誤差相對其他兩種方法較大,重建信噪比相對其他兩種方法最低;ADMM算法(圖6b)重建效果較好,不會對未缺失部分數據產生影響,但缺失數據的振幅尚未得到完全恢復;SR-ADMM算法(圖6c)重建的效果最好,對數據細節的恢復更完備,同相軸清晰可見,重建數據振幅與原始地震數據振幅一致,且該算法重建信噪比最高,重建均方誤差最小。

表1 模擬數據重建信噪比和均方誤差

圖6 三種算法重建結果與原始數據對比(a)OMP;(b)ADMM;(c)SR-ADMM
由于SR-ADMM算法中平衡因子是根據地震數據自適應地選取,所選平衡因子使數據重建誤差最小,從而達到提高重建精度的目的。而且,基于該算法的壓縮感知數據重建模型防止了數據的過擬合現象,一定程度上也提高了數據重建精度。因此,針對隨機缺失模擬地震數據,采用SR-ADMM算法重建的信噪比高于OMP和ADMM算法,即SR-ADMM算法的重建精度高于經典ADMM算法。
選取的M油田實際單炮地震數據(圖7)共有170道,采樣點數為1300,時間采樣率為1ms;其缺失率為40%的單炮地震數據如圖8a所示。分別使用OMP、ADMM及SR-ADMM三種算法對缺失數據(圖8a)進行重建(圖8b~圖8d)。ADMM算法中平衡因子據多次試驗后設定為λ=0.17,SR-ADMM算法的平衡因子據實際地震數據自適應地選取。三種算法對應的重建信噪比和均方誤差如表2所示。

圖7 實際地震數據單炮記錄

表2 實際數據重建信噪比和均方誤差
結合表2和圖8b~圖8d可看出:OMP算法(圖8b)重建效果欠理想,同相軸連續性較差,重建信噪比較低,均方誤差較大;ADMM算法(圖8c)重建效果比OMP算法稍好,大致實現了缺失數據的重建,但重建后的數據附近出現抖動,重建信噪比較高,均方誤差中等;SR-ADMM算法(圖8d)重建效果很好,同相軸更光滑,重建振幅連續性良好,在三種算法中重建信噪比最高,均方誤差最小。

圖8 實際缺失地震數據及三種算法重建后單炮記錄(a)實際缺失數據;(b)、(c)、(d)對應OMP、ADMM、SR-ADMM算法
為了更直觀地對比三種優化算法的重建效果,將隨機采樣原始數據中缺失較多的第59~第109道進行局部放大(圖9a),其中第59~第67道、第72~第74道、第76道、第78~第87道、第90~第95道、第100道及第104道為缺失數據。分別采用OMP(圖9b)、ADMM(圖9c)和SR-ADMM(圖9d)三種算法得到相應重建結果的局部放大圖。

圖9 原始數據及OMP、ADMM、SR-ADMM三種算法重建結果局部放大(a)原始地震數據局部放大;(b)、(c)、(d)對應OMP、ADMM、SR-ADMM三種算法重建結果的局部放大
對比圖9a與圖9b可看出部分缺失數據實現了重建(圖9b藍色箭頭所指),部分缺失數據未得到恢復(紅色箭頭所指),且重建數據振幅與原始數據有差距;對比圖9a與圖9c可見,缺失數據雖得到了恢復(圖9c藍色箭頭),但重建后的振幅與原始數據存在誤差;對比圖9a與圖9d發現,所有缺失數據都得到有效恢復(圖9d藍色箭頭),且重建數據振幅與原始數據一致。
因此,采用SR-ADMM算法可以很好地重建實際地震數據復雜的細節,重建性能較好。
上述模擬和實際地震單炮記錄驗證了本文方法的有效性,進一步將本文方法應用于中國石油大慶油田實際三維數據。從該三維地震數據的一條測線中選取200道(圖10),時間采樣率為1ms,采樣點為500個。圖11a為隨機缺失30%的該三維數據實例。ADMM算法中平衡因子據多次實驗后設定為λ=0.20,SR-ADMM算法的平衡因子基于實際地震數據自適應地選取。

圖10 三維地震數據實例
同樣分別采用OMP、ADMM及SR-ADMM三種算法對隨機采樣實際三維地震數據進行重建,得到相應的重建結果(圖11b~圖11d)。可見OMP算法的數據重建(圖11b)效果一般,未能很好地實現缺失數據的恢復重建,同相軸欠連續;ADMM算法重建數據(圖11c)效果相對較好,同相軸連續性提高,但存在同相軸抖動現象;SR-ADMM算法很好地實現了三維數據重建(圖11d),缺失位置的數據與相鄰位置同相軸的走向相符合,且很連續。

圖11 實際缺失三維地震數據及三種算法重建結果(a)實際缺失數據;(b)、(c)、(d)對應OMP、ADMM、SR-ADMM算法
為便于比較三種算法重建效果,將三種重建結果的第50~第60道與原始數據進行放大對比(圖12)。三種算法的重建信噪比及均方誤差如表3所示。
綜合分析表3和圖12所示局部放大圖,可以看出OMP算法恢復地震數據(圖12a)的振幅與原始地震數據(黑色)相差較大,且重建信噪比較低,均方誤差較大;ADMM算法重建(圖12b)的振幅與原始地震數據更契合,其重建信噪比高于OMP算法,均方誤差較小;而SR-ADMM算法重建后(圖12c)振幅與原始地震數據幾乎一致,重建信噪比最高,均方誤差最小。顯然,SR-ADMM算法在實際三維地震數據重建中也有效,且可提高實際三維地震數據的重建精度。

表3 三維實際數據重建信噪比和均方誤差

圖12 三種算法重建結果(藍色)與原始數據(黑色)的波形對比(a)OMP;(b)ADMM;(c)SR-ADMM
本文提出基于壓縮感知的SR-ADMM算法地震數據重建方法。該算法在原始ADMM算法基礎上加入正則項正則化迭代中子問題,并且自適應地選取平衡因子,利用該算法解決最優化問題實現缺失地震數據的重建。仿真驗證表明利用本文算法能有效地實現缺失地震數據,且其重建精度高于常規壓縮感知地震數據重建優化算法。
需要指出的是,數據的稀疏性是壓縮感知數據重建的前提,對重建效果有很大影響,數據越稀疏,重建效果越好。因此,如何更好的對地震數據進行稀疏表示是下一步的研究方向。