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

張 帥, 楊潤海, 高爾根

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

0 引言

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

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

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

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

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

1 壓縮感知?dú)鈽屝盘?hào)重建方法(CS)

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

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

1.1 方法原理

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

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

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

(1)

(2)

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

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

Θ=ΨTy

(3)

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

x(n)=φΘ=φΨTy

(4)

令α=ΨTy,因此得到:

x(n)=φα

(5)

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

(6)

輸入:測量矩陣φ,測量信號(hào)向量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 信號(hào)重建處理流程

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

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

2 數(shù)值模擬

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

(12)

(13)

2.1 一維數(shù)值模擬

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

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

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

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

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

2.2 二維數(shù)值模擬

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

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

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

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

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

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

3 實(shí)際資料處理

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

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

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

4 結(jié)論與討論

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

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

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

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

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