一種基于稀疏窗S變換的分頻-重構波阻抗反演方法*

冉喜陽 周懷來 張益明 李雷豪 楊吉鑫

(1.成都理工大學地球物理學院 四川成都 610059； 2.成都理工大學油氣藏地質及開發工程國家重點實驗室 四川成都 610059；3.四川省地震局 四川成都 610041； 4.中海油研究總院有限責任公司 北京 100028)

在地震勘探中，不同頻率的地震反射波特征反映了不同的地層巖性、物性及含流體性質[1],故相比常規地震資料，分頻技術的應用能夠從地震數據中獲取更多的有效信息，有利于對目標體進行更全面的分析解釋，進而提高儲層預測精度,降低解釋多解性。地震分頻技術的基礎是時頻分析方法，目前常見的時頻分析方法包括短時傅里葉變換、Gabor變換、Wigner Ville時頻分布、小波變換、S變換等[2]。

為了從地震數據中獲得更多的有效信息，研究人員對分頻技術做了很多的研究工作[3-12]。前人利用各種時頻分析方法對測井數據和地震數據進行分頻處理，發現不同頻率尺度地層特征在此基礎上進行波阻抗反演獲得了較高分辨率的反演結果，但在應用過程中也存在一些問題。王振卿 等[5]對塔中西部臺內數據用廣義S變換進行譜分解，使用帶通濾波處理得到15～25 Hz優勢頻段，對其進行波阻抗反演，雖然提高了識別精度，應用效果良好，但是在濾波過程中也濾掉了大量的有效波，忽視了部分細節。

另外，在常規的時頻分析方法中，由于短時Fourier變換的窗函數固定，不能適應非平穩信號中不同頻率成分特性的分析；S變換通過改進時間窗函數的尺度因子，具備多尺度分析能力，但S變換中基本小波是固定的，使其在應用中受到了限制；高靜懷 等[12]在S變換的基礎上，增加4個可調屬性參數得到廣義S變換，有效地增大了其應用范圍，方法效果尤為依賴對窗參數的調節。

稀疏窗S變換在保證高分辨率的同時還可以降低時間成本[13]。本文在前人研究基礎上，提出一種基于稀疏窗S變換的分頻-重構波阻抗反演方法，將地震數據分割成N個不同頻帶的數據體，并利用分頻井曲線約束進行反演，以提高對薄層的識別精度。

1 基本原理

1.1 稀疏窗S變換原理

稀疏窗S變換基于S變換，通過對窗口參數進行稀疏優化，實現自適應地調整非平穩信號的頻率突變，具有與非優化算法相同的計算復雜度。傳統濾波方法減少了頻譜帶寬，雖然使道分析結果變得穩定，卻無法獲得多個薄層的高分辨圖像。稀疏窗S變換作為一個強大的時頻分析方法能實現高分辨率復雜道分析，可以解決上述沖突。S變換的頻率域表達式為

(1)

用稀疏概念來解決這個問題，就是通過定義窗函數來控制它們的“影響區域”，在“影響區域”內定位強振幅頻率分量，并在整個時頻域內混疊處理低振幅分量。為了以有效的2D方式解決這個問題，Hamid Sattari[13]將式(1)以矩陣形式重寫為如下自適應形式：

(2)

基于振幅譜構建As[l][l,m]以區分不同頻率分量，需要滿足以下2個前提條件：

一是單位分解化，也就是其所有列加起來等于一個列單位向量，其矩陣形式為

(3)

為了使時頻分解可逆，要將其進行歸一化處理。而單位分解在窗參數優化中不僅僅是強制性的前提條件，同時也將稀疏性引入了該方法。

二是窗函數對稱性，表達式為

a(s[l],m)=a(s[-l],m)

(4)

該條件保持了S變換分析窗正負頻率的對稱性。

對于s[l]的求取，由于需要根據相應的振幅譜作調整以獲得自適應窗，在其求取過程中應用了輸入信號頻譜平滑后的結果。Hamid Sattari[11]構建的窗參數優化公式如下：

(5)

(6)

(7)

(8)

最后將得到a(s[l],m)分別置于As[l][l,m]的m列上即可得到自適應窗函數矩陣，再根據式(2)即可得到稀疏自適應S變換的時頻譜。

1.2 基于稀疏窗S變換的分頻原理

基于稀疏窗S變換的時頻濾波器可表示為

(9)

基于稀疏窗S變換的分頻,實際上是將通過稀疏窗S變換獲得的時頻譜與時頻濾波器相乘的過程。利用原始信號的頻譜圖加以控制，通過分頻處理使有用的信息被保留，不必要的信息被去除,獲得地震道時頻譜并進行濾波后得到所需的頻率分量。然而，這些數據仍然只是時頻數據，還需進行反變換才能獲得最終時間域的分頻結果[14]。使用反變換為

(10)

式(10)中：hf(t)是分頻后的地震信號；FFT-1是傅里葉反變換。先將二維的時頻譜乘以頻率域窗函數，接著對時間軸積分，再進行傅里葉逆變換，最終獲得了時域的分頻結果。

1.3 基于稀疏窗S變換的分頻-重構波阻抗反演流程

基于稀疏窗S變換的分頻-重構波阻抗反演流程包含疊后地震資料頻譜分析及分頻處理、測井資料預處理、波阻抗反演、數據融合等主要環節，最終獲得重構波阻抗反演數據體。利用上述分頻方法，設計了基于稀疏窗S變換的分頻-重構反演技術流程，如圖1所示。

圖1 分頻-重構波阻抗反演流程圖Fig.1 Frequency-divided impedance inversion flow chart

首先，對原始疊后地震數據進行頻譜分析，再根據其頻譜圖選取適當的頻率范圍將原始地震數據分為3段，通過分頻處理得到低、中、高頻地震剖面。然后，利用分頻井曲線約束分別對3個分頻地震數據進行波阻抗反演，之后將反演得到的各數據體再次進行頻譜分析，取3個反演結果中相對應的頻率部分進行數據重構，最終得到重構波阻抗數據體。

2 模型試算

為了驗證稀疏窗S變換的優越性，本文設計了2個模擬實驗。圖2為利用短時傅里葉變換、S變換、稀疏窗S變換3種時頻分析方法分別對Chirp信號進行的頻譜分析圖，可以看到短時傅里葉變換和S變換由于窗參數的固定，對不同頻率成分的能量聚集性較弱，隨著頻率變化，對于非平穩信號的時頻刻畫能力有所降低；而稀疏窗函數由于對窗參數進行了稀疏的優化，對強振幅的頻率分量作了相應的窗函數調節，得到的時頻譜分辨率較高、能量聚集性較強，且能很好的適應該信號。

圖2 Chirp信號仿真模擬Fig.2 Chirp signal simulation

圖3為合成記錄仿真模擬圖，其中第1列為設計的反射系數，第2列為反射系數與45 Hz雷克子波褶積的合成記錄。由圖3可以看出，稀疏窗S變換在時間分辨率變小的情況下依然能較準確的分離信號，而短時傅里葉變換在時間間隔較小的情況下信號的時頻譜出現耦合，無法分辨兩個相鄰信號，對識別造成了影響。這表明，稀疏窗S變換確實能有效提高時頻譜的分辨率。

圖3 合成記錄仿真模擬Fig.3 Synthetic record simulation

3 方法驗證

不同頻率的測井曲線對儲層的敏感性不同，大致來說，低頻曲線反映宏觀趨勢，高頻曲線反映薄層細節。為了驗證該分頻方法的可行性，本文選取含氣B井中的P-wave曲線進行試處理，圖4為曲線時頻譜，圖5為曲線分解圖。

圖4精確地反映了P-wave曲線的時頻分布規律，據此圖可將井曲線分成低、中、高頻3段；從圖5中原始井曲線與3條曲線的對比可以看出，聲波曲線對地下巖性變化較為敏感；低頻曲線反映尺度較大，主要刻畫地層組層速度變化；高頻曲線反映尺度較小，重點突出薄層內物性變化。3條分頻曲線與巖性剖面吻合較好，說明該分頻方法具有可行性。

利用上述分頻井曲線約束井旁地震道做單道分頻反演，以驗證該方法的正確性(圖6)。圖6中，第一列曲線為原始P-wave曲線，曲線1為利用本文方法得到的分頻反演曲線，曲線2為常規反演曲線。從原始井曲線與2條反演曲線的對比可以看出，2條反演曲線的大體趨勢都能較好地吻合原始井曲線。由于分頻反演較為充分的利用了地震數據的高頻成分，故在細節反映上分頻反演得出的結果明顯優于常規反演方法得出的結果，說明該方法的正確性以及優越性。

圖4 B井P-wave稀疏窗S變換時頻譜Fig.4 P-wave Sattari time spectrum of Well B

圖5 B井P-wave曲線分解圖Fig.5 P-wave curve decomposition diagram of Well B

圖6 B井旁道分頻反演圖Fig.6 Frequency-divided inversion of seismic traces of Well B

4 應用實例

為了驗證分頻-重構波阻抗反演的實用性，選取中國東南部L工區過A—F井的連井地震剖面進行分頻-重構波阻抗反演測試。圖7為原始地震剖面及其頻譜分析圖，可以看出該地區目的儲層較發育，但地震資料分辨率較低制約了砂體接觸關系的確定，儲層連通性存在不確定性，局部砂體較薄，常規波阻抗反演方法受到了限制。根據分頻-重構波阻抗反演方法進行反演，圖8～10分別為分頻后得到的低、中、高頻地震剖面，圖11為傳統波阻抗反演剖面，圖12為分頻-重構波阻抗反演剖面。由圖11與圖12對比可見，2種反演方法均能識別出有效的氣層，但由于常規反演方法更多受主頻控制，得到的反演結果分辨率較低，難以區分氣層和氣-水過渡帶；而分頻-重構反演合理地利用了地震數據中相對低頻與相對高頻的成分，使得反演結果與井吻合度更高，有效地提高了分辨率。此外，區域中C井附近斷層較發育，但傳統波阻抗反演難以獲得有效相關信息，對斷層反映極其不敏感，而分頻-重構反演的剖面可以在目標井附近顯示出明顯的地層錯斷，說明其對斷層的識別更有效；并且傳統的波阻抗反演難以反映出地層的薄厚變化，而分頻重構波阻抗反演卻能更細膩準確地描述地層。以上3個優勢顯示了該反演方法有較好的應用潛力。

圖7 中國東南部L工區原始地震剖面及其頻譜分析圖Fig.7 Original seismic profile and its spectrum analysis in L work area,Southeast China

圖8 中國東南部L工區低頻地震剖面及其頻譜分析圖Fig.8 Low-frequency seismic profile and its spectrum analysis in L work area,Southeast China

圖9 中國東南部L工區中頻地震剖面及其頻譜分析圖Fig.9 Mid-frequency seismic profile and its spectrum analysis in L work area,Southeast China

圖10 中國東南部L工區高頻地震剖面及其頻譜分析圖Fig.10 High-frequency seismic profile and its spectrum analysis in L work area,Southeast China

圖11 中國東南部L工區常規波阻抗反演剖面Fig.11 Conventional impedance inversion profile in L work area,Southeast China

圖12 中國東南部L工區分頻-重構波阻抗反演剖面Fig.12 Frequency-divided impedance inversion profile in L work area,Southeast China

5 結論

1) 稀疏窗S變換對窗參數進行了基于振幅譜的稀疏優化，令窗參數根據不同信號的振幅譜作不同適應性窗參數優化，提升了運算速度，該時頻分析方法具有分辨率高和能量聚集性好的優勢。分頻-重構波阻抗反演可以充分利用地震資料中的低頻和高頻信息，獲得的最終剖面具有與井吻合度高、斷層反應清晰、地層刻畫準確等特點。

2) 本文方法的關鍵之處在于如何調節參數控制分頻過程和效果以及如何將低、中、高頻波阻抗數據體進行重構，這都將直接影響到最終反演結果。重構結果高頻成分過重，會產生一些假象；重構結果低頻成分過重，將影響反演的分辨率。如何更精確控制分頻反演中各頻率成分占比，使反演結果更合理，也是作者之后的研究內容。