基于一次波稀疏約束的多道預(yù)測反褶積方法
【技術(shù)領(lǐng)域】
[0001] 本發(fā)明屬于地震勘探技術(shù)中地震信號處理領(lǐng)域，具體涉及一種基于一次波稀疏約 束的多道預(yù)測反褶積方法。
【背景技術(shù)】
[0002] 在海洋地震勘探中，預(yù)測反褶積用來消除水層多次波。多道預(yù)測反褶積方法 比單道預(yù)測反裙積方法能更好地消除多次波（M.T. Taner, "Long period sea-floor multiples and their suppression,'，Geophysical Prospecting, vol. 28, no. I, pp. 30-48 ，F(xiàn)eb. 1980.)。多道預(yù)測反褶積方法采用2D預(yù)測濾波器將原始數(shù)據(jù)的多道進行結(jié)合來預(yù)測 多次波。為避免可能存在的一次波損傷，多道預(yù)測反褶積采用相同的2D預(yù)測濾波器來同時 預(yù)測多道中的多次波。因此，可以在相互重疊的時間和空間數(shù)據(jù)窗口內(nèi)利用多道預(yù)測反褶 積進行多次波壓制，并能有效地適應(yīng)地震數(shù)據(jù)的時間和空間變化特性。多道預(yù)測反褶積方 法的數(shù)學模型為：
[0004]或
表示估計一次波，
表示原始數(shù)據(jù)，X = [X … X K,LX K+l, L …X K+1,L... xK, L …XK,L]表不 2D 預(yù)測濾波器，
[0007] U表示數(shù)據(jù)
的褶積矩陣。
[0008] 在方程⑴中，q為時間采樣點個數(shù)，表示時間方向的預(yù)測步長；2K+1和T。為時間 采樣點個數(shù)，分別表示2D預(yù)測濾波器的時間長度和數(shù)據(jù)窗口的時間長度；2L+1和X。為道 數(shù)，分別表示2D預(yù)測濾波器的空間長度和數(shù)據(jù)窗口的空間長度。
[0009] 在多道預(yù)測反褶積中，2D預(yù)測濾波器的系數(shù)絕對值中的大值集中在沿同 相軸傾角方向的有限范圍內(nèi)（J. Claerbout. (2006, Feb 28). Image estimation by example:Geophysical soundings image construction. [Online]. Available:http:// www. reproducibility. org/RSF/book/gee/)，稱之為濾波器的有限支撐域。在濾波器的系 數(shù)空間中，有限支撐域內(nèi)的濾波器系數(shù)對多道預(yù)測反褶積的多次波預(yù)測起主要作用。通常， 與單道預(yù)測反褶積中的ID預(yù)測濾波器相比，2D預(yù)測濾波器能夠更好地消除多次波，但需要 更多的計算時間。
[0010] 傳統(tǒng)多道預(yù)測反褶積對估計一次波施加能量最小化約束。另外，為避免濾波器估 計的不穩(wěn)定性，濾波器系數(shù)也假設(shè)滿足能量最小化約束。相應(yīng)的優(yōu)化問題為
[0012] 其中，λ為正則化參數(shù)。方程（3)中的2D預(yù)測濾波器可以采用最小二乘算法進 行求解：
[0013] X=(UtIHai)1UtU0 (4)
[0014] 最小二乘算法需要一次波和多次波的正交性假設(shè)。當一次波和多次波相互重疊 或有強一次波同相軸存在時，最小二乘算法會產(chǎn)生殘余多次波或造成一次波的損傷。在 單道預(yù)測反褶積中，對一次波施加稀疏約束，并采用迭代重加權(quán)最小二乘算法或快速迭代 收縮閾值算法求解ID預(yù)測濾波器，能有效地均衡一次波的保護和多次波的壓制。其中， 快速迭代收縮閾值算法采用ID收縮閾值算子促進一次波的非高斯性，相比于迭代重加 權(quán)最小二乘算法，在保持計算精度的同時，能進一步提高計算效率（L.Liu and W. Lu, "Α fast Lllinear estimator and its application on predictive deconvolution,'，IEEE Geosciences and Remote Sensing Letters, vol. 12, no. 5, pp. 1056-1060, May. 2015. )〇 然 而，傳統(tǒng)的多道預(yù)測反褶積方法需要求解濾波器系數(shù)空間的所有濾波器系數(shù)，并對一次波 施加能量最小化約束來求解2D預(yù)測濾波器，不能有效地均衡一次波的保護和多次波的壓 制。

【發(fā)明內(nèi)容】

[0015] 本發(fā)明的目的在于提出一種基于一次波稀疏約束的多道預(yù)測反褶積方法，該方法 首先確定2D預(yù)測濾波器的有限支撐域，降低所求解的濾波器系數(shù)的個數(shù)，然后構(gòu)建對一次 波施加稀疏約束的優(yōu)化問題，并采用快速迭代收縮閾值算法降低優(yōu)化問題求解的計算復(fù)雜 度，同時能有效地均衡一次波的保護和多次波的壓制。
[0016] 為了實現(xiàn)上述目的，本發(fā)明采用如下技術(shù)方案：
[0017] 基于一次波稀疏約束的多道預(yù)測反褶積方法，包括如下步驟：
[0018] a設(shè)置變量初始值，需要設(shè)置初始值的變量包括時間方向的預(yù)測步長q，2D數(shù)據(jù)窗 口的時間長度T。，2D數(shù)據(jù)窗口的空間長度X。，2D預(yù)測濾波器的時間長度? :，2D預(yù)測濾波器 的空間長度F，濾波器閾值P，一次波閾值sa，阻尼因子β和迭代次數(shù)#;
[0019] b輸入一個待處理的2D數(shù)據(jù)窗口中的數(shù)據(jù)u，根據(jù)比u的旅行時小q個時間采樣 點的2D數(shù)據(jù)窗口內(nèi)的數(shù)據(jù)以及數(shù)據(jù)窗口的長度參數(shù)、2D濾波器的長度參數(shù)構(gòu)造褶積矩陣 U，并確定濾波器的有限支撐域為：
[0021] 其中，受表示2D預(yù)測濾波器的粗略估計，兔,表示向量i中下標為（i，j)的元素 ，Ω 表示濾波器的有限支撐域；
[0022] c根據(jù)2D預(yù)測濾波器的有限支撐域，確定多道預(yù)測反褶積的數(shù)學模型為：
[0023] V = U-U0X0 ?
[0024] 其中，v表示一次波，^只包含有限支撐域內(nèi)的濾波器系數(shù)，U。為相應(yīng)的褶積矩 陣；計算逆矩P
I表示單位矩陣；
[0025] d利用步驟c得到的逆矩陣％對地震道集的2D數(shù)據(jù)窗口的數(shù)據(jù)u進行處理；
[0026] e判斷該地震道集中所有數(shù)據(jù)窗口內(nèi)的數(shù)據(jù)u是否全部處理完畢；如果否，返回步 驟b ;如果全部處理完畢，則首先采用2D漢寧窗將每一個2D數(shù)據(jù)窗口中估計的一次波進行 加權(quán)，并融合為一個道集％，然后采用同樣的方式將2D漢寧窗融合為一個道集％，最終的 一次波估計結(jié)果表示為
，其中，/表示逐個元素的相除操作。
[0027] 在步驟d中，利用逆矩陣Gq對地震道集的2D數(shù)據(jù)窗口的數(shù)據(jù)u進行處理的具體 過程為：
[0028] dl設(shè)置迭代數(shù)m = 0,利用步驟c得到的逆矩陣?Ω求取一次波的初始估計值：
[0030] d2令m = m+1,對上一步估計的一次波進行收縮閾值操作：
[0032] 其中，Ta為收縮閾值算子，定義為：
[0033] T a (v) j = (I v j I -s a C) +sgn (Vij j),
[0034] i = 1，2，...，T0 j = 1，2，...，X0
[0035] 其中，表示估計一次波的收縮閾值結(jié)果，v(m)表示估計的一次波，v u表示向量 V中下標為（i，j)的元素，
[0036] d3對收縮閾值結(jié)果進行更新：
[0038] 其中，表示更新后的收縮閾值結(jié)果，序列表達為t(1)= 1，
[0039] d4求取更新后的一次波估計結(jié)果：
[0041] d5判斷迭代次數(shù)m是否達到最大迭代次數(shù)M ;如沒達到，返回步驟d2 ;如果達到， 輸出當前數(shù)據(jù)窗口的一次波估計結(jié)果。
[0042] 本發(fā)明具有如下優(yōu)點：
[0043] 本發(fā)明方法首先確定多道預(yù)測反褶積中2D預(yù)測濾波器的有限支撐域和相應(yīng)的數(shù) 學模型，降低所求解的2D預(yù)測濾波器的系數(shù)個數(shù)，然后構(gòu)建對一次波施加稀疏約束的優(yōu)化 問題，并采用快速迭代收縮閾值算法求解2D預(yù)測濾波器，實現(xiàn)對多次波的壓制。相比于傳 統(tǒng)的多道預(yù)測反褶積方法需要估計濾波器系數(shù)空間的所有濾波器系數(shù)，并對一次波施加能 量最小化約束來求解2D預(yù)測濾波器而言，本發(fā)明方法能減小所求解的濾波器系數(shù)個數(shù)，有 效地均衡一次波的保護和多次波的壓制，同時降低求解優(yōu)化問題的計算復(fù)雜度。
【附圖說明】
[0044] 圖1為本發(fā)明中基于一次波稀疏約束的多道預(yù)測反褶積方法的流程示意圖；
[0045] 圖2a為偏移距等于40米的共偏移距道集圖；
[0046] 圖2b為將圖2a中每一道延遲400毫秒的共偏移距道集圖；
[0047] 圖2c為在圖2a中人為地加入一個傾斜的一次波同相軸的共偏移距道集圖；
[0048] 圖2d為真實的一次波的共偏移距道集圖；
[0049] 圖2e為濾波器有限支撐域的示意圖；
[0050] 圖3a為本發(fā)明中方法估計的一次波圖；
[0051 ] 圖3b為本發(fā)明中方法去除的多次波圖；
[0052] 圖3c為不利用濾波器有限支撐域的基于快速迭代收縮閾值算法的多道預(yù)測反褶 積方法估計的一次波圖；
[0053] 圖3d為不利用濾波器有限支撐域的基于快速迭代收縮閾值算法的多道預(yù)測反褶 積方法去除的多次波圖；
[0054] 圖4a為利用濾波器有限支撐域的基于迭代重加權(quán)最小二乘算法的多道預(yù)測反褶 積方法估計的一次波圖；
[0055] 圖4b為利用濾波器有限支撐域的基于迭代重加權(quán)最小二乘算法的多道預(yù)測反褶 積方法去除的多次波圖；
[0056] 圖4c為不利用濾波器有限支撐域的基于迭代重加權(quán)最小二乘算法的多道預(yù)測反 褶積方法估計的一次波圖；
[0057] 圖4d為不利用濾波器有限支撐域的基于迭代重加權(quán)最小二乘算法的多道預(yù)測反 褶積方法去除的多次波圖；
[0058] 圖5a為傳統(tǒng)的基于最小二乘算法的多道預(yù)測反褶積方法估計的一次波圖；
[0059] 圖5b為傳統(tǒng)的基于最小二乘算法的多道預(yù)測反褶積方法去除的多次波圖；
[0060] 圖5c為傳統(tǒng)的基于快速迭代收縮閾值算法的單道預(yù)測反褶積方法估計的一次波 圖；
[0061] 圖5d為傳統(tǒng)的基于快速迭代收縮閾值算法的單道預(yù)測反褶積方法去除的多次波 圖。
【具體實施方式】
[0062] 本發(fā)明的基本思想是：
[0063] 逐個2D數(shù)據(jù)窗口地進行多次波壓制，首先確定2D預(yù)測濾波器的有限支撐域，然后 構(gòu)建相應(yīng)的褶積矩陣和數(shù)學模型，并構(gòu)建對一次波施加稀疏約束的優(yōu)化問題：
[0065] 其中，u為