應(yīng)用快速不動(dòng)點(diǎn)延續(xù)算法的地震數(shù)據(jù)重建

彭佳明 付麗華 張雪敏

(中國(guó)地質(zhì)大學(xué)(武漢)數(shù)學(xué)與物理學(xué)院，湖北武漢 430074)

0 引言

從缺失的原始數(shù)據(jù)重建完整地震波場(chǎng)是經(jīng)典的反問(wèn)題。顯然，地震數(shù)據(jù)重建可為后續(xù)的偏移成像、全波形反演等處理提供可靠的基礎(chǔ)數(shù)據(jù)[1]。現(xiàn)今，地震數(shù)據(jù)的重建方法主要有以下六類(lèi)。

基于濾波的重建方法是通過(guò)褶積插值濾波器實(shí)現(xiàn)插值，較常用的有基于預(yù)測(cè)誤差濾波方法等[2]。此類(lèi)方法的優(yōu)點(diǎn)是不需反射同相軸的先驗(yàn)信息，不受假頻影響； 缺點(diǎn)是運(yùn)算量較大，且只適用于規(guī)則采樣的地震數(shù)據(jù)。

基于波場(chǎng)延拓算子的重建方法是通過(guò)Kirchhodff積分算子實(shí)現(xiàn)插值，主要有炮域延拓法[3]、炮檢距域延拓法[4-5]和傾角時(shí)差校正法[6]等。該類(lèi)方法的優(yōu)點(diǎn)是靈活度高，能最大程度地利用地下信息； 缺點(diǎn)是計(jì)算復(fù)雜度較高。

基于變換域重建的方法常采用Randon變換[7-8]、Fourier變換[9-10]、Curvelet變換[11-12]等，再在變換域完成地震數(shù)據(jù)的重建。該類(lèi)方法的優(yōu)點(diǎn)是計(jì)算速度快，可處理不同缺失形式的地震數(shù)據(jù)； 缺點(diǎn)是不適用于有限帶寬和含有空間假頻的地震數(shù)據(jù)。

基于相干傾角插值的重建方法主要包括最大相干傾角插值法[13]和多傾角同相軸方法[14]。其優(yōu)點(diǎn)是可自適應(yīng)地做傾角掃描； 缺點(diǎn)是不適用于含噪地震數(shù)據(jù)，對(duì)于高維地震數(shù)據(jù)的計(jì)算復(fù)雜度很高。

基于壓縮感知理論的地震數(shù)據(jù)重建[15-18]，其思想是利用變換域的稀疏性，將重建作為反問(wèn)題求解。其主要方法有正交匹配追蹤算法[19]、字典學(xué)習(xí)法[20-21]、光滑L1范數(shù)法[22]和Bregman迭代算法[23]等。它們的優(yōu)點(diǎn)是計(jì)算速度快，可處理大規(guī)模問(wèn)題； 缺點(diǎn)是受噪聲干擾影響較大。

近年來(lái)，作為壓縮感知的二維推廣形式，“矩陣補(bǔ)全”得到快速發(fā)展，基于“低秩補(bǔ)全”的地震數(shù)據(jù)重建也成為現(xiàn)今主流方法。其優(yōu)點(diǎn)是比現(xiàn)有方法更簡(jiǎn)單、高效； 缺點(diǎn)是須通過(guò)預(yù)變換建立低秩結(jié)構(gòu)。基于低秩補(bǔ)全理論的重建方法的主要原理是：有限個(gè)線(xiàn)性同相軸的無(wú)缺失或無(wú)噪聲污染的地震數(shù)據(jù)在頻率切片上是低秩結(jié)構(gòu)[24]，而數(shù)據(jù)的缺失和噪聲的存在導(dǎo)致地震數(shù)據(jù)構(gòu)造的Hankel矩陣的秩的增加，因此可通過(guò)降秩處理實(shí)現(xiàn)數(shù)據(jù)重建。

經(jīng)典的矩陣秩最小化是一個(gè)“NP難”(高復(fù)雜度)問(wèn)題。對(duì)此，人們嘗試將秩最小化近似轉(zhuǎn)化為核范數(shù)最小化問(wèn)題，再利用凸優(yōu)化法進(jìn)行求解。Trickett等[25]提出基于頻率切片降秩的方法。Oropeza等[26]證明地震時(shí)頻切片形成的塊Hankel矩陣是低秩矩陣。Yang等[27]梳理了矩陣補(bǔ)全與地震數(shù)據(jù)重建的關(guān)系，針對(duì)地震數(shù)據(jù)重建提出核范數(shù)最小化模型，且說(shuō)明了低秩矩陣補(bǔ)全的可靠性。Cai等[28]提出的奇異值閾值(Singular Value Thresholding， SVT)算法，主要是利用交替更新和梯度下降法求解核范數(shù)最小化問(wèn)題，而SVT的加速算法可參見(jiàn)文獻(xiàn)[29-30]。Goldfarb等[31]通過(guò)分離變量后，利用參數(shù)交替迭代更新和軟閾值收縮算子求解，實(shí)現(xiàn)了不動(dòng)點(diǎn)延續(xù)(Fixed Point Continuation，F(xiàn)PC)算法，且提出利用迭代方式更新觀測(cè)數(shù)據(jù)，即FPC與Bregman迭代[32]相結(jié)合的算法。Toh等[33]通過(guò)利用目標(biāo)函數(shù)的二階泰勒逼近函數(shù)，再利用不同變量的相互交替更新，實(shí)現(xiàn)了加速近端梯度(Accelerated Proximal Gradient，APG)算法； 并給出了改進(jìn)的加速算法(Accelerated Proximal Gradient by Linesearch-like technique，APGL)。

FPC算法在每次迭代更新時(shí)都需進(jìn)行奇異值分解(Singular Value Decomposition，SVD)，因此會(huì)降低算法效率。傳統(tǒng)解決方法是利用PROPACK加速包，該方法能對(duì)矩陣進(jìn)行快速SVD。為了進(jìn)一步提高計(jì)算效率，本文提出一種快速不動(dòng)點(diǎn)延續(xù)(Fast Fixed Pointed Continuation，F(xiàn)FPC)算法，利用塊克雷洛夫迭代近似奇異值分解算法(Block Krylov Iteration Methods for Approximate SVD，BKIMASVD)和子空間復(fù)用技術(shù)，在保證重構(gòu)精度的同時(shí)可大幅度縮短運(yùn)行時(shí)間。BKIMASVD是基于隨機(jī)投影形式實(shí)現(xiàn)的，其計(jì)算復(fù)雜度遠(yuǎn)小于PROPACK中所采用的快速算法；子空間復(fù)用技術(shù)則通過(guò)減少迭代次數(shù)，大幅度縮短運(yùn)行時(shí)間。通過(guò)以上兩種技術(shù)的運(yùn)用，相比傳統(tǒng)加速算法，F(xiàn)FPC有更高的計(jì)算效率。

1 基于FPC算法的地震數(shù)據(jù)重建

1.1 基于低秩理論的地震數(shù)據(jù)重建模型

地震數(shù)據(jù)重建的秩最小化模型為

(1)

式中：X∈Rm×n為重建的地震數(shù)據(jù)，其中m、n為重建地震數(shù)據(jù)的維度，且有p=m×n； Rank(X)表示求X的秩； A是一個(gè)Rm×n→Rp的線(xiàn)性映射；b∈Rp為缺失的地震數(shù)據(jù)列向量化結(jié)果。

式(1)的求解是個(gè)“NP難”問(wèn)題，因此通常可將目標(biāo)函數(shù)凸松弛為核范數(shù)最小化進(jìn)行求解。

假設(shè)矩陣X有r個(gè)正奇異值σ1≥σ2≥…≥σr>0，則X的核范數(shù)定義為

則基于核范數(shù)最小化的地震數(shù)據(jù)重建的模型如下

(2)

式中μ為正則化參數(shù)。

1.2 FPC重建算法的步驟

式(2)目標(biāo)函數(shù)關(guān)于X的次梯度為

μ?‖X‖*+g(X)

(3)

式中： ?‖X‖*為核范數(shù)的導(dǎo)數(shù)；g(X)=A*[A(X)-b]，且A*是A的“逆運(yùn)算”。

假設(shè)重構(gòu)矩陣X*為式(2)的最優(yōu)解，則滿(mǎn)足

0∈τμ?‖X*‖*+τg(X*)

=τμ?‖X*‖*+X*-[X*-τg(X*)]

(4)

式中變量τ>0。定義Y*=X*-τg(X*)，式(4)等價(jià)于

0∈τμ?‖X*‖*+X*-Y*

(5)

對(duì)式(5)關(guān)于X*求積分，可知X*也是下式最優(yōu)解

(6)

這樣，式(2)的秩最小化問(wèn)題即可利用文獻(xiàn)[32]中定理3求解。表1給出了具體的FPC算法步驟。

表1 不動(dòng)點(diǎn)延續(xù)算法(FPC)步驟

1.3 FPC算法復(fù)雜度分析

FPC算法中最重要的部分為循環(huán)體中的第2～第4步。假設(shè)循環(huán)終止迭代次數(shù)為M，由于第2步的截?cái)郤VD主要是針對(duì)稀疏矩陣Y操作的，再假設(shè)Y∈Rm×n中非零個(gè)數(shù)為Φ，選取的奇異值個(gè)數(shù)為k，則第2步復(fù)雜度為cSVDk|Φ|=O[mnmin(m,n)]，其中cSVD為標(biāo)準(zhǔn)SVD復(fù)雜度常數(shù)； 利用PROPACK加速包，可將計(jì)算復(fù)雜度降低為O(kmn)。計(jì)算第3步的復(fù)雜度為O[m(k+n)(2k-1)]，第4步的復(fù)雜度為O(1)。

綜上所述，F(xiàn)PC算法的平均復(fù)雜度為

kO[mnmin(m,n)]+kO[m(k+n)(2k-1)]+kO(1)

=M{cSVD[mnmin(m,n)]+

cmul[m(k+n)(2k-1)+1]}

(7)

傳統(tǒng)加速型FPC算法的平均復(fù)雜度為

M{cSVD(kmn)+cmul[m(k+n)(2k-1)+1]}

(8)

式中cmul為矩陣乘法的復(fù)雜度常數(shù)。

2 FFPC算法

為了解決SVD算法復(fù)雜度較高的問(wèn)題，本文提出的FFPC算法采用隨機(jī)投影的分解方法，即將高維矩陣投影到低維矩陣，再進(jìn)行SVD運(yùn)算； 同時(shí)利用子空間復(fù)用的方法減少迭代次數(shù)，降低計(jì)算復(fù)雜度，以進(jìn)一步提高運(yùn)行效率。

2.1 塊克雷洛夫迭代近似奇異值分解算法

BKIMASVD算法的主要思想是基于RandQB過(guò)程[34-35]，將需要進(jìn)行SVD的矩陣A降維分解，通過(guò)對(duì)維度更低的矩陣B進(jìn)行SVD以達(dá)到快速分解的目的。RandQB過(guò)程可表示為

A≈QBA∈Rm×n,Q∈Rm×l,B∈Rl×n

(9)

式中l(wèi)?min(m,n)。

BKIMASVD的具體步驟為：首先對(duì)觀測(cè)矩陣A隨機(jī)采樣得到采樣矩陣AΩ，其中Ω∈Rn×(k+s)為隨機(jī)產(chǎn)生的矩陣； 將采樣矩陣AΩ通過(guò)QR分解(將一個(gè)矩陣分解為一個(gè)正交矩陣與上三角矩陣的乘積)得到正交矩陣H0，在P次循環(huán)迭代過(guò)程中通過(guò)[Hi,～]=QR[A(ATHi-1),0]迭代方法逐步更新子空間Hi,i∈{1,2,…,P}； 再將迭代產(chǎn)生的子空間Hi按列依次排列，擴(kuò)充為更高維度的空間H； 然后對(duì)高維空間H進(jìn)行QR分解得到子空間矩陣Q，這樣便得到低維矩陣B=QTA； 最后對(duì)B進(jìn)行SVD即可提高其計(jì)算效率。迭代參數(shù)P可提升子空間維度，提高重建精度[34]，詳細(xì)原理參見(jiàn)文獻(xiàn)[35-38]。表2為BKIMASVD的具體算法步驟。

表2 BKIMASVD算法步驟

2.2 子空間復(fù)用

在FPC算法中，每個(gè)觀測(cè)矩陣需進(jìn)行多次SVD迭代更新。事實(shí)上，當(dāng)FPC算法執(zhí)行多次迭代后，迭代更新的矩陣Yt會(huì)趨于穩(wěn)定，因此在后續(xù)對(duì)Yt進(jìn)行迭代更新時(shí)，可復(fù)用第P次迭代產(chǎn)生的子空間Q，再對(duì)B=QTYt進(jìn)行SVD[34]。當(dāng)FPC算法中迭代次數(shù)t超過(guò)P后，復(fù)用技術(shù)可省去高維空間H和Q的迭代計(jì)算，從而有效縮短運(yùn)行時(shí)間。

具體操作步驟為：當(dāng)?shù)螖?shù)小于P時(shí)，利用BKIMASVD對(duì)迭代更新的矩陣進(jìn)行分解； 當(dāng)?shù)螖?shù)超過(guò)P后，保留第P次產(chǎn)生的Q，直接從BKIMASVD算法(表2中)的第6步開(kāi)始執(zhí)行，這樣就可保證在高精度范圍內(nèi)提高計(jì)算效率。迭代中可將P設(shè)置較小，只要滿(mǎn)足(k+s)P?min(m,n)就可保證多次循環(huán)迭代時(shí)SVD的計(jì)算量大幅度降低，也就可在確保精度的情況下顯著縮短運(yùn)行時(shí)間。

本文提出的FFPC算法是以BKIMASVD替代SVD，同時(shí)利用子空間復(fù)用技術(shù)減少迭代更新次數(shù)，可極大地提高計(jì)算效率。表3給出了FFPC算法的具體步驟。

表3 FFPC算法步驟

2.3 FFPC算法復(fù)雜度分析

FFPC與FPC算法最大的不同點(diǎn)在第2步，因此需首先分析BKIMASVD的復(fù)雜度。對(duì)于矩陣A∈Rm×n，BKIMASVD的第3步(表2)中的每一次循環(huán)迭代所需計(jì)算量為2cmul(k+s)Nnnz+cQRm(k+s)2，則循環(huán)P次的計(jì)算量為

2cmul(k+s)PNnnz+cQRm(k+s)2P

(10)

生成Q的計(jì)算量為cQRm(kP+sP)2，生成B的計(jì)算量為

cmul(k+s)PNnnz

(11)

對(duì)B進(jìn)行SVD的計(jì)算量為

cSVD[m(kP+sP)min(m,kP+sP)]

(12)

第8步(表2)計(jì)算量為cmul[(k+s)P]2m， 第9步(表2)計(jì)算量為cmul[m(k+n)(2k-1)]。因此，BKIMASVD的復(fù)雜度為

cSVDm(kP+sP)min(m,kP+sP)+

cmulm(k+s)(2k-1)

(13)

式中：Nnnz表示矩陣A中的非零元素的個(gè)數(shù)；cQR表示QR分解的復(fù)雜度常數(shù)。且(k+s)P?min(m,n)，而傳統(tǒng)SVD算法的復(fù)雜度為cSVD[mnmin(m,n)]=O[mnmin(m,n)]，所以利用該算法可實(shí)現(xiàn)快速SVD。

同樣，假設(shè)循環(huán)終止迭代次數(shù)為M，因第2步(表3)的SVD主要針對(duì)稀疏矩陣Yt，假設(shè)Yt∈Rm×n中非零個(gè)數(shù)為Φ，選取奇異值個(gè)數(shù)為k，子空間復(fù)用個(gè)數(shù)為P，則第2步(表2)主要復(fù)雜度為

(kP+sP)2m+cSVDm(kP+sP)×

min(m,kP+sP)+cmulm(k+s)(2k-1)]

(14)

當(dāng)次數(shù)超過(guò)P以后，由于復(fù)用前面計(jì)算的Q，所以只需計(jì)算B的分解，那么循環(huán)所需時(shí)耗為

(M-P)[cmul(k+s)PNnnz+

cSVDm(kP+sP)min(m,kP+sP)+

cmul(kP+sP)2m+cmulm(k+s)(2k-1)]

(15)

則FFPC算法的復(fù)雜度為

3cmulkP2Nnnz+(cQR+cQRP)m(kP+sP)2+

M[cSVDm(kP+sP)min(m,kP+sP)+

cmulm(kP+sP)2+cmulm(k+n)(2k-1)]

(16)

當(dāng)矩陣維度m、n較大、P很小(P?M)，且滿(mǎn)足條件(k+s)P?min(m,n)時(shí)， FFPC算法的復(fù)雜度近似為

M[cSVDm(kP+sP)min(m,kP+sP)+

cmulm(kP+sP)2+cmulm(k+n)(2k-1)]

=M[(cSVD+cmul)m(kP+sP)2+

cmulm(k+n)(2k-1)]

(17)

而傳統(tǒng)加速型FPC算法的復(fù)雜度可近似為

M[cSVDkmn+cmulm(k+n)(2k-1)]

(18)

當(dāng)P、k、s滿(mǎn)足以上條件時(shí)，對(duì)比FFPC與FPC算法的復(fù)雜度，可知FFPC算法相對(duì)FPC具有更高的計(jì)算效率。

3 基于FFPC算法的地震數(shù)據(jù)重建

理論證明了地震數(shù)據(jù)在頻率切片上具有低秩結(jié)構(gòu)特性，但地震數(shù)據(jù)的缺失會(huì)導(dǎo)致其構(gòu)造Hankel矩陣的秩增加，于是缺失位置信息就可通過(guò)降秩進(jìn)行重建。針對(duì)地震數(shù)據(jù)，圖1給出了基于FPC和FFPC算法的重建流程圖。具體操作過(guò)程如下。

首先對(duì)缺失觀測(cè)地震數(shù)據(jù)做一維傅里葉變換，得到實(shí)部和虛部?jī)刹糠郑僭陬l率域的每個(gè)切片上構(gòu)造Hankel矩陣； 對(duì)于每個(gè)Hankel矩陣，分別利用FPC和FFPC算法進(jìn)行降秩重建； 針對(duì)重建的Hankel矩陣用反對(duì)角取平均，即可得重建的頻率切片； 最后將重建的實(shí)部和虛部頻率切片通過(guò)逆傅里葉變換到時(shí)間域，即完成了地震數(shù)據(jù)重建(圖1)。

圖1 基于FPC和FFPC算法的地震數(shù)據(jù)重建流程

4 實(shí)驗(yàn)

測(cè)試實(shí)驗(yàn)選用兩個(gè)主要衡量指標(biāo)，即信噪比和運(yùn)算耗時(shí)。信噪比的計(jì)算公式為

(19)

式中：Ps是信號(hào)能量；Pn是噪聲能量；I為原始信號(hào)；In為噪聲信號(hào)。此處運(yùn)算耗時(shí)為重建過(guò)程所消耗的CPU 時(shí)間。采用重建信噪比相對(duì)誤差表征不同算法重建信噪比的差異，其計(jì)算表達(dá)式為

(20)

4.1 仿真地震數(shù)據(jù)實(shí)驗(yàn)

仿真數(shù)據(jù)由Ricker子波產(chǎn)生，它包含兩條交叉的線(xiàn)性同相軸，仿真地震數(shù)據(jù)(圖2a)尺寸為256(道)×256(個(gè))。對(duì)比原始仿真地震數(shù)據(jù)(圖2a)、隨機(jī)缺失50%數(shù)據(jù)(圖2b)及其對(duì)應(yīng)的頻率—波數(shù)圖(圖3a、圖3b)，可見(jiàn)圖3b上有很明顯的假頻； 從基于FFPC與FPC算法的重建結(jié)果(圖2c、圖2d)及其對(duì)應(yīng)頻率—波數(shù)圖(圖3c、圖3d)上，清晰可見(jiàn)抗假頻重構(gòu)特征。當(dāng)缺失率為50%時(shí)(表4)，F(xiàn)FPC與FPC的重建信噪比相對(duì)誤差為0.38%，但FFPC算法耗時(shí)為FPC算法的1/1.66。顯然，針對(duì)仿真地震數(shù)據(jù)，在確保高精度重建效果的同時(shí)， FFPC算法比FPC具有更高的計(jì)算效率。

為了更詳實(shí)地對(duì)比不同缺失率下FFPC和FPC算法的重建效果，分別測(cè)試了缺失率為10%～80%情況下基于兩算法的重建信噪比和耗時(shí)(表4、圖4和圖5)。可知在固定缺失率情況下，F(xiàn)FPC與FPC算法重建信噪比誤差范圍為0.38%～7.00%，而FFPC算法的重建耗時(shí)為FPC算法的1/3.53～1/1.59。FFPC和FPC算法的重建信噪比都隨缺失率的增加而降低； 重建耗時(shí)也隨缺失率的增加總體上逐漸降低，且FFPC算法提升的倍率為先降后升(以60%缺失率為界)。當(dāng)缺失率≤50%時(shí)，F(xiàn)PC算法耗時(shí)較長(zhǎng)； 而利用FFPC算法能在確保高精度重建的同時(shí)，顯著縮短運(yùn)行時(shí)間，提高計(jì)算效率。

4.2 疊前地震數(shù)據(jù)實(shí)驗(yàn)

選取實(shí)際疊前地震數(shù)據(jù)(256(道)×256(個(gè)))進(jìn)行實(shí)驗(yàn)。對(duì)比原始疊前數(shù)據(jù)(圖6a)、隨機(jī)缺失50%數(shù)據(jù)(圖6b)及其對(duì)應(yīng)頻率—波數(shù)圖(圖7a、圖7b)、基于FFPC和FPC算法的重建結(jié)果(圖6c、圖6d)及其對(duì)應(yīng)頻率—波數(shù)圖(圖7c、圖7d)，可見(jiàn)基于FFPC的重建效果與FPC基本相當(dāng)，都可較好地對(duì)隨機(jī)缺失數(shù)據(jù)進(jìn)行插值。當(dāng)缺失率為50%時(shí)(表5)，F(xiàn)FPC與FPC的重建信噪比相對(duì)誤差為0.58%； 但FPC算法的耗時(shí)是FFPC的3.83倍，表明針對(duì)疊前數(shù)據(jù)，在確保高精度重建效果的同時(shí)，F(xiàn)FPC計(jì)算效率更高。

圖2 缺失率為50%的仿真地震數(shù)據(jù)重建效果對(duì)比(a)原始數(shù)據(jù)； (b)隨機(jī)缺失50%的數(shù)據(jù)； (c)基于FFPC重建的數(shù)據(jù)； (d)基于FPC重建的數(shù)據(jù)

圖3 缺失率為50%時(shí)仿真地震數(shù)據(jù)頻率—波數(shù)對(duì)比(a)原始數(shù)據(jù)； (b)隨機(jī)缺失50%的數(shù)據(jù)； (c)基于FFPC重建的數(shù)據(jù)； (d)基于FPC重建的數(shù)據(jù)

圖4 不同缺失率下仿真數(shù)據(jù)FPC與FFPC算法重建耗時(shí)

圖5 不同缺失率下仿真數(shù)據(jù)FPC與FFPC算法重建信噪比

圖6 缺失率為50%的疊前地震數(shù)據(jù)重建效果對(duì)比(a)原始數(shù)據(jù)； (b)隨機(jī)缺失50%的數(shù)據(jù)； (c)基于FFPC重建的數(shù)據(jù)； (d)基于FPC重建的數(shù)據(jù)

表4 不同缺失率下仿真數(shù)據(jù)FPC與FFPC算法性能對(duì)比

表5 不同缺失率下疊前數(shù)據(jù)FPC與FFPC算法性能對(duì)比

表5給出缺失率為10%～80%的疊前地震數(shù)據(jù)的重建信噪比和耗時(shí)(圖8、圖9)。可見(jiàn)在固定缺失率時(shí)，F(xiàn)FPC與FPC算法重建信噪比的相對(duì)誤差范圍是0.37%～3.80%，即二者重建的精度接近； 而FPC算法的重建耗時(shí)卻是FFPC算法的2.37～5.88倍。隨著缺失率的增加，F(xiàn)FPC和FPC算法的重建信噪比都降低，二者的重建耗時(shí)也逐漸降低。當(dāng)缺失率≤50%時(shí)，F(xiàn)PC算法耗時(shí)很長(zhǎng)，而FFPC算法在確保高精度重建的前提下能顯著縮短運(yùn)行時(shí)間，具有更高計(jì)算效率。

圖7 缺失率為50%時(shí)疊前地震數(shù)據(jù)頻率—波數(shù)對(duì)比(a)原始數(shù)據(jù)； (b)隨機(jī)缺失50%的數(shù)據(jù)； (c)基于FFPC重建的數(shù)據(jù)； (d)基于FPC重建的數(shù)據(jù)

圖8 不同缺失率下疊前數(shù)據(jù)FPC和FFPC算法重建耗時(shí)

圖9 不同缺失率下疊前數(shù)據(jù)FPC和FFPC算法重建信噪比

4.3 疊后地震數(shù)據(jù)實(shí)驗(yàn)

選取相同尺度和特征參數(shù)的實(shí)際疊后地震數(shù)據(jù)做進(jìn)一步測(cè)試。得到并對(duì)比原始疊后數(shù)據(jù)(圖10a)、隨機(jī)缺失50%的數(shù)據(jù)(圖10b)、基于FFPC(圖10c)和FPC(圖10d)的重建結(jié)果及其對(duì)應(yīng)的頻率—波數(shù)圖(圖11a～圖11d)，可見(jiàn)基于FFPC和FPC的重建結(jié)果對(duì)抑制假頻都起了較好作用。當(dāng)缺失率為50%時(shí)(表6)，F(xiàn)FPC與FPC重建信噪比誤差為0.15%，但FFPC算法的重建耗時(shí)是FPC的1/1.75。

表6給出缺失率為10%～80%的實(shí)際疊后地震數(shù)據(jù)的重建信噪比和耗時(shí)(圖12、圖13)。可見(jiàn)在固定缺失率時(shí)，基于FFPC與FPC算法的重建信噪比的相對(duì)誤差范圍是0.15%～5.20%； 而FPC算法的重建耗時(shí)是FFPC的1.75～2.99倍。FFPC和FPC算法的重建信噪比都隨缺失率的增加而降低； 重建耗時(shí)也都隨缺失率的增加而逐漸降低，且FFPC算法提升的倍率為先降后升。當(dāng)缺失率≤50%時(shí)，F(xiàn)PC算法耗時(shí)很長(zhǎng)，但利用FFPC算法卻能在確保高精度重建效果的同時(shí)，有效縮短運(yùn)行時(shí)間，提高計(jì)算效率。

因此，基于BKIMASVD算法與子空間復(fù)用的FFPC算法比傳統(tǒng)加速型的FPC算法的效率至少可提升1.75倍。

圖10 缺失率為50%的疊后地震數(shù)據(jù)重建效果對(duì)比(a)原始數(shù)據(jù)； (b)隨機(jī)缺失50%的數(shù)據(jù)； (c)基于FFPC重建的數(shù)據(jù)； (d)基于FPC重建的數(shù)據(jù)

圖11 缺失率為50%的疊后地震數(shù)據(jù)的頻率—波數(shù)對(duì)比(a)原始數(shù)據(jù)； (b)隨機(jī)缺失50%的數(shù)據(jù)； (c)基于FFPC重建的數(shù)據(jù)； (d)基于FPC重建的數(shù)據(jù)

表6 不同缺失率下疊后數(shù)據(jù)FPC與FFPC算法性能對(duì)比

圖12 不同缺失率下疊后數(shù)據(jù)FPC與FFPC算法重建耗時(shí)

圖13 不同缺失率下疊后數(shù)據(jù)FPC 與FFPC算法重建信噪比

5 結(jié)論與展望

基于核范數(shù)最小化的FPC算法在每次迭代過(guò)程中需進(jìn)行SVD，導(dǎo)致算法的計(jì)算復(fù)雜度較高； 傳統(tǒng)的加速方法為直接利用PROPACK包，但仍然耗時(shí)較長(zhǎng)。本文提出利用BKIMASVD算法和子空間復(fù)用的FFPC算法，并分析了該算法的計(jì)算復(fù)雜度。仿真數(shù)據(jù)和實(shí)際地震數(shù)據(jù)的實(shí)驗(yàn)表明：基于FFPC算法與基于FPC算法的地震數(shù)據(jù)重建精度基本相當(dāng)，都可較好地對(duì)隨機(jī)缺失地震數(shù)據(jù)進(jìn)行插值，實(shí)現(xiàn)抗假頻重構(gòu)，但本文提出的FFPC算法進(jìn)一步提高了計(jì)算效率。

實(shí)驗(yàn)過(guò)程中還發(fā)現(xiàn)：對(duì)于m≤n類(lèi)型的矩陣，F(xiàn)FPC算法提速效果更明顯。由于算法所涉及的參數(shù)較多，它們的選取和設(shè)定都會(huì)一定程度地影響重建效果。其中冪迭代次數(shù)P可通過(guò)前后兩次迭代的重構(gòu)矩陣誤差滿(mǎn)足一定精度而自適應(yīng)地改變。