基于壓縮感知的信號重建方法及在氣槍震源信號處理中的應(yīng)用

張 帥, 楊潤海, 高爾根

(1. 云南省地震局, 云南 昆明 650224; 2. 安徽建筑大學(xué), 安徽 合肥 230000)

0 引言

在地震孕育和發(fā)生的過程中往往伴隨著地下介質(zhì)狀態(tài)的變化[1]。地震波穿過其中往往能夠帶回介質(zhì)變化信息。為此,2012年全球第一個陸地大容量氣槍源發(fā)射站在大理賓川建成[2],利用氣槍震源具有能量轉(zhuǎn)換效率高,信號高度重復(fù)性等優(yōu)勢,通過重復(fù)測量地震波走時來監(jiān)測地球內(nèi)部介質(zhì)波速微弱變化,為研究地震孕育和發(fā)生過程及尋找地震發(fā)生前兆提供可能[3]。

但在氣槍源探測過程中,受各種干擾因素的影響,部分信號常?；烊腚S機干擾信號或者數(shù)據(jù)缺失。傳統(tǒng)濾波方法雖然可以對其識別,但無法有效濾除[4],常用的處理手段是人為剔除這些受隨機信號干擾的數(shù)據(jù)道(這一處理過程類似人工勘探中對壞道進行的道編輯處理),但會導(dǎo)致數(shù)據(jù)隨機缺失,影響后續(xù)走時計算、結(jié)構(gòu)成像等處理研究。為了得到連續(xù)完整數(shù)據(jù),需要一種高精度非線性的插值方法對缺失數(shù)據(jù)進行插值重建。傳統(tǒng)方法大都受 Nyquist 采樣理論的限制,對于不滿足采樣定理的超稀疏地震數(shù)據(jù),插值效果不佳[5]。而壓縮感知理論表明:即使是隨機缺失的數(shù)據(jù),也有可能恢復(fù)出滿足一定精度要求的完整數(shù)據(jù)[6]。

壓縮感知方法發(fā)展起步較晚,第一篇有關(guān)壓縮感知的完整論述是由Cande’s等[7]于2006年才正式發(fā)表。隨著近幾年的不斷發(fā)展,其已經(jīng)成功應(yīng)用于多個地震學(xué)和地球物理反演問題中去[8],包括地震勘探領(lǐng)域的信號恢復(fù)和重建[9-11],大地震破裂和能量釋放過程研究[12-13],地震信號分析和處理[14]等領(lǐng)域。

壓縮感知理論假設(shè):若一個信號在某一個變換域中是稀疏的,那么通過一個與稀疏變換基函數(shù)不相關(guān)的測量矩陣將該信號由高維空間投影至低維空間,可以得到一個遠小于信號長度的觀測信號[15],再通過某種重構(gòu)算法就可以將低維信號重建到原始高維信號。

本文基于傅里葉變換域構(gòu)建一種基于壓縮感知理論的氣槍源缺失信號重建方法(以下簡稱CS),并與傳統(tǒng)的三次樣條插值方法(以下簡稱Spline)進行對比分析。數(shù)值模擬和實際數(shù)據(jù)處理表明,該方法能夠有效地重建因受干擾缺失的氣槍源信號。

1 壓縮感知氣槍信號重建方法(CS)

重構(gòu)過程涉及稀疏變換,采樣方法(測量矩陣),重構(gòu)算法三個方面,大部分研究都是基于這三個方面開展[16],①稀疏變換,如ZHANG等[17]、Herrmann等[18-19]分別基于傅里葉變換和曲波變換對缺失地震數(shù)據(jù)進行恢復(fù)重建研究,顯示基于曲波變換的重建效果更好。②采樣方法(測量矩陣),規(guī)則采樣方法容易引入假頻,而隨機采樣方法可以將混淆真實頻率的相干噪聲轉(zhuǎn)化成容易濾除的非相干噪聲[20],避免引入假頻(如圖1)。Hennenfent等[21]、唐剛等[22]將jitter隨機采樣方法引入到地震數(shù)據(jù)重建中,能夠很好地控制采樣間隔。③重構(gòu)算法,實際是通過不斷的稀疏促進策略解去解決矩陣的最優(yōu)化問題。HALE等[23]將FPC重構(gòu)算法應(yīng)用于數(shù)據(jù)重建中,證明了其收斂性能優(yōu)于傳統(tǒng)算法。MA等[24]、LIU等[25]將不動點迭代與Bregman迭代相結(jié)合求解矩陣最小化問題。有效信號在稀疏域中系數(shù)很稀疏,噪聲的系數(shù)很稠密,因此,在通過最小化的非零系數(shù)對原系數(shù)進行估計最小化估計過程中,同時起到了對噪聲的壓制作用。

圖1 規(guī)則采樣與隨機采樣頻譜對比示意圖Fig.1 Frequency spectra comparison between regular sampling and random sampling

1.1 方法原理

依據(jù)壓縮感知理論:若一個信號在某一個變換域中是稀疏的,那么就可以通過一個遠小于信號長度的觀測信號,利用一定的重構(gòu)算法可以低維信號重建到高維信號。本文根據(jù)壓縮感知理論這一特點,基于傅里葉稀疏變換域構(gòu)建了一種基于壓縮理論的缺失氣槍信號重建方法。其方法原理及重構(gòu)過程如下:

假設(shè)不完整氣槍信號x(n),它的重構(gòu)問題可以通過以下的模型進行表示

x(n)=Rvec[y(n)]

(1)

(2)

式中:x(n)表示通過測量矩陣得到的不完整的地震數(shù)據(jù);vec表示二維的向量化算子;y(n)表示二維完整地震數(shù)據(jù)。R∈RPXMN(P?MN)是根據(jù)缺失數(shù)據(jù)的定義的單位測量隨機矩陣;O為缺失數(shù)據(jù)位置,它的值為0;I為未缺失位置數(shù)據(jù)的位置,它的值為1[11]。

根據(jù)CS理論(Compressive Sensing),如果y(n)在某一個變換域Ψ中是稀疏的。那么其變換系數(shù)為

Θ=ΨTy

(3)

Θ是Ψ的等價或者稀疏表示。根據(jù)缺失的位置設(shè)計一個平穩(wěn)與變換基Ψ不相關(guān)M*N的維的測量矩陣φ對Θ進行觀測得到觀測數(shù)據(jù):

x(n)=φΘ=φΨTy

(4)

令α=ΨTy,因此得到:

x(n)=φα

(5)

式(5)為欠定方程,方程解多個,但是α具有稀疏性,我們可以通過不斷促進α的稀疏性進行反演求解。這就為實現(xiàn)缺失數(shù)據(jù)的精確重構(gòu)成為可能。其中的稀疏促進策略如下[18]:

(6)

輸入:測量矩陣φ,測量信號向量y,

(1) 初始化殘差r0=y;初始支撐集F0=?;迭代次數(shù)k=1;終止最大迭代次數(shù)P=15;

(2) 找到索引:

λk={j:|gk(j)|>tσk}

(7)

式中:

gk=φTrk-1

(8)

gk(j)為向量gk中第k元素。

(9)

式中:σk為正常噪聲水平;t為閾值參數(shù)。

(3) 由最小二乘得到:

(10)

(4) 更新殘差:

(11)

1.2 信號重建處理流程

基于傅里葉變換的壓縮感知缺失氣槍信號重建方法(CS)處理流程如下。

圖2 重構(gòu)流程圖Fig.2 Reconstruction flow chart

2 數(shù)值模擬

為了驗證方法的有效性,利用信噪比值Rsn及均方根誤差Mse來衡量缺失地震數(shù)據(jù)的重建效果

(12)

(13)

2.1 一維數(shù)值模擬

首先模擬一維氣槍信號對上述方法進行測試。該數(shù)據(jù)采樣頻率為1 000 Hz,采樣時間為0.15 s,共150個采樣點,如圖3(a)所示。對模擬信號進行20%隨機缺失,如圖3(b)所示。利用本文方法(CS)及三次樣條插值方法(Spline)分別對缺失數(shù)據(jù)進行重構(gòu),并進行效果對比,如圖3(c)所示。

圖3 兩種方法重建效果對比Fig.3 Comparison between reconstruction effect of two methods

考慮到噪聲對重構(gòu)的影響,本文對模擬信號加入一定量的隨機噪聲如圖3(d)所示。同樣進行20%的數(shù)據(jù)缺失如圖3(e)所示。利用上述兩種方法分別進行重建處理;重建結(jié)果如圖3(f)所示。分別對模擬加噪聲與不加噪聲氣槍信號進行了缺失重建效果分析,并與傳統(tǒng)的三次樣條插值方法處理結(jié)果進行對比,從信噪比(SNR)、均方根誤差、互相關(guān)三個方面對重建效果進行對比分析,結(jié)果如表1,表2所列。表明CS方法重構(gòu)后的波形誤差更小、波形一致性更強、信噪比更高,重建效果好于傳統(tǒng)的三次樣條插值方法。

表1 不加噪聲重構(gòu)效果對比分析

表2 加噪聲重構(gòu)效果對比分析

2.2 二維數(shù)值模擬

模擬二維氣槍信號對本文方法做進一步測試。該二維信號是將某一維模擬信號按照類正弦曲線排列得到,采樣頻率為100 Hz,采樣時間5 s,共500個采樣點,100道槍信號,如圖4(a)所示。二維信號的缺失重建問題需轉(zhuǎn)換到一維中處理,也就是從原來的時間域采樣變換到空間域采樣。

圖4 模擬二維信號缺失重建Fig.4 Simulation of missing and reconstruction of two-dimensional signal

同樣考慮噪聲對重構(gòu)的影響,模擬信號中加入一定量的隨機噪聲,如圖4(b)所示。將模擬信號隨機加入20%的干擾信號,如圖4(c)所示。通過設(shè)計測量矩陣(采集到的設(shè)置為1,未采集到的設(shè)置為0),對受隨機干擾的數(shù)據(jù)進行相應(yīng)的缺失,如圖4(d)所示,再利用本文方法重構(gòu)出滿足一定精度的有效信號。

利用本文方法重建后的效果如圖5(a)所示,觀察到重構(gòu)后的氣槍信號同相軸清晰且連續(xù)性相同,受干擾缺失的數(shù)據(jù)得到很好的恢復(fù),對隨機噪聲的壓制較好,重構(gòu)精度較高。而三次樣條插值的重構(gòu)結(jié)果如圖5(b)所示,觀察到由于對噪聲的壓制效果較差,導(dǎo)致重構(gòu)精度較低,重構(gòu)效果不佳。對兩種方法的重構(gòu)效果做進一步的信噪比及誤差分析如圖6所示。可以觀察到本文方法重建信號的信噪比值遠遠高于三次樣條插值重構(gòu)信號的信噪比值,同樣本文方法重構(gòu)前后信號的誤差也整體小于三次樣條插值方法重構(gòu)前后信號的誤差。也說明本文方法重構(gòu)效果要好于傳統(tǒng)的三次樣條插值方法。但同樣也觀察到部分連續(xù)缺失的信號,重建后信號的信噪比值要相對較小、誤差相對較大,說明連續(xù)缺失對重構(gòu)精度有一定的影響。

圖5 兩種方法重建效果對比Fig.5 Comparison between reconstruction effect of two methods

圖6 兩種方法重建效果分析Fig.6 Analysis of reconstruction effect of two methods

3 實際資料處理

選取距離賓川氣槍震源發(fā)射站大概27 km的高興臺站數(shù)據(jù)進行處理。該數(shù)據(jù)有1 000道槍信號,采樣頻率為100 Hz,采樣時間12 s。由于臺站距離相對較遠及受各種干擾的影響,導(dǎo)致部分有效信號受隨機信號干擾嚴重如圖7(a)所示。并且隨機干擾信號分布具有隨機性。因此,將所有槍信號進行線性疊加得到疊加信號,再將疊加信號與每個槍信號做互相關(guān),得到每個槍信號與疊加信號的互相關(guān)系數(shù),找出相關(guān)系數(shù)低于0.7的槍信號,根據(jù)計算,互相關(guān)系數(shù)小于0.7的槍信號有211道,依據(jù)它的分布設(shè)計出相應(yīng)的測量矩陣進行采樣,將相關(guān)系數(shù)低于0.7的槍信號進行缺失如圖7(b)所示。利用本文方法對缺失信號進行重構(gòu)。

圖7 構(gòu)建測量矩陣對隨機干擾信號進行缺失 Fig.7 Missing random interference signal by constructing a measurement matrix

重建結(jié)果如圖8(a)所示。重建后的信號同相軸連續(xù)性相同、振幅一致性準確,且能量分布基本一致。對重構(gòu)前后信號進行誤差分析,如圖8(b)所示。未缺失的信號重建前后的誤差很小,從能量分布上來看基本維持在零值附近,說明重建前后對未缺失的信號影響很小。受隨機干擾影響而進行缺失重建的信號,基本上得到完整的重建,重構(gòu)效果較好。

4 結(jié)論與討論

本研究基于傅里葉變換域構(gòu)建了一種基于壓縮感知理論的氣槍源信號重建方法,解決部分有效氣槍源信號因受隨機干擾信號的影響,而需進行缺失重建的問題。進行數(shù)值模擬與實際資料的處理。并與傳統(tǒng)的三次樣條插值方法進行了對比。對重建效果從信噪比、均方根誤差、互相關(guān)系數(shù)三個方面進行了對比分析。結(jié)果表明:本文方法重構(gòu)前后的一維信號形態(tài)振幅一致性更強,信噪比更高、誤差小、互相關(guān)系數(shù)高。二維信號缺失重建問題轉(zhuǎn)換到一維中處理,變成從時間域采樣到空間域采樣,本文的方法重構(gòu)前后的信號同相軸清晰,連續(xù)性相同,振幅一致性準確,對噪聲壓制較好,隨機干擾得到缺失,并較好的重建出滿足精度要求的有效信號。而三次樣條插值方法重構(gòu)前后信噪比較低、誤差較大,對噪聲壓制較差,重構(gòu)精度較低。對實際資料處理結(jié)果表明:受隨機干擾信號影響的有效信號得到很好的重構(gòu),重構(gòu)后的信號具有同相軸的連續(xù)性強,空間振幅的一致性準確,能量分布一致,重構(gòu)精度高。

圖8 重建結(jié)果Fig.8 Reconstruction results

對于隨機缺失的地震信號,傳統(tǒng)的插值方法方法大都受 Nyquist 采樣理論的限制,插值效果不佳,本文方法對隨機欠采樣缺失的氣槍信號重構(gòu)效果具有明顯的優(yōu)勢。

在重構(gòu)的過程中隨機噪聲對重構(gòu)精度有一定影響,針對信噪比低的信號,也可以通過與稀疏域去噪聲相結(jié)合,可以更好地去除噪聲和提高重構(gòu)精度。

致謝:感謝陳颙院士工作站2018-ET11小隊各位導(dǎo)師在大理培訓(xùn)期間對本文的指導(dǎo)！