秦　晅,宋維琪

(中國石油大學(xué)(華東)地球科學(xué)與技術(shù)學(xué)院,山東青島266580)

?

基于同步壓縮變換微地震弱信號提取方法研究

秦晅,宋維琪

(中國石油大學(xué)(華東)地球科學(xué)與技術(shù)學(xué)院,山東青島266580)

摘要:針對微地震資料的信噪比低,無法清晰識別P波和S波的問題,根據(jù)微地震信號具有隨機性、非平穩(wěn)性的特點,研究了基于同步壓縮變換(synchrosqueezing transform,SST)微地震弱信號提取方法。首先利用SST對信號進行自適應(yīng)閾值去噪,然后在有效信號的頻率中心附近進行SST系數(shù)的積分抽取,再利用抽取的有效信號進行SST重構(gòu)實現(xiàn)弱信號的提取。應(yīng)用于合成的含不同強度噪聲的非平穩(wěn)信號模型以及實際微地震單道記錄的處理結(jié)果表明,該方法具有較好的抗噪能力和較高的信號提取精度。將該方法應(yīng)用于實際井中微地震數(shù)據(jù)的試處理和分析,并與常規(guī)低通濾波結(jié)果進行了對比,表明該方法能夠較好地將弱有效信號從噪聲中提取出來,具有較好的實用價值。

關(guān)鍵詞：微地震;同步壓縮變換;弱信號提取;自適應(yīng)閾值去噪

在油氣田開發(fā)過程中,微地震監(jiān)測是獲得水力壓裂裂縫分布的一種較為有效的方法。由于微地震波能量微弱,微地震資料信噪比較低,有時噪聲甚至完全淹沒了有用的微地震信號,從而無法清晰地識別出P波和S波,影響了后續(xù)的初至拾取、震源定位及裂縫解釋工作[1]。因此,需要對微地震信號進行濾波去噪,將有效信號從噪聲中分離出來。

常用的微地震弱信號提取方法中,傳統(tǒng)的傅里葉變換雖然能提高信噪比,但是對于具有非平穩(wěn)信號特點的微地震信號來說效果不佳[2-3];連續(xù)小波變換(continuous wavelet transform,CWT)在低頻處具有較高的頻率分辨率,在高頻處具有較高的時間分辨率,能夠在時頻域很好地區(qū)分非平穩(wěn)信號的突變部分,但結(jié)果受Heisenberg不確定性原理和小波基函數(shù)的影響[4-5];經(jīng)驗?zāi)B(tài)分解(empirical mode decomposition,EMD)是一種信號自適應(yīng)分解方法,不依賴于基函數(shù)的選取,但這是一種基于經(jīng)驗的方法，并沒有嚴格的數(shù)學(xué)物理基礎(chǔ),同時其穩(wěn)定性受模態(tài)分量分解過程的影響,并且分解時常出現(xiàn)模態(tài)分量混疊現(xiàn)象,造成結(jié)果不準確[6-8];偏振特性分析方法利用有效信號與隨機信號偏振度不同的特性來區(qū)分微地震信號[9],但其本身沒有一定的尺度,無法單獨進行有效信號的檢測;獨立分量分析(independent component analysis,ICA)是基于高階統(tǒng)計特性的用于含噪信號的盲源分離信號處理方法[10-11],但是方法本身需要估計原信號的近似值,應(yīng)用于微地震數(shù)據(jù)資料處理時,受低信噪比的影響,處理結(jié)果不準確[12]。

基于以上分析,本文針對微地震信號隨機性、非平穩(wěn)性的特點,研究了同步壓縮變換(synchrosqueezing transform,SST)微地震弱信號提取方法。SST屬于一種時頻重排算法,不同于原始重排算法的是，其在提高時頻分辨率的同時,也支持信號的重構(gòu)[13-15]。本文首先介紹了SST的基本原理,然后利用合成的不同噪聲強度非平穩(wěn)信號模型和實際微地震資料進行驗證,說明了該方法能夠很好地實現(xiàn)弱信號提取,將微地震的有效信號從噪聲中分離出來。

1基本原理

1.1同步壓縮變換基本理論

一個微地震信號s(t)可以表示成一系列不同頻率的諧波之和[14],即：

(1)

式中：Ak(t)表示瞬時振幅;θk(t)表示第k個諧波成分的瞬時相位;η(t)表示噪聲或干擾。瞬時頻率fk(t)可由瞬時相位的導(dǎo)數(shù)求得：

(2)

首先對微地震信號s(t)進行小波變換：

(3)

式中：Ws(a,b)表示小波系數(shù)譜;ψ*表示母小波函數(shù)的共軛復(fù)數(shù);a表示尺度因子;b表示時間平移因子。

雖然實際得到的小波系數(shù)譜Ws(a,b)經(jīng)常發(fā)生能量擴散,不能很好地聚焦,降低了脊線的可讀性,使得時頻譜變得模糊,但其相位不受尺度變化的影響[16]。因此可以利用小波系數(shù)Ws(a,b)得到的相位來計算其瞬時頻率ωs(a,b)：

(4)

利用計算得到的瞬時頻率,可以建立(a,b)→(ωs(a,b),b)之間的映射關(guān)系,同步壓縮變換再對時間—尺度平面的能量進行重新分配,將其轉(zhuǎn)化為時間—頻率平面,這就是小波變換與重排瞬時頻率結(jié)合的壓縮變換的基本思想。離散情況下,尺度坐標表示為Δak=ak-ak-1,頻率坐標表示為Δω=ωl-ωl-1,那么壓縮變換Ts(ω,b)確定在ωl的中心位置,其中ωl的范圍為[ωl-Δω/2,ωl+Δω/2]。因此壓縮變換的離散公式如下：

(5)

由于同步壓縮變換是對小波變換的復(fù)數(shù)譜僅沿著頻率軸方向進行重排,因此是可逆的,其重構(gòu)的信號s(t)可表示為：

(6)

其中,Cφ為依賴選取母小波的常數(shù),可表示為：

(7)

式中：Ψ*表示母小波的傅里葉變換;ξ表示母小波的主頻。

1.2同步壓縮變換參數(shù)選取

median[Ws(a1:nv,b)]|}/0.6745

(8)

式中：Ws(a1:nv,b)表示尺度因子長度上小波系數(shù);0.6745表示歸一化高斯分布標準偏差系數(shù),利用計算的ση就可以自適應(yīng)估算最優(yōu)的閾值γ：

(9)

式中：n表示微地震信號的采樣點數(shù)。

1.3小波閾值去噪

小波閾值去噪方法包括硬閾值和軟閾值方法,硬閾值函數(shù)表達式為：

(10)

軟閾值函數(shù)表達式為：

(11)

由于硬閾值能夠保留較大的小波系數(shù),也可以保持信號有效部分的特征;而軟閾值改變了原信號,不能完全保持信號的特征[7],所以本文采用了硬閾值方法進行小波閾值去噪。

1.4同步壓縮變換有效信號提取

利用同步壓縮變換提取有效信號的分量主要分為三步：首先根據(jù)時頻分布的脊線確定每個有效信號分量頻帶范圍lk(tm)；然后在有效信號頻帶范圍內(nèi)對Ts(ωl,b)系數(shù)進行積分抽?。蛔詈罄贸槿〉挠行盘朤s(ωl,b)系數(shù)進行同步壓縮變換重構(gòu),實現(xiàn)sk(tm)有效信號提取。可將(6)式改寫為：

(12)

式中：sk(tm)表示第k個提取的有效信號;Cφ表示依賴選取母小波的常數(shù);Re表示對信號取實部;lk(tm)表示有效信號頻帶范圍;Ts(ωl,b)表示同步壓縮變換系數(shù)。

綜上所述,基于同步壓縮變換微地震弱信號提取方法的實現(xiàn)步驟為：

1) 先對微地震信號s(t)進行同步壓縮變換,計算得到系數(shù)Ts(ωl,b);

2) 再利用公式(9)計算得到的γ進行自適應(yīng)閾值去噪;

3) 在有效信號頻帶范圍內(nèi)對同步壓縮變換Ts(ωl,b)系數(shù)進行積分抽取;

4) 利用抽取的有效信號Ts(ωl,b)系數(shù)進行同步壓縮變換重構(gòu),實現(xiàn)微地震弱信號提取。

2仿真分析

2.1模擬非平穩(wěn)信號分析

為了驗證同步壓縮變換弱信號提取的可行性,本文模擬了一個非平穩(wěn)信號f(t)進行分析。該信號采樣間隔為5ms,采樣點數(shù)為2000,信號由3組非平穩(wěn)信號組成,具體如下：

(13)

x1(t)=[1+0.5cos(t)]cos(4πt)

(14)

x2(t)=2e-0.1tcos{2π[3t+0.25sin(1.4t)]}

(15)

x3(t)=[1+0.5cos(2.5t)]cos[2π(5t+2t1.3)]

(16)

為了驗證不同噪聲水平對SST弱信號提取的影響,定義了原始信號與加入噪聲的能量比值為信噪比(SNR),計算公式為：

(17)

圖1a和圖1b分別為模擬的非平穩(wěn)原信號與信噪比為0.75的含噪混合信號的波形圖。利用SST方法對含噪信號進行提取,結(jié)果見圖2。圖2a和圖2b 分別為含噪信號SST結(jié)果和自適應(yīng)閾值去噪后SST結(jié)果,可以看出,SST對于時頻譜的脊線刻畫準確,提高了時頻譜的分辨率,并且自適應(yīng)閾值去噪也能很好地壓制隨機噪聲。從圖2c到圖2f可以看出,SST較好地重構(gòu)與提取了信號,與原信號相比,提取的信號除了邊界以外,其余部分吻合較好。

圖1　模擬的原信號(a)和含噪混合信號(b)

圖2　信噪比為0.75時的含噪信號SST提取結(jié)果a 含噪混合信號SST結(jié)果; b 自適應(yīng)閾值去噪后SST結(jié)果; c 自適應(yīng)閾值去噪后SST重構(gòu)信號與原信號的對比; d 自適應(yīng)閾值去噪后SST提取信號x1(t)與原信號x1(t)的對比; e 自適應(yīng)閾值去噪后SST提取信號x2(t)與原信號x2(t)的對比; f 自適應(yīng)閾值去噪后SST提取信號x3(t)與原信號x3(t)的對比

不斷改變含噪信號的信噪比,觀察SST方法提取的效果,并計算重構(gòu)與提取的信號與原信號的相關(guān)系數(shù),結(jié)果如表1所示。由表1可以看出,隨著噪聲的增強,SST提取的信號與原信號的相關(guān)系數(shù)雖然呈現(xiàn)下降趨勢,但是相關(guān)系數(shù)都在0.92以上,說明了SST具有較好的抗噪能力,在一定的噪聲強度內(nèi),仍然具有較高的信號提取精度。

2.2實際微地震單道記錄分析

為了驗證SST對實際微地震數(shù)據(jù)提取效果,對圖3a所示的實際微地震單道記錄進行了分析,該信號的采樣間隔為0.5ms,抽取采樣點數(shù)為600。圖3b和圖3c分別為短時傅里葉變換(SIFT)結(jié)果與同步壓縮變換(SST)結(jié)果,可以看出,SIFT結(jié)果的時頻譜模糊,不能很好地分辨出P波和S波,而SST結(jié)果分辨率明顯提高了很多,能夠很好地分辨出P波和S波。圖3d和圖3e分別為SST提取的P波和S波與原P波和S波的對比結(jié)果,可以看出,SST提取的P波和S波與原P波和S波的相關(guān)系數(shù)分別為0.9925和0.9984,說明了SST對實際微地震數(shù)據(jù)信號提取的有效性。

表1　不同噪聲水平下SST提取信號與原信號的相關(guān)系數(shù)對比

在圖3a所示的實際微地震數(shù)據(jù)中加入噪聲,并計算SST提取的P波和S波與原P波和S波的相關(guān)系數(shù),結(jié)果如表2所示。由表2可以看出,隨著噪聲的增強,P波與S波的相關(guān)系數(shù)雖然呈現(xiàn)下降的趨勢,但是相關(guān)系數(shù)都能保持在0.90以上,說明了SST具有較好的抗噪能力,在一定的噪聲強度內(nèi),仍然具有較高的信號提取精度。從P波與S波相關(guān)系數(shù)對比來看,S波提取的精度明顯比P波高,這主要是由于P波的信號能量相對較弱,當加入較強噪聲后,信號會被噪聲淹沒,造成了提取精度的降低。

圖3　實際微地震單道記錄SST提取結(jié)果a 實際微地震單道記錄; b SIFT結(jié)果; c SST結(jié)果; d SST提取P波與原P波對比; e SST提取S波與原S波對比

信號信噪比SST提取信號相關(guān)系數(shù)P波S波0.50.90350.98101.50.92440.98364.00.93490.991410.00.95600.9968

3應(yīng)用實例分析

本文截取了某油田一口壓裂井的一段井中微地震數(shù)據(jù)進行分析,原始數(shù)據(jù)如圖4a所示,該數(shù)據(jù)采用16個三分量檢波器接收,采樣間隔為0.25ms,采樣點數(shù)為1000。圖4a中1～16,17～32,33～48道分別為數(shù)據(jù)的Z,X,Y分量,由于數(shù)據(jù)的信噪比較低,噪聲幾乎淹沒了有用的微地震信號,無法清晰地識別出P波和S波。圖4b為對原始數(shù)據(jù)進行低通濾波處理后的結(jié)果,可以看出,雖然高頻噪聲被濾除,但是微地震的有效信號同相軸仍然不夠清晰。圖4c為利用SST提取弱信號的結(jié)果,可以看出,微地震的噪聲得到很好的壓制,并且較好地將有效信號從噪聲中提取出來。

圖4　實際微地震數(shù)據(jù)弱信號提取結(jié)果對比a 原始數(shù)據(jù); b 低通濾波結(jié)果; c SST提取結(jié)果

4結(jié)論

本文針對微地震資料的信噪比低,無法清晰地識別出P波和S波的問題,根據(jù)同步壓縮變換在提高信號時頻分辨率的同時,也能夠支持信號重構(gòu)的特點,研究了同步壓縮變換微地震弱信號提取方法。合成的不同噪聲強度非平穩(wěn)信號模型和實際微地震單道記錄的測試表明本文方法具有較好的抗噪能力與較高的信號提取精度。利用本文方法對實際井中微地震數(shù)據(jù)進行了分析處理,并與常規(guī)的低通濾波處理結(jié)果進行對比,說明本文方法具有較好的實用價值,能夠?qū)⒂行盘枏脑肼曋刑崛〕鰜怼Ｓ捎谕綁嚎s變換具有高分辨率時頻分布的特點,對后續(xù)的初至拾取、震源定位及裂縫解釋的研究也具有良好的應(yīng)用前景。

參考文獻

[1]宋維琪,何可,郭全仕,等.微地震資料自適應(yīng)濾波方法研究[J].石油物探,2013,52(3):229-233

SONG W Q,HE K,GUO Q S,et al.Micro-seismic data adaptive filtering method[J].Geophysical Prospecting for Petroleum,2013,52(3):229-233

[2]ALVANITOPOULOS P F,PAPAVASILEIOU M,ANDREADIS I,et al.Seismic intensity feature construction based on the Hilbert-Huang transform[J].IEEE Transactions on Instrumentation and Measurement,2012,61(2):326-337

[3]GACI S.The use of wavelet-based denoising techniques to enhance the first-arrival picking on seismictraces[J].IEEE Transactions on Geoscience and Remote Sensing,2014,52(8):4558-4563

[4]SINHA S,POUTH P S,ANNO P D,et al.Spectral decomposition of seismic data with continuous-wavelet transforms[J].Geophysics,2005,70(6):19-25

[5]BEENAMOL M,PRABAVATHY S,MOHANALIN J.Wavelet based seismic signal de-noising using Shannon and Tsallis entropy[J].Computers & Mathematics with Applications,2012,64(11):3580-3593

[6]BATTISTA B M,KNAPP C,MCGEE T,et al.Application of the empirical mode decomposition and Hilbert Huang transform to seismic reflection data[J].Geophysics,2007,72(2):31-45

[7]王姣,李振春,王德營.基于CEEMD的地震數(shù)據(jù)小波閾值去噪方法研究[J].石油物探,2014,53(2):164-172

WANG J,LI Z C,WANG D Y.A method for wavelet threshold denoising of seismic data based on CEEMD[J].Geophysical Prospecting for Petroleum,2014,53(2):164-172

[8]李月,彭蛟龍,馬海濤,等.過渡內(nèi)蘊模態(tài)函數(shù)對經(jīng)驗?zāi)B(tài)分解去噪結(jié)果的影響研究及改進算法[J].地球物理學(xué)報,2013,56(2):626-634

LI Y,PENG J L,MA H T,et al.Study of the influence of transition IMF on EMD do-noising and the improved algorithm[J].Chinese Journal Geophysics,2013,56(2):626-634

[9]吳治濤,駱循,李仕雄.聯(lián)合小波變換與偏振分析自動拾取微地震P波到時[J].地球物理學(xué)進展,2012,7(1):131-0136

WU Z T,LUO X,LI S X.United wavelet transform and polarization analysis automatically identify micro seismic P-arrival[J].Progress in Geophysics,2012,27(1):131-136

[10]王維強,楊國權(quán).基于EMD和ICA的地震信號去噪技術(shù)研究[J].石油物探,2012,51(1):19-29

WANG W Q,YANG G Q.Random noise denoising technique based on EMD and ICA[J].Geophysical Prospecting for Petroleum,2012,51(1):19-29

[11]VAN DER BAAN M.PP/PS Wavefield separation by independent component analysis[J].Geophysical Journal of International,2006,166(1):339-348

[12]賈瑞生,趙同彬,孫紅梅,等.基于經(jīng)驗?zāi)B(tài)分解及獨立成分分析的微震信號降噪方法[J].地球物理學(xué)報,2015,58(3):1013-1023

JIA R S,ZHAO T B,SUN H M,et al.Micro-seismic signal denoising method based on empirical mode decomposition and independent component analysis[J].Chinese Journal Geophysics,2015,58(3):1013-1023

[14]HERRERA R H,TARY J B,VAN DER BAAN M,et al.Body wave separation in the time-frequency domain[J].IEEE Geoscience and Remote Sensing Letters,2015,12(2):364-368

[15]陳小旺,馮志鵬,LIANG Ming.基于迭代廣義同步壓縮變換的時變工況行星齒輪箱故障診斷[J].機械工程學(xué)報,2015,51(1):131-137

CHEN X W,FENG Z P,LIANG M.Planetary gearbox fault diagnosis under time-variant conditions based on iterative generalized synchrosqueezing transform[J].Journal of Mechanical Engineering,2015,51(1):131-137

[16]尚帥,韓立國,胡瑋,等.壓縮小波變換地震譜分解方法應(yīng)用研究[J].石油物探,2015,54(1):51-55

SHANG S,HAN L G,HU W,et al.Applied research of synchrosqueezing wavelet transform in seismic spectral decomposition method[J].Geophysical Prospecting for Petroleum,2015,54(1):51-55

[17]HERRERA R H,HAN J,VAN DER BAAN M.Applications of the synchrosqueezing transform in seismic time-frequency analysis[J].Geophysics,2014,79(3):55-64

[18]楊心超,朱海波,崔樹果,等.P波初動震源機制解在水力壓裂微地震監(jiān)測中的應(yīng)用[J].石油物探,2015,54(1):43-50

YANG X C,ZHU H B,CUI S G,et al.Application of P-wave first-motion focal mechanism solutions in microseismic monitoring for hydraulic fracturing[J].Geophysical Prospecting for Petroleum,2015,54(1):43-50

[19]王小杰,欒錫武.基于小波分頻技術(shù)的地層Q值提取方法研究[J].石油物探,2015,54(3):260-266

WANG X J,LUAN X W.Q value extraction method based on wavelet frequency division technology[J].Geophysical Prospecting for Petroleum,2015,54(3):260-266

[20]汪祥莉,王斌,王文波,等.混沌干擾中基于同步擠壓小波變換的諧波信號提取方法[J].物理學(xué)報,2015,64(10):100201

WANG X L,WANG B,WANG W B,et al.Harmonic signal extraction from chaotic interference based on synchrosqueezed wavelet transform[J].Acta Physica Sinica,2015,64(10):100201

[21]盛冠群,李振春,王維波,等.基于小波分解與高階統(tǒng)計量的微地震初至拾取方法研究[J].石油物探,2015,54(4):388-395

SHENG G Q,LI Z C,WANG W B,et al.A new automatic detection method of microseismic events based on wavelet decomposition and high-order statistics[J].Geophysical Prospecting for Petroleum,2015,54(4):388-395

(編輯：朱文杰)

Weak signal extraction method of microseismic data based on synchrosqueezing transform

QIN Xuan,SONG Weiqi

(SchoolofGeosciences,ChinaUniversityofPetroleum,Qingdao266580,China)

Abstract:It is very difficult to clearly identify the P-wave and S-wave due to the low SNR of microseismic data.According to the characteristics of randomness and non-stationary of microseismic signals,we proposed a method to extract the weak signal from microseismic data based on synchrosqueezing transform (SST).Firstly,we utilize the SST to conduct the adaptive threshold denoising.Then,the SST coefficients are extracted near the center frequency of effective signal by integrating.Finally,we carry out the SST reconstruction using the extracted effective signal SST coefficients to implement the weak signal extraction.The test results of the synthesis non-stationary signal models with different noise intensity and the actual microseismic single-channel record show that this method has better noise immunity and higher signal extraction accuracy.Through the processing and analysis of the actual borehole microseismic data and the comparison with conventional low-pass filtering,the proposed method has good capability and is better to extract the weak effective signal from the noise.

Keywords:microseismic,synchrosqueezing transform,weak signal extraction,adaptive threshold denoising

文章編號：1000-1441(2016)01-0060-07

DOI:10.3969/j.issn.1000-1441.2016.01.008

中圖分類號：P631

文獻標識碼：A

基金項目：國家高技術(shù)研究發(fā)展計劃(863計劃)(2011AA060303)項目資助。

作者簡介：秦晅(1988—),男,碩士在讀,主要從事微地震信號處理方面研究工作。

收稿日期：2015-06-19;改回日期：2015-11-19。

秦晅,宋維琪.基于同步壓縮變換微地震弱信號提取方法研究[J].石油物探,2016,55(1):-66

QIN Xuan,SONG Weiqi.Weak signal extraction method of microseismic data based on synchrosqueezing transform[J].Geophysical Prospecting for Petroleum,2016,55(1):-66

This research is financially supported by the National High-tech R&D Program (863 Program) (Grant No.2011AA060303).