疊前深度偏移地震記錄直接模擬方法

吳 涵 廉西猛 孫成禹* 芮擁軍 蔡瑞乾 鄧小凡

（①中國石油大學（華東）地球科學與技術學院地球物理系，山東青島266580；②中國石化勝利油田分公司物探研究院，山東東營257022）

0 引言

疊前深度偏移是目前油氣地震勘探資料處理中的關鍵技術。獲取疊前深度偏移地震記錄要經過炮記錄采集或正演、常規處理和偏移歸位處理等諸多環節。而地震正演與偏移成像存在計算量大、運行周期長等問題，且常規處理和偏移流程較為繁瑣。為提高正演和偏移的效率，研究人員做了許多工作：王德利等［1］對黏彈介質中三維波動方程進行了有限差分MPI并行正演及驗證；Micikevicius［2］針對GPU的有限差分實現給出了快速算法；龍桂華等［3］應用GPU實現了三維地震波交錯網格有限差分模擬的加速計算；李博等［4］、劉守偉等［5］、Foltinek等［6］研究了基于GPU加速的高階有限差分疊前逆時偏移，取得了顯著的成果；李振華等［7］利用無記憶擬牛頓—模擬退火法對偏移算子方程進行求解，提出了一種正則化偏移成像的全局優化快速算法；屠寧［8］對基于壓縮感知的快速最小二乘逆時偏移進行了研究，降低了最小二乘逆時偏移中波動方程的計算量。

盡管如此，獲取疊前深度偏移剖面流程繁瑣的根本問題依舊未能得到解決。為簡化這一復雜流程、快速獲得偏移剖面，人們往往應用一維褶積進行偏移剖面的模擬，但對于二維、三維復雜速度模型模擬效果不甚理想［9］。Hamran等［10-11］、Gelius等［12-14］、Wu等［15］、Zhu等［16］、朱小三等［17］在廣義散射層析成像的研究中詳細介紹了波場的成像條件及成像分辨率算子，實現了高精度的散射層析成像反演；馬在田［18］研究了地震偏移成像分辨率與觀測系統、成像點、震源等的關系，給出地震成像分辨率的定量計算公式，闡述了觀測地震道的分辨率和成像分辨率的空間變化的特征。

Lecomte等［19-22］分析了疊前深度偏移成像分辨率算子，提出用散射波數算子表征地震成像分辨率函數，并以此實現局部小區域的疊前深度直接模擬成像。在此基礎上，本文發展了一種基于深度域廣義卷積技術的疊前深度偏移地震記錄直接模擬方法，在已知速度模型和野外采集參數的前提下，無需經過波場正演和偏移處理，便可以直接快速生成偏移剖面，得到理論成像結果。根據理論成像結果，可以檢測觀測系統的優劣和地震數據的質量，為優化觀測系統參數和評價成像質量提供參考，以滿足數據采集和處理中質量控制的需求［23-25］。

1 方法原理

對地下介質正演和成像過程可以用算子的形式分別表示為

式中：F是模型的正演算子，用于生成模型S的正演數據U；B是成像算子，用于生成U的成像結果I；r表示模型中散射點位置；r′表示模型中散射點成像位置；rg為檢波點位置；rs為震源位置（圖1）。由式（2）可知，通過算子可實現由模型S直接獲取成像結果I。

圖1 地震散射波場傳播示意圖

標量聲波方程可表示為

式中：速度擾動函數O（r）一定程度上可近似看作反射率；C0（r）為背景速度場；k0＝ω/C0為背景波數；Psc、Pin分別代表散射波場和入射波場；v為散射體范圍。對于一個非均勻且無衰減的背景場，可以使用動態射線追蹤得到其格林函數G（rg，r）［26］，結合Stamnes［27］給出的格林函數近似解，并引入一個偏移算子l（rg，rs，r′），則此時散射點的成像結果可表示為

式中：D（r，r′）即是直接成像算子，本文將其稱為點擴散算子；a和τ為

其中γ表征振幅項。

結合Beylkin等［28-29］提出的理論，在成像點的鄰域內將振幅項和相位項Talor展開，結合程函方程并使點擴散算子趨近于δ函數，同時為了得到攜帶地震波響應的疊前深度偏移剖面，需將地震子波作用在點擴散算子上，則有

式中：s（ω）為地震子波的頻譜；，其中ks為沿著震源點到成像點的射線路徑上的波數向量，kg為沿著成像點到檢波器的射線路徑上的波數向量，且有。由此，可以使用地震子波信息、觀測系統參數和射線路徑信息在波數域構建點擴散算子，與地下速度模型的速度擾動或反射率信息相互作用，即可得到理論成像結果。同時，還可以將擾動點近似看成是一個小區域內的照明點，即擾動點處的點擴散算子可以對擾動點臨近區域進行局部模擬成像。

2 點擴散算子分析

由上文可知，將點擴散算子與速度擾動函數褶積即可得到模擬成像結果。速度擾動函數可以近似看作反射率，由地下速度模型直接得出。空間域的褶積運算相當于波數域的乘積運算，點擴散算子可以在波數域由震源信息、觀測系統信息及其射線路徑信息構建。為說明具體算法和效果，本文以主頻為25Hz的Ricker子波作為震源，起始炮點位置為50m，終止炮點位置為4600m，炮點距為50m，每炮共461道接收，道間距為10m，計算了Marmousi模型（圖2）中四個點處的點擴散算子（圖3），以此觀察點擴散算子的形態。由圖3可以看出，點擴散算子隨著深度的增加，地震波頻譜帶寬變窄、波長變長，表明成像分辨率逐漸降低。圖4、圖5展示了（1500m，2000m）附近局部區域（1000m×1000m）的模擬成像過程，可以看出，在給定速度模型信息、震源信息、觀測系統信息的情況下，構建的點擴散算子可以對地下構造進行準確模擬成像。

圖2 Marmousi模型及四個散射點示意圖

圖3 Marmousi模型不同位置處的波數（左）和空間（右）域點擴散算子

圖4 Marmousi模型（1500m，2000m）附近的速度分布（a）、反射率模型（b）及反射率模型的波數域振幅譜（c）

圖5 Marmousi模型（1500m，2000m）附近的直接模擬成像波數域結果（a）及空間域結果（b）

3 模型試算

通過簡單凹槽模型和Marmousi模型，驗證疊前深度偏移地震記錄直接模擬方法的準確性及其在指導觀測系統設計方面的高效性。

3.1 簡單凹槽模型

建立圖6所示的凹槽模型，尺寸為200×300個網格，網格間距Δx＝Δz＝5m。采用主頻為30Hz的Ricker子波作為震源進行試算，起始炮點位置為5m，終止炮點位置為1000m，炮點距為50m（共20炮），每炮共200道接收，道間距5m。疊前深度偏移地震記錄直接模擬結果（圖7）能夠準確反映地下層位位置、形態以及對應的波形信息。

3.2 Marmousi模型

為了檢驗算法的適用性和準確性，使用圖2所示的Marmousi速度模型進行試算。該模型尺寸為461×284個網格，網格間距Δx＝Δz＝10m。采用主頻為25Hz的Ricker子波作為震源，起始炮點位置為50m，終止炮點位置為4600m，炮點距為50m（共92炮），每炮共461道接收，道間距為10m。圖8為疊前深度偏移地震記錄直接模擬結果，能夠準確反映地下地質構造信息。與地下反射系數深度域褶積波形相似度很高（圖9），相關系數為0.9494，驗證了疊前深度偏移地震記錄直接模擬方法的準確性。

圖10為相同觀測參數下Marmousi模型不同方法的疊前深度偏移結果（正演數據由空間八階、時間二階的波動方程有限差分正演方法計算）。與圖8對比可以看出：逆時偏移成像結果在深部細節上，尤其是深部低速目標區以及中部層狀背斜構造處，成像效果較差（圖10a深部方框所示），且淺部存在明顯的噪聲干擾（圖10a淺部方框所示）；最小二乘逆時偏移能夠對深部細節進行準確成像，但深部與淺部反射率信息對應性較差（圖10b方框所示）；Kirchhoff偏移在中部層狀背斜構造和深部低速層附近的成像效果明顯較差，不能夠完整準確描述地下真實構造（圖10c方框所示）。由此可見，疊前深度偏移地震記錄直接模擬方法能夠為成像質量提供評價的參考標準。

圖6 凹槽速度模型

圖7 凹槽模型疊前深度偏移記錄直接模擬結果（a）及其波形顯示（b）

圖8 Marmousi模型疊前深度偏移地震記錄直接模擬結果

圖9 Marmousi模型直接模擬結果第200道歸一化波形與反射系數深度域褶積的對比

圖10 不同疊前深度偏移方法成像結果對比

圖11 Marmousi模型不同炮點距低速目標層模擬成像結果

圖12 Marmousi模型不同覆蓋次數模擬成像剖面第340道與對應反射系數深度域褶積的波形相關性分析

采用主頻為25Hz的Ricker子波作為震源，起始炮點位置為0，炮點距分別為100、50、30、20、10m，炮點右側接收，最小炮檢距為10m，最大炮檢距為1000m，道間距為10m，共100道接收，計算不同覆蓋次數（分別為5、10、17、25、50）的直接模擬成像結果。圖11展示了不同觀測系統參數下Marmousi模型深度2500m附近低速目標區域疊前深度偏移地震記錄直接模擬結果，覆蓋次數越高的觀測系統成像效果越好。炮點距為20m（覆蓋次數為25）時成像準確（圖11d），成像剖面第340道波形與對應反射系數深度域褶積波形相關系數為0.8190（圖12）；炮點距為30m（覆蓋次數為17）時成像結果（圖11c）也可識別出目標層及其附近構造，第340道波形與對應反射系數深度域褶積波形相關系數為0.7878（圖12）。兼顧經濟、效率兩方面要求，建議采用20～30m的炮點距進行地震采集。

4 結束語

本文研究了一種在已知速度模型和野外采集參數的前提下直接模擬疊前深度偏移記錄的方法。該方法可以簡單高效地獲得疊前深度偏移的理論成像結果，省去了繁瑣的炮記錄正演、常規處理及偏移等一系列過程。其模擬成像精度能夠準確反映地質構造信息及反射率信息，可為評價觀測系統和成像質量提供參考標準，以滿足采集觀測系統設計、成像處理和儲層分析等領域的研究需求。