彭佳明 付麗華 張雪敏
(中國(guó)地質(zhì)大學(xué)(武漢)數(shù)學(xué)與物理學(xué)院,湖北武漢 430074)
從缺失的原始數(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)?;诘椭妊a(bǔ)全理論的重建方法的主要原理是:有限個(gè)線性同相軸的無(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ì)算效率。
地震數(shù)據(jù)重建的秩最小化模型為
(1)
式中:X∈Rm×n為重建的地震數(shù)據(jù),其中m、n為重建地震數(shù)據(jù)的維度,且有p=m×n; Rank(X)表示求X的秩; A是一個(gè)Rm×n→Rp的線性映射;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ù)。
式(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)解,則滿足
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)步驟
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ù)。
為了解決SVD算法復(fù)雜度較高的問(wèn)題,本文提出的FFPC算法采用隨機(jī)投影的分解方法,即將高維矩陣投影到低維矩陣,再進(jìn)行SVD運(yùn)算; 同時(shí)利用子空間復(fù)用的方法減少迭代次數(shù),降低計(jì)算復(fù)雜度,以進(jìn)一步提高運(yùn)行效率。
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算法步驟
在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è)置較小,只要滿足(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算法步驟
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),且滿足條件(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滿足以上條件時(shí),對(duì)比FFPC與FPC算法的復(fù)雜度,可知FFPC算法相對(duì)FPC具有更高的計(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ī)刹糠?,再在頻率域的每個(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ù)重建流程
測(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)
仿真數(shù)據(jù)由Ricker子波產(chǎ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)??芍诠潭ㄈ笔是闆r下,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ì)算效率。
選取實(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算法重建信噪比
選取相同尺度和特征參數(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算法重建信噪比
基于核范數(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)矩陣誤差滿足一定精度而自適應(yīng)地改變。