時空變自適應反Q濾波

郭志偉，曹思遠*，袁殿

1 中國石油大學(北京)油氣資源與探測國家重點實驗室，北京 102249

2 中國石油大學(北京)地球物理學院，北京 102249

3 北京億華通科技股份有限公司，北京 100084

1 引言

目前，全球油氣勘探呈現著明顯的全面(常規和非常規)、多樣(各種類型油氣資源)、更深(深層和超深層)的特點[1]。隨著地震勘探要求的日益提高，高分辨率的地震資料處理對后期的油藏描述和開發顯得至關重要。地震波在地層介質傳播過程中的吸收衰減是影響地震資料分辨率的主要原因，主要表現為振幅衰減、相位畸變、頻率降低等特征[2]。地層介質的吸收衰減作用是造成地震波能量衰減的一個重要因素，頻率越高，衰減越快，嚴重影響了地震資料的中深層成像能力和儲層預測的精度[3]。反Q濾波是對吸收衰減進行有效補償的方法，能較好地補償振幅和校正相位，具有較準確的物理機制[4]。在反Q濾波的實現算法方面，不同學者進行了大量的研究[5-14]。Robinson[15-16]提出了頻率域內插值的相位反Q濾波及其改進算法。Hale[17-18]根據Futterman模型提出了用泰勒級數展開近似高頻補償的反Q濾波。Bickel等[19]提出了基于復函數平面波的反Q濾波算法。McCarley[20]提出了基于自回歸模型的常Q濾波方法。Varela等[21]提出了Hale反Q濾波算法的改進方法。裴江云等[22]通過麥克勞林級數展式推導了基于Kjartansson模型的高階近似反Q濾波公式，其一階近似式與Hale的反Q濾波式一致。Bano[23]將常Q模型相位反Q濾波拓展到層Q值模型。在補償過程中，由于低信噪比和機器誤差等因素，高頻成分出現補償溢出的現象。針對該缺點，人們提出了增益補償函數，其作用是實現各頻率振幅合理有效的補償，通常在某個頻率點之后逐漸減小補償系數(即遵循低頻充分補償、中頻欠補償、高頻不補償的原則)。Wang[24-28]基于波場向下延拓提出了穩定高效的Gabor域反Q濾波，采用的增益函數只需一個控制參數(穩定因子法)，即可實現自適應的穩定補償處理。李合群等[29]在對比不同增益控制方法的基礎上，提出了選擇性更靈活的補償算子，取得了較好的補償效果。張固瀾等[30]針對固定增益函數容易造成深層資料高頻補償不足的缺點，提出了時變增益函數，較好地恢復地震有效頻帶內的能量。吳吉忠等[31]提出了一種基于等效Q值的反Q濾波算法，進行補償高頻能量耗散和相位畸變校正。SHI等[32]將增益限與穩定因子設為隨時間和地層Q值連續變化的參數，提出變增益限的反Q濾波方法，該方法能夠補償深層能量衰減，提高分辨率。趙巖[33]提出了一種變穩定因子的反Q濾波方法，其中穩定因子是地震波的傳播距離、地震記錄的頻率以及品質因子Q的函數。WEI等[34]引入整形正則矩陣來實現振幅補償的穩定計算，利用增益控制因子和噪聲壓制系數來壓制噪聲，改進了穩定化反Q濾波方法。XUE等[35]提出一種自適應穩定的反Q濾波方法，該方法利用同步壓縮變換對地震信號主頻段有效成分進行重構，同時壓制主頻段外環境噪聲，能有效地補償振幅衰減，抑制環境噪聲，提高地震資料分辨率。

由于反Q濾波處理過程中需要兼顧信噪比，一般情況下，地震資料的信噪比越高，則期望補償越充分。穩定因子法中的增益函數只有一個控制參量，其值越大，對高頻的補償越充分。基于上述特點和處理目標，通過建立增益參數和信噪比的映射關系，對資料不同信噪比的區域進行不同程度的補償處理，即自適應增益反Q濾波。

2 基本原理

2.1 反Q濾波

根據地層吸收衰減等價于偏移的思想，地表波場的向下延拓可以表示為：

其中，U(r,ω)是角頻率為ω傳播距離為r的平面波波場，k(ω)代表波數，i是虛數單位。地震剖面圖記錄了從震源到反射點然后再回到地表的反射波。在波數的定義中引入地層品質因子Q：

v(ω)為與頻率有關的相速度，將復波數k代入式(1)得：

用旅行時增量Δτ=Δ(ωr)代替距離增量，得到大地Q濾波的表達式：

ν(ωr)參考頻率ωr對應的相速度。根據Wang提出的改進的Kolsky頻散模型，相速度ν(ω)近似為：

ωh為調諧頻率，一般取值為地震帶寬的最高頻率，可以在反Q濾波時提供相對準確的相位校正。其中γ=(πQr)-1。注意：在頻散式(5)中只考慮正頻率，將之代入式(4)得：

式(6)為基于波場延拓的Q濾波的基本公式，根據該式得到反Q濾波式：

式中有兩個指數算子，前者為振幅補償算子，后者為相位補償算子。反Q濾波是在頻率域進行的，最后對所有平面波求和得到補償后的時間域地震信號：

2.2 自適應增益函數

直接對振幅進行補償時會放大隨機噪聲與機器噪聲，產生極大的不穩定性。為解決振幅補償值溢出和高頻噪音增大問題，采用穩定因子法構建增益函數。

衰減函數為β(τ,ω)：

增益函數為Λ(τ,ω)：

其中，σ2為穩定因子，Glim為增益參數(一般取10～100)，值越高，補償越充分。根據式(10)，當ω非常大時，Λ(τ,ω)趨向于1，表明增益函數對該高頻段不補償也不壓制。

我們接收到的地震記錄一般包括有效地震反射波和各種各樣的噪音，不同時刻的信噪比是不同的，即記錄的信噪比是時間和空間的函數。以信噪比函數作為自適應增益參數的選取準則，建立Glim與信噪SNR(τ)之間的映射關系F，則有

式中，穩定因子和增益參數與地震資料的信噪比有關，是時變參量，增益函數Λ~(τ,ω)具有自適應特征。當資料信噪比在時空域變化較大時，有必要采用式(11)來確定自適應穩定的增益參數進行反Q濾波。

映射關系F可根據信噪比SNR(τ)的實際特征設計，本文給出一種簡單的線性定義：

式中，待定參量Gmin、Gmax分別為最小和最大自適應增益參數，SNRmin、SNRmax分別為SNR(τ)的最小和最大值。

3 模型測試

3.1 無噪模型

設計5個反射系數界面的模型，反射系數值為1，分別位于 200 ms、350 ms、500 ms、650 ms和 800 ms，地面初始地震子波為50 Hz主頻的Ricker子波，地下介質Q值為100。圖1a是基于波場延拓正演的無噪衰減記錄(不含幾何擴散、反射/透射損失等固有衰減)，隨著傳播距離的增加，地震子波的延續長度逐漸增加、振幅衰減、相位改變。圖1b是利用式(10)反Q濾波的結果，增益參數Glim取50。可以看到，在無噪情況下，中深層地震子波的波長、振幅和相位都得到較好的恢復和校正。圖2是反Q濾波前(2a)、后(2b)時頻譜的對比，吸收衰減降低了中深層記錄的主頻，反Q濾波較好地恢復了中深層記錄的中高頻成分。

3.2 含噪模型

對無噪模型添加不同程度的隨機噪音，考慮在不同信噪比條件下，對比固定增益參數和自適應增益參數的反Q濾波效果。

圖3a是含噪記錄一，在無噪模型的基礎上添加了0.75%的隨機噪音，圖3b是固定增益參數反Q濾波結果，圖3c是自適應增益參數反Q濾波結果，其中，固定增益參數(見圖4a)取30，自適應增益參數(見圖4b)最大值取30，最小值取5。對比圖3a和圖3b，采用固定增益參數的反Q濾波較好地恢復了深層的弱反射信息(如800 ms處的同相軸)，同時擴大了隨機噪音，大幅降低了信噪比水平；對比圖3b和圖3c，采用自適應增益參數的反Q濾波在恢復弱反射信息的同時，有效地控制了噪音的放大，深層的信噪比得到較好的改善。圖5是反Q濾波前后的時頻譜對比，圖5a是原始含噪記錄一的時頻譜，圖5b是固定增益參數反Q濾波記錄的時頻譜，圖5c是自適應增益參數反Q濾波記錄的時頻譜。可以看到，固定增益參數的反Q濾波對噪音中高頻成分的放大效果較為明顯，自適應增益參數的反Q濾波則較好地控制了低信噪比區域的噪音(見紅色圈中部分)，提高了時頻譜的整體信噪比水平。

圖1 無噪模型反Q濾波剖面對比(a)反Q濾波前；(b)反Q濾波后Fig. 1 Comparison of inverse Q filtering pro file of the noise-free model. (a)before inverse Q filtering; (b)after inverse Q filtering

圖2 無噪模型反Q濾波時頻譜對比(a)反Q濾波前；(b)反Q濾波后Fig. 2 Spectrum comparison of the noise-free model with inverse Q filtering. (a)before inverse Q filtering; (b)after inverse Q filtering

圖6a是含噪記錄二，在無噪模型的基礎上添加了5%的隨機噪音，圖6b是固定增益參數反Q濾波結果，圖6c是自適應增益參數反Q濾波結果，其中，固定增益參數(見圖7a)取20，自適應增益參數(見圖7b)最大值取20，最小值取1。圖6a中的含噪記錄，隨著傳播距離的增加，信噪比越來越低，500 ms和650 ms處的有效反射波較弱，同相軸信息不明顯，800 ms處的反射信息淹沒在隨機噪音中，幾乎看不到同相軸。圖6b是采用固定增益參數的反Q濾波結果，與反Q濾波前剖面相比，信噪比進一步降低，500 ms和650 ms處的反射信息淹沒在隨機噪音中，分辨能力還不如反Q濾波前的剖面。圖6c是采用自適應增益參數的反Q濾波結果，500 ms和650 ms處的反射信息得到較好的增強，信噪比得到較好的控制，800 ms處的弱反射同相軸得到一定的恢復，分辨能力得到一定的改善。圖8是反Q濾波前后的時頻譜對比，圖8a是原始含噪記錄二的時頻譜，圖8b是固定增益參數反Q濾波記錄的時頻譜，圖8c是自適應增益參數反Q濾波記錄的時頻譜。在原始含噪記錄的時頻譜中(圖8a)，650 ms以下的隨機噪音與反射波的能量在同一量級；如果采用固定增益參數進行反Q濾波，對噪音的增強效果大于對有效反射波的恢復效果(圖8b中紅圈)，導致剖面中(圖6b)500 ms以下的同相軸淹沒在隨機噪音中，不可分辨；如果采用自適應增益參數進行反Q濾波，650 ms和800 ms處的有效反射振幅得到較好的恢復，有效反射的能量團在時頻譜中仍然占主要成分(圖8c紅圈)，清晰可辨。

圖3 含噪模型一(0.75%隨機噪音)反Q濾波剖面對比(a)原始含噪記錄；(b)固定增益反Q濾波；(c)自適應增益反Q濾波Fig. 3 Comparison of inverse Qfiltering profile in noise model 1 (0.75% random noise) (a)original noise recording; (b)fixed gain inverse Qfiltering; (c)adaptive gain inverse Qfiltering

圖4 模型一增益參數(a)固定增益參數；(b)自適應增益參數Fig. 4 Model 1 gain parameter (a)fixed gain parameter; (b)adaptive gain parameter

圖5 含噪模型一(0.75%隨機噪音)反Q濾波時頻譜對比(a)原始含噪記錄；(b)固定增益反Q濾波；(c)自適應增益反Q濾波Fig. 5 Spectrum comparison of inverse Qfiltering with noise model 1 (0.75% random noise) (a)original noise records;(b)fixed gain inverse Q filtering; (c)adaptive gain inverse Qfiltering

圖6 含噪模型二(5%隨機噪音)反Q濾波剖面對比(a)原始含噪記錄；(b)固定增益反Q濾波；(c)自適應增益反Q濾波Fig. 6 Comparison of inverse Q filtering profile of noise model 2 (5% random noise) (a)original noise records; (b)fixed gain inverse Qfiltering; (c)adaptive gain inverse Qfiltering

圖7 模型二增益參數((a)固定增益參數；(b)自適應增益參數)Fig. 7 Model 2 gain parameter ((a)fixed gain parameter; (b)adaptive gain parameter)

圖8 含噪模型二(5%隨機噪音)反Q濾波時頻譜對比(a)原始含噪記錄；(b)固定增益反Q濾波；(c)自適應增益反Q濾波Fig. 8 Spectrum comparison of inverse Qfiltering with noise model 2 (5% random noise) (a)original noise records; (b)fixed gain inverse Qfiltering; (c)adaptive gain inverse Qfiltering

4 資料處理

實際資料來自于沙漠工區，由于近地表吸收衰減嚴重，導致記錄的主頻較低(約為15～20 Hz)，嚴重影響了中深層的成像精度。圖9a是原始疊加剖面，目的層4 s附近，圖9b是采用固定增益參數的反Q濾波結果，圖9c是采用自適應增益參數的反Q濾波結果，圖9d是圖9b與圖9c所示的剖面之差。圖10為固定增益參數與自適應增益參數對比，其中，固定增益參數為40(圖10a)；自適應增益參數(圖10b)最小值為10，最大值為60，主要分布范圍在25～45。對比圖9a、圖9b和圖9c，反Q濾波后剖面的同相軸變細，垂向分辨率得到提高(見黑色箭頭)；對比圖9b和圖9c，兩個剖面的分辨率相近，但圖9c的信噪比高于圖9b(如黑圈部分)，圖9d顯示了兩者之差，主要是凌亂的反射和隨機噪音，如果把這部分信息加到圖9c中，則只會降低資料的信噪比，而對分辨率的提升卻有限，這表明采用自適應增益參數的反Q濾波能較好地控制低信噪比區域的噪音擴大，有效改善反Q濾波剖面的信噪比。

圖11是疊加記錄反Q濾波前后的時頻譜對比，圖11a是原始疊加記錄的時頻譜，圖11b是固定增益參數反Q濾波記錄的時頻譜，圖11c是自適應增益參數反Q濾波記錄的時頻譜。對比發現，反Q濾波前疊加記錄的時頻譜(圖11a)能量主要分布在10～30 Hz，經反Q濾波處理后，時頻譜(圖11b和圖11c)能量主要集中在15～40 Hz之間，中高頻成分得到較好的恢復；圖11c的高頻成分比圖11b弱，這是由于該段記錄的信噪比相對較低，采用了較小的增益參數(小于固定增益參數40)進行反Q濾波。值得注意的是，在圖11b中，4000 ms處的主頻(約為40 Hz)明顯高于4500 ms處的主頻(約為32 Hz)，在圖11c中，4000 ms和4500 ms處的主頻基本一致(約為30 Hz)，這說明自適應增益參數能較好地調節不同信噪比區域的補償效果。

圖9 疊后資料反Q濾波剖面對比(a)原始疊加記錄；(b)固定增益反Q濾波；(c)自適應增益反Q濾波；(d)剖面b與剖面c之差Fig. 9 Comparison of inverse Qfiltering profile of post-stack data (a)original overlay records; (b)fixed gain inverse Qfiltering;(c)adaptive gain inverse Qfiltering; (d)the difference between section b and section c

5 結論

穩定的反Q濾波處理為了兼顧分辨率和信噪比，需要不斷調試增益參數，將原始資料的信噪比作為先驗信息，用于約束增益參數的選取，可實現地震資料自適應的反Q濾波處理。采用自適應增益參數的反Q濾波，根據資料不同區域的信噪比水平，能自動調節補償效果，有效控制低信噪比區域的補償，避免噪音被過度放大。需要注意的是，為避免反Q濾波剖面的突變，自適應增益參數應進行適當的二維平滑處理。

圖10 實際資料增益參數。(a)固定增益；(b)自適應增益Fig. 10 Actual data gain parameter. (a)fixed gain; (b)adaptive gain

圖11 疊加記錄反Q濾波時頻譜對比。(a)原始疊加記錄；(b)固定增益反Q濾波；(c)自適應增益反Q濾波Fig. 11 Spectrum comparison inverse Qfiltering is recorded in superposition records. (a)original superposition records;(b)fixed gain inverse Q filtering; (c)adaptive gain inverse Qfiltering