一種頻率域三維不規(guī)則地震數(shù)據(jù)重建方法
【技術(shù)領(lǐng)域】
[0001] 本發(fā)明涉及的是一種缺失道地震數(shù)據(jù)重建的方法，具體是一種頻率域三維不規(guī)則 地震數(shù)據(jù)重建方法。 技術(shù)背景
[0002] 在野外數(shù)據(jù)采集過程中，由于采集設(shè)備、項(xiàng)目成本以及地形條件的限制，地震數(shù)據(jù) 沿空間方向上通常低于奈奎斯特采樣定理要求進(jìn)行不規(guī)則采樣，從而導(dǎo)致采集到的地震數(shù) 據(jù)不規(guī)則、不完整，出現(xiàn)空間假頻，難以滿足后續(xù)處理的要求，例如速度分析、自由表面多次 波消除、偏移歸位等，在這種情況下，需要發(fā)展較好的地震數(shù)據(jù)重建方法，重建出缺失的地 震道，得到滿足處理要求的地震數(shù)據(jù)，并且該方法也可以指導(dǎo)野外數(shù)據(jù)采集，壓縮地震數(shù)據(jù) 采集量。
[0003]基于曲波變換的數(shù)據(jù)重建方法不需要地下結(jié)構(gòu)的先驗(yàn)信息，能夠重建出規(guī)則缺失 和不規(guī)則缺失的地震道，且計算速度快，精度高，得到了廣泛的應(yīng)用。薛念等采用閾值迭代 法，實(shí)現(xiàn)了基于曲波變換的二維地震數(shù)據(jù)重建(基于Curvelet變換的地震數(shù)據(jù)插值和去噪 [碩士論文]，西南交通大學(xué)，2010)，劉國昌等采用凸集投影算法，在曲波域?qū)崿F(xiàn)了基于P0CS 算法的二維地震數(shù)據(jù)重建(基于curvelet變換的缺失地震數(shù)據(jù)插值方法[J].石油地球物理 勘探，2011，46(2)，237-245。），該方法更能夠重建非線性同相軸，更加適合具有各項(xiàng)異性特 征的地震波場，但沒有深入研究基于曲波變換的三維數(shù)據(jù)重建方法，而且對閾值參數(shù)選擇 也沒有進(jìn)一步研究，僅僅采用指數(shù)閾值參數(shù)公式，其重建后的效果還是有限。而后張華等在 此基礎(chǔ)上實(shí)現(xiàn)了基于曲波變換的三維地震數(shù)據(jù)重建(基于jitter采樣和曲波變換的三維地 震數(shù)據(jù)重建[J].地球物理學(xué)報，2013,56(5): 1637-1649。），通過對時間切片進(jìn)行逐次重建， 取得較好的重建效果，并且探討了閾值參數(shù)選擇策略。王本鋒等提出了不同閾值下的曲波 重建方法(P0CS聯(lián)合改進(jìn)的Jitter采樣理論曲波域地震數(shù)據(jù)重建方法[J].地球物理學(xué)報， 2015,50(1):20-27。），然而該方法在時間域進(jìn)行重建，且只采用硬閾值參數(shù)進(jìn)行處理，當(dāng)重 建數(shù)據(jù)同時包含高和低能量振幅同相軸時，由于P0CS算法在前幾次迭代過程中，低于閾值 參數(shù)以下的同相軸系數(shù)被濾除，因此這些方法在重建低振幅有效波同相軸的效果有限，不 能完全恢復(fù)出地震信號的低能量同相軸及提高局部異常區(qū)的信噪比，并且迭代次數(shù)多，重 建效率低，制約了該方法的進(jìn)一步應(yīng)用。

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

[0004]本發(fā)明的目的是為了能夠快速高精度的重建出缺失道三維地震數(shù)據(jù)，而提供了一 種頻率域三維不規(guī)則地震數(shù)據(jù)重建方法。
[0005]本發(fā)明技術(shù)方案：
[0006]本發(fā)明提出頻率域三維不規(guī)則地震數(shù)據(jù)重建方法，首先利用傅立葉變換將時間域 三維地震數(shù)據(jù)變換到頻率域，然后采用凸集投影算法思想和引入能夠刻畫地震數(shù)據(jù)局部化 特征的曲波變換。在迭代過程中提出新的按指數(shù)規(guī)律衰減的閾值參數(shù)，采用軟閾值算子，對 每個頻率切片單獨(dú)進(jìn)行重建，降低了迭代次數(shù)和提高了重建精度，從而達(dá)到重建三維不規(guī) 則地震數(shù)據(jù)的目的。
[0007] 一種頻率域三維不規(guī)則地震數(shù)據(jù)重建方法，地震數(shù)據(jù)的重建問題描述為由一組不 完整數(shù)據(jù)通過線性算子的作用恢復(fù)出完整的數(shù)據(jù)，假設(shè)如下線性正演模型：
[0008] y〇bs=Mf
[0009] 這里yobseR"代表采集的地震數(shù)據(jù);feRN，且N2n，表示待重建的無假頻數(shù)據(jù);Me RnXN表示對角矩陣，其元素1和0分別表示已知地震道和未知地震道，假設(shè)系數(shù)X是f在曲波域 C中的稀疏表示，則上述方程為：
[0010]
[0011] 這里上標(biāo)Η代表共輒轉(zhuǎn)置矩陣，當(dāng)從采集的數(shù)據(jù)yQbs中重建無假頻數(shù)據(jù)f時，由于X 是稀疏的，從而可以采用稀疏算法求解該欠定方程。
[0012] 在稀疏促進(jìn)反演后，重建信號由；7 = 決定，同時
[0013]
[0014] 在這個表達(dá)式中，i代表估計值，Li范數(shù)定義為K ，x[i]是向量X中第i個 元素，通過求解上述方程，原始無假頻的數(shù)據(jù)就重建出來；
[0015]在地震數(shù)據(jù)重建過程中，采用凸集投影算法，具體重建步驟如下：
[0016]步驟1:在時間域中輸入具有缺失道的三維地震數(shù)據(jù)，然后采用傅立葉變換，將重 建數(shù)據(jù)從時間域變換到頻率域；
[0017]步驟2:采用曲波變換對頻率域三維地震數(shù)據(jù)進(jìn)行稀疏變換，得到曲波域中曲波系 數(shù)，根據(jù)曲波系數(shù)的大小選擇合適的閾值參數(shù)公式M(i=l，2,3，· · ·，Ν，其中N為迭代次 數(shù))及其公式；
[0018] 步驟3:在曲波域中，采用軟閾值算子進(jìn)行處理，也即大于閾值人,的曲波系數(shù)減去 一個閾值大小的值，小于閾值的曲波系數(shù)加上一個閾值大小的值，而其它的曲波系數(shù)置 零；
[0019] 步驟4:將閾值化后的曲波系數(shù)做反曲波變換得到時間域地震數(shù)據(jù)，然后再用缺失 地震數(shù)據(jù)的未缺失地震道填充到反變換后的地震數(shù)據(jù)中去；
[0020] 步驟5:最后將得到的地震數(shù)據(jù)代入步驟(2)，重新進(jìn)行迭代，直到運(yùn)行Ν次結(jié)束，再 對迭代Ν次后的地震數(shù)據(jù)做反傅立葉變換即得到最終的重建結(jié)果。
[0021] 進(jìn)一步，所述曲波變換的定義為：
[0022]
[0023] 式中：表示曲波函數(shù)，j，l，k分別表示尺度，方向和位置參數(shù)，f(x)為地震數(shù) 據(jù)，其頻率域定義式為：
[0024]
[0025] 經(jīng)過變換后得到的曲波系數(shù)，可用(：{」}{1}仏1氺2)表示其結(jié)構(gòu)，其中」表示尺度，1 表示方向，（lu，k2)表示j尺度1方向上的矩陣系數(shù)。
[0026] 進(jìn)一步，所述軟閾值算子，其表達(dá)式如下：
[0027]
[0028] Tk表示軟閾值算子，Sk-i表示第k-Ι次迭代得到的重建數(shù)據(jù)的曲波系數(shù)，滿足 L=G7(^Ft和F-\表示關(guān)于時間變量t的正反傅立葉變換，74._,表示第k_l次迭代得到的 時間域重建數(shù)據(jù)，其中表示N維閾值集合，· · ·，λΝ}，且滿足λρλ；^ > · · · >λΝ，Ν表示最大迭代次數(shù)。
[0029] 進(jìn)一步，閾值參數(shù)公式，其表達(dá)式如下：
[0030]
[0031]其中Max為|CFy|的最大值，即曲波變換系數(shù)絕對值的最大值。Ν為總的迭代次數(shù)，ε為接近零的小值，與數(shù)據(jù)中噪聲的能量有關(guān)。
[0032]本發(fā)明的優(yōu)點(diǎn)：本發(fā)明通過利用傅立葉變換將時間域三維地震數(shù)據(jù)變換到頻率 域，并采用凸集投影算法和能夠刻畫地震數(shù)據(jù)局部化特征的曲波變換進(jìn)行數(shù)據(jù)重建，從而 避免了低能量有效波同相軸系數(shù)被濾除，提高了缺失地震數(shù)據(jù)重建精度，同時在迭代過程 中提出新的按指數(shù)規(guī)律衰減的閾值參數(shù)，采用軟閾值算子，對每個頻率切片單獨(dú)進(jìn)行重建， 可以減少迭代次數(shù)，降低算法的計算復(fù)雜度，大幅度地提高了計算效率，節(jié)省了運(yùn)算時間。
【附圖說明】
[0033]圖1表示本發(fā)明實(shí)施例中頻率域三維不規(guī)則地震數(shù)據(jù)重建方法流程圖。
[0034] 圖2不同閾值參數(shù)重建信噪比曲線對比圖。
[0035] 圖3為二維隨機(jī)欠采樣及其重建結(jié)果對比圖。
[0036]圖4為重建前后二維頻譜對比圖。
[0037]圖5為信噪比與采樣率關(guān)系曲線對比圖。
[0038]圖6為含噪數(shù)據(jù)重建過程圖。
【具體實(shí)施方式】
[0039]以下實(shí)施案例用于說明本發(fā)明，但不用來限制本發(fā)明的范圍。
[0040] 實(shí)施例1
[0041] 實(shí)現(xiàn)該方法的步驟主要包括，重建方程的構(gòu)建，時間域變換到頻率域，凸集投影重 建算法，閾值參數(shù)處理等。具體步驟如下：
[0042] 步驟1:重建方程的構(gòu)建。地震數(shù)據(jù)的重建問題可以描述為由一組不完整數(shù)據(jù)通過 線性算子的作用恢復(fù)出完整的數(shù)據(jù)，假設(shè)如下線性正演模型
[0043] y〇bs=Mf
[0044] 這里yobseR"代表采集的地震數(shù)據(jù);feRN，且N2n，表示待重建的無假頻數(shù)據(jù);Me RnXN表示對角矩陣，其元素1和0分別表示已知地震道和未知地震道，假設(shè)系數(shù)X是f在曲波域 C中的稀疏表示，則上述方程可以寫成
[0045] y〇bs=Axiiα=μ€η
[0046]這里上標(biāo)H代表共輒轉(zhuǎn)置矩陣，當(dāng)從采集到的數(shù)據(jù)yQbs中重建無假頻數(shù)據(jù)f時，由于 X是稀疏的，從而使得該欠定方程有解。
[0047] 在稀疏促進(jìn)反演后，重建信號由7 = (:〃 ^決定，同時
[0048]
r
[0049] 在這個表達(dá)式中，^代表估計值，Li范數(shù)定義為|χ| ，x[i]是向量x中第i 個元素，通過求解上述欠定方程，原始完整而規(guī)則的地震數(shù)據(jù)就可以重建出來。
[0050] 步驟2:時間域變換到頻率域。本發(fā)明采用傅立葉變換對時間域地震數(shù)據(jù)進(jìn)行計 算，其離散傅立葉變換公式如下：
[0051]
[0052] 式中：m、n均為整數(shù)。m= 0，1，"·Μ-1;n= 0，1，"·Ν-1 為頻譜分析的時窗長 度。
[0053]步驟3:凸集投影重建算法。針對傳統(tǒng)重建算法重建精度低和計算速度慢的缺點(diǎn)， 提出在頻率域中進(jìn)行重建的算法，并采用軟閾值算子和新的按指數(shù)規(guī)律衰減的閾值參數(shù)公 式，其算法迭代表達(dá)式為：
[0054]
[0055]其中，fk表示第k次迭代得到的時間域重建數(shù)據(jù)，表示原始采集到數(shù)據(jù)y〇bs(t，x， y)，滿足/n，F(xiàn)t和F、表示關(guān)于時間變量t的正反傅立葉變換，c和(Γ1表示正反曲波變 換，其曲波變換的定義為：
[0056]
[0057] 式中：表示曲波函數(shù)，j，l，k分別表示尺度，方向和位置參數(shù)，f(x)為地震數(shù) 據(jù)，其頻率域定義式為：
[0058]
[0059]經(jīng)過變換后得到的曲波系數(shù)，可用(：{」}{1}仏1氺2)表示其結(jié)構(gòu)，其中」表示尺度，1 表示方向，（lu，k2)表示j尺度1方向上的矩陣系數(shù)。
[0060] Tk表示軟閾值算子，其元素滿足：
[0061]
[0062]Sk-i表示第k-Ι次迭代得到的重建數(shù)據(jù)