沙漠區地震數據高分辨率處理關鍵方法及其在尼日爾Agedem 地區的應用

蘇 勤，曾華會，徐興榮，王德英，4，孟會杰

（1.電子科技大學資源與環境學院，成都 611731；2.中國石油勘探開發研究院西北分院，蘭州 730020；3.中國石油超深層復雜油氣藏勘探開發技術研發中心，新疆庫爾勒 841000；4.中國石油大學（華東）地球科學與技術學院，山東青島 257061）

0 引言

提高地震數據分辨率是一項系統工程，涉及采集、處理和解釋三大環節。近年來，雖然地震數據采集技術與工程設計方法快速發展，但由于受地震基礎理論和硬件發展限制，很難在采集階段進一步減少地震數據的能量損失。利用地震資料解釋中的分頻、交匯等技術，雖然能夠檢測出與巖性相關的衰減信號，但其嚴重依賴地震資料所保留的巖性相應信息。因此，地震資料的處理質量對提高數據分辨率和解決地質問題意義重大。

沙漠近地表地層多為疏松沙礫，其對地震波能量的吸收衰減嚴重，導致地震波在穿過表層后出現振幅衰減和相位畸變問題，因此沙漠工區資料多具有主頻低、頻帶窄的特點。如何有效且保真地提高資料主頻、拓寬資料頻帶、恢復衰減信號一直是地震資料處理人員亟待解決的難題［1-3］。反褶積作為諸多提高分辨率方法中最常用的方法之一，已經得到了快速發展，其中根據不同的地震子波估計方法，已經發展出了譜模擬反褶積和預測反褶積等多種方法。由于實際地震資料復雜多變，對其子波的估計技術存在瓶頸，因此反褶積方法對于消除沙漠表層吸收衰減影響的能力有限。針對這一問題，眾多學者提出了假設條件簡單、結果穩定的品質因子Q補償方法，如Wang［4］提出的穩定補償方法，在一定程度上解決了高頻補償過量的問題，但在面對雙復雜地區資料時，其假設條件難以保證補償結果的保真性。因此，為了提高假設的準確性，更多學者開始從地層品質因子Q入手，通過對估算Q和反Q濾波穩定補償算法的研究，從地質屬性的角度來解決補償假象的問題。從此，Q值的估算方法得到了快速發展，如王曉濤等［5］認為Q是決定能量衰減程度與速度的關鍵因素，因此從反射波的傳播衰減規律角度估算Q值；蔣立等［6］驗證了地表一致性反褶積方法與近地表Q補償方法的適用條件及穩定性，證明戈壁區適合采用地表一致性反褶積方法，沙漠區適合采用近地表Q補償方法；宋吉杰等［7］結合檢波點的實測Q值和折射波法估算的相對Q值構建了高精度Q場，提出了基于Q值融合估算方法，進一步促進了Q值的準確估算；李偉娜等［8］提出了一種基于雙線性回歸的Q估算法，通過模型與實際資料的測試和應用，獲得了與近地表模型分層一致且地質意義明確的估計Q場。以上Q值估算方法雖在一定程度上提高了估計精度，但并沒有對Q值估算的真實性進行研究，針對這一問題，利用雙測井約束估算Q值，可為結果的保真性提供重要依據。

壓縮感知作為一種新型采樣理論，現已在數學科學、醫學成像和天文成像中取得重要應用成果。國內外研究人員對這一理論進行了大量的研究，并將其應用于地震勘探領域。韓立國等［9］利用壓縮感知理論，在頻率域將擴展頻率數據的低頻能量替換原始數據等范圍內的低頻能量，得到拓頻后的頻率域信息，再進行反傅里葉變換進而得到低頻補償數據；張瑩［10］將擴展頻率數據的低頻和高頻能量與原始數據的主頻能量相結合，形成了最終的全頻補償地震數據；Sun 等［11］利用信息約束來加強低頻補償對噪聲的可靠性，但該方法僅限于低頻。以上3 種主要方法在拓展低頻和高頻方面均取得了良好的效果，但壓縮感知理論作為數學算法，其拓頻方法為數據驅動，并沒有對保真性進行保護。針對這一問題，利用井約束建立基于壓縮感知的拓頻基函數，在保證數據驅動算法優勢的同時，也利用井數據保護了其準確性。

針對沙漠區提出一種基于Q補償、Q層析建模成像和壓縮感知的聯合高保真、高分辨率處理（雙高處理）方法，首先利用約束估計的相對Q場對疊前記錄進行穩定Q補償，再對補償后的數據進行Q層析建模成像來增加成像的縱向分辨率，最后利用壓縮感知拓頻方法重建反射系數，進行疊后補償處理，并將該方法在尼日爾Agedem 地區進行野外數據應用，以期為提高薄儲層識別能力及提高油氣勘探效率提供依據。

1 近地表Q 疊前吸收補償方法

1.1 近地表Q疊前吸收補償方法理論

針對地震波傳播的吸收衰減，近地表Q補償方法可消除疏松地表帶來的影響。在整個聯合高分辨率處理過程中，首先將Q補償方法應用于疊前記錄的高分辨率處理，其中如何估計準確的Q值是關鍵［12］。為使結果具有準確性，借助微測井測量來求取近地表的實測Q值，具體分為如下幾個步驟。

（1）通過雙井微測井測量數據求取近地表的實測Q值，其結果對整個工區Q值的準確求取起決定性作用［13］。利用雙井微測井數據的直達波來進行Q值的求取，能很好地反映近地表低降速帶對地震波能量和頻率的影響。通過分析微測井數據井底檢波器和地面檢波器的頻率屬性和振幅屬性，結合質心頻移法和峰值頻移法各自的優點，采用改進的峰值頻移法，分析井底檢波器峰值頻率fm到地面檢波點拾取主頻f的移動，進而反演調查點的表層Q值，峰值頻移法的計算公式為

式中：Q實測為近地表實測Q值；fm為井底檢波器峰值頻率，Hz；f為地面檢波點拾取主頻，Hz；t為表層旅行時，s。

（2）計算表層旅行時t。根據拾取的折射波初至信息，利用層析反演等近地表速度反演方法，計算出工區的近地表速度、厚度、表層旅行時等信息。

（3）計算相對振幅衰減系數R。表層速度小或旅行時大，則表現為較強的吸收衰減，這種吸收衰減和地表具有一致性，為了體現這種相對差異的比例關系，定義R為相對振幅衰減系數。首先統計原始地震數據（折射波或反射波）每一道的振幅值，再在給定的時窗內計算每個地震道的均方根振幅值，最后分別統計炮點、檢波點的振幅能量，計算出近地表因素對原始地震數據振幅的影響，求取期望補償平均能量水平，進而得到相對振幅衰減系數R［14］。

（4）結合相對振幅衰減系數和表層旅行時，利用譜比法計算Q值：

式中：R為相對振幅衰減系數；Fscale為比例調節因子（目的是使RFscale<1，從而計算出的Q>0）；A0(f)為震源激發的地震波振幅；A(f)為衰減之后的地震波振幅；Q相對為近地表相對Q值。

（5）實測Q值和相對Q值的標定。在同一位置點，將微測井計算的實測Q值作為約束條件，標定基于疊前地震數據求取的相對Q值，根據兩者的差異，建立兩者的比值關系M，即

根據標定M值的區間范圍，進行Q值擬合處理。對于標定后的Q值與微測井計算的Q值線性誤差太大的點，進行人工剔除，直至擬合后的Q值與微測井計算的Q值非常接近為止，最后得到全工區Q場。

得到工區估計Q場后，按如下方法對疊前地震記錄進行補償。

將一維波動方程進行傅里葉變換：

式中：U（z，ω）為波函數的譜，其中z為波的傳播距離，m；ω為角頻率，Hz；kz為空間復波數。

考慮地震波吸收效應得到穩定的表層吸收補償算法，可對疊前數據在頻率域進行表層Q空變補償：

式中：U(τ+Δτ，ω) 表示經過振幅和相位補償后的頻率域地震數據，其中Δτ 為表層旅行時增量；U（τ，ω）為未經補償的頻率域數據，其中τ為近地表旅行時，s；ωh為參考頻率，Hz，是與地震波頻帶的最高頻率有關的一個調整參數。

式（6）是由波動方程推導的補償算式的基礎，等式右側由振幅補償項和相位補償項構成完整的Q補償。Wang［4］提出了穩定的Q補償算法，對振幅項做了如下改進：

式中：Λ(ω) 為穩定的振幅補償量；β(ω) 為原始振幅量；ρ2為穩定因子。

式中：Glim是增益限制，dB，是一個可調節參數；ωm為中心頻率，Hz，與地震波頻帶的最高頻率有關。

上述算法的穩定性體現在可以通過增益限制控制補償的頻帶范圍，防止高頻噪音的過度補償［16］，同時還可調整相位，解決地層吸收造成的地震波能量損失和頻散問題。

1.2 近地表Q疊前吸收補償測試

為驗證近地表Q的疊前補償效果，實驗首先選取尼日爾Agedem 地區疊前數據進行測試。從圖1 的單炮記錄展示中可以看出，處理后的數據具有更高的信噪比、分辨率和一致性，且同相軸更光滑連續，細小結構刻畫得更為清晰。由于低降速帶的影響，炮集頻率較低，2～3 s 處目的區域有效軸無法呈現，而處理后的數據有效信息能量增強，有效信息呈現更為清晰明顯。從頻譜中可以看出近地表Q補償后數據頻帶更寬（頻帶寬度從35 Hz 拓寬至約70 Hz），中—高頻能量得到明顯提高（主頻從30 Hz 提升至約45 Hz）。

圖1 尼日爾Agedem地區近地表Q補償疊前效果對比Fig.1 Comparison of near-surface Q compensation pre-stack effect in Agedem area，Niger

2 基于數據驅動的Q 層析建模方法

在Q疊前深度偏移實際應用的過程中，首先需要建立準確的深度域速度模型，再通過Q層析或其他方法得到與深度域速度模型相對應的深度域Q模型，最后進行Q疊前深度偏移，因此Q場的準確求取是開展Q疊前偏移的關鍵。在常規旅行時層析反演過程中，需要建立觀測到的旅行時誤差與地震射線所穿過介質的慢度的關系方程，再根據網格層析反演的方法求取地下介質每個點的Q值，而Q層析反演需建立觀測到的等效Q場（地震波傳播過程中Q效應的累計測量）與地下某個度量間隔內的Q值的對應關系，等效Q場定義為

式中：t1為射線總走時，s；Qeff為等效Q場；Qi為沿射線路徑上覆地層第i層介質的Q值；n為射線傳播路徑上的介質總層數；ti是沿射線路徑第i層介質的垂直旅行時，s。

使用初始速度模型、初始Q場共成像點道集進行射線追蹤，對于每一條射線，沿著射線路徑以及利用初始Q模型的層間Q值，可以累計得到等效Q場的模擬值，從而得到

定義Q層析的目標函數，求解大型稀疏線性方程組，即

式中：L為與吸收旅行時相關的反演大型稀疏線性方程組；S為平滑因子；t*=，為衰減旅行時，s，表示衰減和速度沿著傳播路徑的積分；μ為阻尼系數。

通過層析反演，得到對于網格單元（i，j，k）中初始Q模型的更新量，通過多次更新得到地下三維Q場。

3 Q 疊前深度偏移方法

在常規偏移成像過程中，振幅補償項一般只補償球面擴散能量損失，而不考慮介質吸收效應對不同頻率的能量影響，另外反Q濾波多為一維算法，僅補償地震波在垂直方向上的振幅損失和相位，而不考慮地震波實際傳播路徑的差異。

傳統積分法深度偏移沒有補償黏彈性介質吸收引起的幅值衰減和頻散校正，導致偏移成像分辨率較低。Q疊前深度偏移方法考慮了地震波在傳播過程中的介質吸收影響，消除了地震波傳播過程中由于介質非彈性引起的地震波吸收衰減和頻散問題，得到了垂向分辨率更高的成像結果，同時對吸收造成的能量衰減也進行了一定程度的補償。

用品質因子Q來描述地下空間介質，一般假定品質因子Q與頻率無關，即常Q模型。復速度c()

ω與實速度vreal()ω表達式為

如果考慮地震波的速度頻散公式，即

式中：ωref為參考頻率，Hz；j為虛數單位。

當頻率趨向于ωref時，vreal接近于常數，所以可用主頻ω0來替代高階頻率ωref，而ω0對應的vreal剛好是通過速度估計方法得到的，式（15）可以變換為

將式（15）代入式（13）可得

式中：v為常規的偏移速度，m/s；Q表示與頻率無關的恒Q值。

對任一地震道d（t），假設常規處理流程中反褶積處理消除了震源子波的影響，應用反褶積成像條件即可得到單道數據的成像結果，即脈沖響應：

式中：x，y，z為成像點的三維坐標；D(ω)為任一地震道d(t)的傅里葉變換；rs為炮點到成像點的直線距離，m；rg為成像點到檢波點的直線距離，m；τs為炮點沿射線路徑到成像點的走時，s；τg為成像點沿射線路徑到檢波點的走時，s；Qs為炮點處與射線路徑相關的等效Q 值；Qg為檢波點處與射線路徑相關的等效Q值；為成像權系數，補償了地震波的球面擴散影響。

4 基于井約束的壓縮感知疊后拓頻方法

4.1 基于井約束的壓縮感知疊后拓頻方法理論

壓縮感知作為一種新型采樣理論，現已在地震資料處理中獲得廣泛應用［15］。該理論克服了傳統方法中地質參數預構的弊端，基于數據驅動的優點讓其在尋找剖面潛力、拓寬頻帶寬度上限方面發揮了重要作用，可這也對壓縮感知拓頻方法的保真性帶來了一定影響［16］。為消除該影響，提出基于井約束的壓縮感知疊后拓頻方法，無論在拓頻效果還是在結果保真度方面都獲得了較好的效果。

壓縮感知理論認為一個不滿足奈奎斯特采樣定理的信號若是稀疏的，則其可以通過某些數學算法恢復為全采樣數據［17］。利用如下線性正演模型說明：

式中：y為原始地震數據；M表示隨機采樣矩陣；x為欠采樣數據；e代表隨機噪聲。

根據壓縮感知理論可知，若x為稀疏，則其可以被重建為接近y的數據，但由于地震數據普遍復雜多變，無法滿足稀疏條件，因此需要借助數學變換對其進行稀疏表示［18］。假設數據x'是數據x在稀疏域F的稀疏表示，則式（18）可以改寫為

式中：FH為稀疏變換F的共軛轉置矩陣。

若x'中有k個非零元素，且k<

式中：是x'的重建結果；η為自定義恢復誤差。

通過求解式（20）中的范數‖x'‖1，找到可以滿足誤差σ的最稀疏解，然后通過逆稀疏變化得到恢復數據y'：

從地震記錄的褶積模型入手來進行壓縮感知拓頻的研究，假設w為地震子波，r為地下反射系數，g為含噪地震數據，n為隨機噪聲，含噪地震數據g可表示為

將式（22）轉換到傅里葉域可表示為

式中：G，W，R和N分別為g，w，r和n的稀疏域表示。

理想的褶積模型中地下反射系數以脈沖形式，即持續時間趨向于零、帶寬趨向于無窮，經過地層吸收、干擾等因素影響得到的一個隨機序列，其原本頻譜應為全頻寬。由于地震子波具有濾波作用（這里子波充當壓縮感知理論中的觀測矩陣），褶積得到的地震數據就會損失一定的高、低頻信息，其頻譜也就展現為非全頻帶形式。在壓縮感知理論框架下，首先需要構造滿足RIP（有限等距性質）條件的觀測矩陣（前人驗證隨機矩陣具有再利用信號稀疏性和普適性），再利用信號稀疏性，通過求解稀疏正則化的凸優化問題重建缺失數據［19］。由于反射系數在時間域中具有稀疏性，因此將時間域作為稀疏域滿足壓縮感知重建理論應用的前提條件?；趬嚎s感知拓頻方法也就可以理解為重建被地震子波過濾的反射系數［20］。

在壓縮感知理論框架的指導下，將式（23）的欠定方程問題轉化為稀疏規則算子約束下的最小誤差能量問題：

式中：‖· ‖2為L2范數約束的誤差函數；‖· ‖1為L1范數約束的準確解；σ為反射系數的保真性誤差估計。

式（24）中第1 項為最小誤差能量約束，在求解過程中控制重建誤差，確保拓頻保真性；第2 項為稀疏規則算子約束以不斷獲取稀疏性更好的反射系數r。

為保證保真性，將井數據引入式（24）來約束反演過程，提高反射系數的準確性：

式中：rw為井位反射系數；r'w為根據井信息估計的反射系數；λ為正則化參數，用于調節兩范數之間的權重以確保稀疏性與保真性之間的平衡。

使用迭代閾值法求解式（25），迭代公式為

式中：(WF)T為WF的共軛轉置；Rw為提取相應井位信息的采樣算子。

通過如下迭代閾值法求解式（26）［21］：

初始化：k=1，賦初始值r0，外循環迭代次數Lin，內循環迭代次數Lout，預設誤差ε；

開始迭代：r0=r0；

外循環：判斷‖G-WFrk‖2>ε和k≤Lout；

迭代繼續：r0=rk-1；

計算拓頻后頻譜R'：R'=Fr'；

計算高頻補償后頻譜GH'：GH'=L(G)+H(R')；

求得高頻補償數據：gH'=F-1(GH')。

在上述算法中：Tλ為閾值函數；α代表算法參數，α必須大于(WF)TWF的最大特征值；λ用來調整反射系數的稀疏性；L和H分別為提取低頻和高頻算子。

4.2 基于井約束的壓縮感知疊后拓頻測試

為驗證井約束稀疏反演的保真優勢，同樣選用尼日爾Agedem 地區部分偏移剖面進行保真性對比試驗。如圖2 所示，相比于無井約束拓頻，有井約束結果與井合成記錄更為匹配，20～40 ms 層位位置更加準確，在180 ms 處，無井約束結果的層位假象問題被解決，總體補償效果具有更好的保真性。

圖2 尼日爾Agedem地區近地表Q補償疊前處理疊加效果對比Fig.2 Comparison of stacking effects of near-surface Q compensation pre-stack processing in Agedem area，Niger

5 野外實際測試

為了驗證上述組合高分辨率處理方法在沙漠區的應用效果，選用尼日爾Agedem 沙漠區采集的野外地震資料進行應用。Agadem 區塊位于尼日爾東南部的Termit 盆地，面積約27 516 km2，其典型地形條件是起伏的大沙丘，連片工區內地形南北高、中部低，全區范圍內沙丘呈南北向條帶狀分布，海拔為308～404 m，沙丘高度變化大，最高約70 m（圖3）。近地表低降速帶大致分為三層結構，全部是黃沙，低速帶厚度較大，速度約為300 m/s；降速帶厚度變化較大，主要受地表沙丘的影響，速度約為630 m/s。高大沙丘對地震波高頻信號吸收衰減作用明顯，提高分辨率是本次處理的一項重要工作。

圖3 尼日爾Agadem地區地表高程圖Fig.3 Surface elevation map in Agedem area，Niger

首先，通過分析地震波振幅、主頻和表層地震波傳播時間等屬性，發現表層模型和初至波反演得到的表層相對Q值與屬性相關性高，反映了地表對地震波吸收的空間變化（圖4），因此采用了近地表Q補償，以消除近地表吸收造成的吸收衰減，為精細解釋和儲層預測提供可靠資料，保證疊加成像。從圖5 中可以看出，處理后剖面淺—中層能量增強，深層有效信息增加，能量的補償有效克服了低降速帶的影響，使得層位清晰、斷層明顯、中深層成像準確。

圖4 尼日爾Agadem區塊地震波平面屬性Fig.4 Plane attributes of seismic wave in Agedem area，Niger

圖5 尼日爾Agadem區塊近地表Q補償前后疊加剖面及頻譜Fig.5 Stacked sections and spectrum before and after near-surface Q compensation in Agedem area，Niger

其次，針對底辟構造埋深大、反射能量弱、頻帶窄的問題，在此Q補償基礎上開展Q疊前深度偏移，消除地震波傳播過程中由于介質非彈性因素導致的能量衰減和頻散。按波場傳播路徑完成振幅補償和相位校正，消除地層吸收效應，得到縱向分辨率更高的成像結果，同時對吸收引起的能量衰減也進行一定程度的補償。從圖6 可以看出，處理后的剖面分辨率顯著提高，振幅和相位一致性也得到明顯改善，細節成像刻畫更為清晰，成像精度更高；強軸之間的目的儲層區域分辨率明顯提高，由于沙漠區低降速帶影響所導致的能量衰減嚴重問題得到了很好的解決。

圖6 尼日爾Agadem區塊Kirchhoff積分法疊前深度偏移（a）與Q疊前深度偏移（b）效果對比Fig.6 Comparison of Kirchhoff integration method（a）and Q-depth migration（b）for pre-stack depth migration in Agedem area，Niger

最后，在深度偏移成果上，依據區內已有測井資料，開展基于壓縮感知理論的拓頻處理，不僅進一步壓縮了子波，還加強了全區子波空間的一致性，同時，還提高了井-震分辨率吻合度。利用新資料，開展了多井精細標定，其中DB1 井標定后，全井段相關系數大于0.92（圖7），為儲層巖性描述和小斷層刻畫提供了更可靠的數據。

通過上述綜合高分辨率處理方法，進一步提高了地震資料的縱向分辨率。成果資料顯示，斷層結構清晰、斷面波干脆、小斷層錯段明顯容易識別，為準確落實和精細刻畫三維區塊內各級斷層提供了可靠資料，滿足了區域地質研究、區帶綜合評價、構造精細研究和儲量規模預測的地質需求。

6 結論

（1）近地表Q補償時采用雙井微測井測量數據求取近地表的實測Q值，恢復了地震波的高頻成分，高頻能量提升明顯，頻帶寬度有所改善，成像質量大幅提高。

（2）聯合應用近地表Q補償、Q層析及Q疊前深度偏移方法對尼日爾Agedem 地區地震資料進行處理后，地震數據的分辨率比常規中深層Q補償方法和地表一致性反褶積方法高，目的儲層區域分辨率明顯提高，由于沙漠區低降速帶影響所導致的能量衰減嚴重問題得到很好的解決；在該方法的基礎上結合區內已有測井，開展基于壓縮感知理論的拓頻處理，可進一步壓縮子波，全區子波空間一致性更好，提高了縱向分辨率，井-震吻合度更高。