廣義S變換在背景噪聲數(shù)據(jù)處理中的應(yīng)用①
李紅光1, 孫剛1, 邢藝蘭2
(1.中國地震應(yīng)急搜救中心,北京 100049; 2.北華航天工業(yè)學(xué)院,河北 廊坊 065000)
摘要:廣義S變換解決了標(biāo)準(zhǔn)S變換中基本小波形態(tài)固定的缺陷,在分析非平穩(wěn)信號時能更好地刻畫其時頻特性。編寫廣義S變換的計算程序,并用它對地震背景噪聲數(shù)據(jù)進(jìn)行去噪處理。結(jié)果顯示:廣義S變換處理過的數(shù)據(jù)信噪比大大提高,可以為數(shù)據(jù)反演提供更精確的數(shù)據(jù)。
關(guān)鍵詞:廣義S變換; 背景噪聲; 信噪比; 提高
收稿日期:①2014-12-25
基金項目:國家重大科學(xué)儀器設(shè)備開發(fā)專項(2012YQ16018506)
作者簡介:李紅光,男,高級工程師,主要從事地震工程、地震信號處理等工作。E-mail:fengzhe122@126.com。
中圖分類號:P631文獻(xiàn)標(biāo)志碼:A
DOI:10.3969/j.issn.1000-0844.2015.03.0867
Application of Generalized S-transform in Ambient
Seismic Noise Data Processing
LI Hong-guang1, SUN Gang1, XING Yi-lan2
(1.NationalEarthquakeResponseSupportService,Beijing100049,China;
2.NorthChinaInstituteofAerospaceEngineering,Langfang065000,Hebei,China)
Abstract:In this study, we used the generalized S-transform method to improve the signal-to-noise ratio of ambient seismic noise data, and improved the efficiency of geophysical inversion. The S-transform method is a time-frequency representation of local spectral phase properties. A key feature of the S-transform method is that it can uniquely combine a frequency dependent resolution of the time-frequency space and absolutely reference local phase information. As such, the phase in a local spectrum setting can be defined, and this makes possible the production of many desirable wave characteristics. The scaling property of the Gaussian window is reminiscent of the scaling property of continuous wavelets, because one Fourier frequency wavelength is always equal to one standard deviation of the window. In this study, we introduce a variant of the original S-transform, which replaces f with λfp. In this variant, one standard deviation of the Gaussian window contains λfp wavelengths of the Fourier sinusoid at all frequencies. The generalized S-transform solves the defect in the analysis of non-stationary signals caused by the fixed form of the basic wavelet in the S-transform, and the time-frequency characteristics can better describe the non-stationary signals. We wrote a computer program for the generalized S-transform, and for ambient seismic noise we used data denoising. The results showed that the signal-to-noise ratio is greatly improved through the generalized S-transform denoising, and yields more accurate data for the inversion process. There have been few studies on the denoised processing of background noise data using the generalized S-transform, and the resolution of our processing results is better than that achieved by the ordinary S-transform.
Key words: generalized S-transform; ambient noise; signal-to-noise ratio; improvement
0引言
在分析非平穩(wěn)信號時,傳統(tǒng)的傅里葉變換只能將信號從時域?qū)?yīng)到一維頻域,不能有效地反映非平穩(wěn)信號的頻率隨時間的變化,因此很難分析信號的時頻變化特性。時頻分析方法是將一維的時域信號變換到二維時頻平面上,然后在時頻域平面上區(qū)分并提取信號分量。常用的時頻分析方法有短時傅里葉變換、Wigner-Vile分布和小波變換等[1]。
1996年美國地球物理學(xué)家Stockwell等[2]提出一種加時窗傅里葉變換方法(S變換),它吸收并發(fā)展了短時傅里葉變換和小波變換[3]。信號S變換的分辨率與頻率有關(guān),且其結(jié)果具有無損可逆性[4]。但是S變換中的基本小波函數(shù)是固定的,這使其在實際應(yīng)用中受到了限制。高靜懷等[5]提出根據(jù)“實際需要”恰當(dāng)?shù)剡x擇或構(gòu)造基本小波函數(shù)的廣義S變換方法,并將其應(yīng)用在地震數(shù)據(jù)的時頻分析工作中。Pinnegar[6]給出了既可調(diào)節(jié)窗口標(biāo)準(zhǔn)差又非對稱的另外一種廣義S變換方法。Pinnegar[7]在研究信號頻率隨時間變化的問題時,提出了多尺度和非零相位復(fù)數(shù)窗的廣義S變換方法。
實際應(yīng)用中S變換是先把信號從時間域轉(zhuǎn)換到時頻域,然后再轉(zhuǎn)換回頻率域,最后轉(zhuǎn)回時間域,因此具有無損可逆性,在實際應(yīng)用中具有非常重要的意義。被國內(nèi)外廣泛應(yīng)用于地震勘探[8-9]、地震震相和主頻分析[10-11]、人工地震探測[3,12-17]、地震波的衰減分析[18]、地震層特征分析[19]和油氣探測與預(yù)測[20-22]等領(lǐng)域,尤其在時變?yōu)V波中去除噪聲具有明顯的效果。
地震臺站記錄一直把噪聲當(dāng)作干擾,在地震數(shù)據(jù)處理中會用各種方法將其去除。近年來理論研究表明:通過對漫射波場(背景噪聲、散射尾波)做互相關(guān)計算可以提取出臺站間的格林函數(shù)[23-24],這促使背景噪聲互相關(guān)技術(shù)開始被研究者從多個角度進(jìn)行研究[25-27]。目前該技術(shù)的數(shù)據(jù)處理方法一般是對每天提取出的地震背景噪聲互相關(guān)函數(shù)(簡稱NCF)直接疊加,其信噪比較低[28]。為了得到比較穩(wěn)定的NCF數(shù)據(jù),通常需要幾個月甚至更長時間的觀測數(shù)據(jù)。因此應(yīng)用較好的時頻分析方法對NCF數(shù)據(jù)進(jìn)行降噪處理,可進(jìn)一步提高NCF的信噪比,同時可進(jìn)一步減少其疊加次數(shù)。本文給出了一種廣義S變換,對從首都圈地區(qū)地震臺站的觀測資料中得到的NCF進(jìn)行降噪處理,以提高其信噪比。
1原理方法
1.1S變換與反變換
函數(shù)h(t)的時間域S變換可以表示為:
其中:h(t)是信號時間函數(shù);τ是時頻平面內(nèi)的時間;f是頻率。
根據(jù)傅里葉變換與卷積定理得到頻率域S變換公式:
式(2)中參數(shù)意義同上,H(f+α)為信號h(t)的傅里葉譜。
S變換的良好性質(zhì)之一在于它是完全可逆的。式(3)為頻率反S變換:
1.2廣義S變換
S變換中的小波形態(tài)是固定的,這使其在分析不同特征的信號時受到了一定限制。本文給出一種廣義S變換的表達(dá)形式:
其中:λ>0,p>0。
當(dāng)λ=1,p=1時,上式即為標(biāo)準(zhǔn)S變換。
廣義S變換的反變換與式(3)為同一形式。廣義S變換通過引入兩個參數(shù)λ和p改造標(biāo)準(zhǔn)S變換的高斯窗函數(shù),在實際應(yīng)用中,它能根據(jù)信號的頻率分布特點和時頻分析的目的靈活地調(diào)節(jié)高斯窗函數(shù)隨頻率的變化趨勢。不僅可以根據(jù)需要調(diào)整時窗寬度隨信號頻率變化的速度,還能使時頻平面上時窗峰值點的連線呈現(xiàn)多種變化特征,使廣義S變換能夠更好地應(yīng)用到具體信號的分析和處理中。
圖1給出了不同λ、p值對應(yīng)的高斯窗函數(shù)。從圖中明顯可以看出:λ和p的增大都使窗函數(shù)變寬,而減小會使窗函數(shù)變窄。根據(jù)這個特征,可通過試驗不同的λ、p值對信號進(jìn)行廣義S變換,得出最佳的參數(shù)。
2背景噪聲數(shù)據(jù)處理
2.1去噪步驟
采用廣義S變換去噪的具體步驟為:
圖1 不同的λ、p值對高斯窗函數(shù)的影響 Fig.1 Influence of different value of λ and p on the Gauss window
(1) 對背景噪聲原始數(shù)據(jù)進(jìn)行廣義S變換,得到時頻譜。
(2) 對時頻譜進(jìn)行閥值濾波。閥值濾波是指先尋找有效信號時頻窗內(nèi)的最大譜值A(chǔ)max,然后根據(jù)信號本身時頻特征選擇合適的閥值比例系數(shù)λ,得到截斷閥值A(chǔ)c=λ×Amax,將時頻窗內(nèi)所有小于截斷閥值A(chǔ)c的譜值都設(shè)為0。
(3) 對濾波后的譜值數(shù)據(jù)進(jìn)行反變換,得到去噪后的數(shù)據(jù)。
2.2實例驗證
選取一條由中國地震局地球物理研究所流動觀測臺提供的背景噪聲數(shù)據(jù)記錄。該數(shù)據(jù)記錄的是由華北某地兩臺站間的背景噪聲提取的格林函數(shù)(圖2),并對此記錄進(jìn)行廣義S變換去噪處理。
圖2 華北某兩個臺間提取的NCF Fig.2 NCF obtained from a double stations in North China
圖3 原始數(shù)據(jù)的群速度頻散曲線 Fig.3 Group velocity dispersion curves of the original data
利用CPS330(ComputerProgramsforSeismologyVersion3.30)軟件包中的do_mft子程序提取面波群速度頻散信息。該程序是基于Dziewonski提出的多重濾波法(MultipleFilterTechnique)編制的,提供了人機(jī)交互式的界面,得到的頻散曲線結(jié)果較為可靠。原始數(shù)據(jù)的群速度頻散曲線如圖3所示。
從圖3可以看出,群速度頻散曲線的速度在2~5km/s之間;周期在6~30s的部分比較連續(xù)、平滑,周期大于30s時能量比較低,提取的頻散不可信,周期小于6s時能量比較分散,使得頻散無法追蹤。由于有效信號的速度在2~5km/s之間,根據(jù)兩個臺站間的距離h可以得到有效信號的時間窗,對時間窗外的信號歸零(圖4)。
圖4 根據(jù)有效速度進(jìn)行預(yù)處理后的信號 Fig.4 Signal pretreated by the effective speed
對圖4中的信號進(jìn)行廣義S變換(取λ=2.0,p=1.2),得到它的廣義S變換頻譜[圖5(a)]。
圖5(a)顯示此信號的能量比較分散,信噪比不高,選取不同的閾值對比信號記錄進(jìn)行濾波。取不同的閥值比例系數(shù)進(jìn)行試驗,經(jīng)過對比最終選擇了去噪效果比較理想的閥值比例系數(shù)λ=0.15。從圖5(b)看出,原始信號記錄的廣義S變換頻譜中,有效信號區(qū)域周圍的噪聲被所設(shè)定的閾值濾去,從而得到了信噪比較高的信號。圖5(c)也顯示,有效信號時間窗內(nèi)的波形變得更加平滑和連續(xù),時間窗外的微小抖動是由吉布斯現(xiàn)象引起的。
仍然采用CPS330軟件包中的do_mft程序?qū)θピ牒蟮男盘柼崛∶娌ㄈ核俣阮l散信息(圖6)。
從圖6(b)可以看出:原始信號經(jīng)廣義S變換去噪后,群速度能夠連續(xù)追蹤至4s,長周期部分能量也比去噪前更集中了,但頻散曲線在25s左右分成了兩支,這跟閥值濾波可能有關(guān)系。圖6(a)中左圖傅里葉功率譜在去噪后,短周期部分比去噪前連續(xù)了很多;雖然長周期部分能量明顯降低,但功率譜變得連續(xù)并且能量集中。
圖5 有效信號的S譜和去噪波形 Fig.5 S spectrum and denoising waveform of the effective signal
圖6 去噪后的信號群速度頻散曲線 Fig.6 Group velocity dispersion curves of the denoising signal
3結(jié)論
給出一種用廣義S變換進(jìn)行閾值濾波的方法。實例去噪的結(jié)果表明:這種去噪方法不僅更直觀地得出信號時頻特性,而且濾波效果也優(yōu)于傳統(tǒng)的一維濾波;對背景噪聲數(shù)據(jù)進(jìn)行處理后得到了清晰平滑的頻散曲線,對信號的短周期頻散測量尤其有效。
本文得出下幾點結(jié)論:
(1) 廣義S變換是一種很好的時頻二維變換方法,可以將數(shù)據(jù)中的信號和噪聲在時-頻域中分離。
(2) 不同的數(shù)據(jù)應(yīng)該用不同的閥值,合適的閾值可以濾去與有效信號混疊在一起的噪聲。
致謝:本研究得到房立華博士的指導(dǎo),并提供了背景噪聲觀測數(shù)據(jù)。
參考文獻(xiàn)(References)
[1]張賢達(dá),保錚.非平穩(wěn)信號分析與處理[M].北京:國防工業(yè)出版社,2001.
ZHANGXian-da,BAOZheng.NonstationarySignalAnalysisandProcessing[M].Bejing:NationalDefenceIndustryPress,2001.(inChinese)
[2]StockwellRG,MansinhaL,LoweRP,etal.LocalizationoftheComplexSpectrum:TheSTransform[J].IEEETransactionsofSignalProcessing,1996,44(4):998-1001.
[3]趙淑紅,朱光明.S變換時頻濾波去噪方法[J].石油地球物理勘探,2007,42(4):402-406.
ZHAOShu-hong,ZHUGuang-ming.StransformTime-frequencyFilteringDenoisingMethod[J].OilGeophysicalProspecting,2007,42(4):402-406.(inChinese)
[4]趙淑紅,王璇.S變換時頻濾波與其它濾波方法的比較[J].西北地震學(xué)報,2007,29(3):224-229.
ZHAOShu-hong,WANGXuan.ComparisonbetweentheTime-frequencyFilteringofSTransformMethodandOtherFilterMethod[J].NorthwesternSeismologicalJournal,2007,29(3):224-229.(inChinese)
[5]高靜懷,陳文超,李幼銘,等.廣義S變換與薄互層地震響應(yīng)分析[J].地球物理學(xué)報,2003,46(4):526-532.
GAOJing-huai,CHENWen-chao,LIYoumingetal.GeneralizedSTransformandSeismicResponseAnalysisofThinInterbeds[J].ChineseJournalofGeophsics,2003,46(4):526-532.(inChinese)
[6]PinnegarCR,MansinhaL.TheS-transformwithWindowsofArbitraryandVaryingShape[J].Geophysics,2003,68(1):381-385.
[7]PinnegarCR,MansinhaL.Time-localFourierAnalysiswithaScalable,Phase-modulatedAnalyzingFunction:theS-transformwithaComplexWindow[J].2004,84(1):1167-1176.
[8]鄒文,陳愛萍,顧漢明.聯(lián)合時頻分析技術(shù)在地震勘探中的應(yīng)用[J].勘探地球物理進(jìn)展,2004,27(4):246-250.
ZOUWen,CHENAi-ping,GUHan-ming.JointTime-frequencyAnalysisandItsApplicationinSeismicProspecting[J].ProgressinExplorationGeophysics,2004,27(4):246-250.(inChinese)
[9]鄧攻.S變換在地震時頻分析中的對比研究[J].研究與探討,2011(4):73-76.
DENGGong.ComparativeStudyonSChangeinAnalysisofTime-frequencyModulation[J].ResearchandDiscussion,2011(4):73-76.(inChinese)
[10]鄒文,陳愛萍,賀振華,等.基于S變換的地震相分析技術(shù)[J].石油物探,2006,45(1):48-51.
ZOUWen,CHENAi-ping,HEZhen-hua,etal.TheTechniqueofSeismicPhaseAnalysisBasedonStransform[J].GeophysicalProspectingforPetroleum,2006,45(1):48-51.(inChinese)
[11]熊曉軍,賀振華,黃德濟(jì).基于廣義S變換的主頻率分析[J].工程地球物理學(xué)報,2007, 4(6):525-528.
XIONGXiao-jun,HEZhen-hua,HUANGDe-ji.DominatingFrequencyAnalysisBasedonGeneralizedSTransform[J].ChineseJournalofEngineeringGeophysics,2007, 4(6):525-528.(inChinese)
[12]馬見青,李慶春.S變換和TT變換聯(lián)合壓制地震面波[J].物探化探計算技術(shù),2011, 33(3 ):252-257.
MAJian-qing,LIQing-chun.JointSTransformandTTTransformSuppressionofSurfaceWave[J].ComputingTechniquesforGeophysicalandGeochemicalExploration,2011,33(3 ):252-257.(inChinese)
[13]楊海濤,朱仕軍,楊愛國,等.S變換時變?yōu)V波在去噪處理中的應(yīng)用研究[J].西南石油大學(xué)學(xué)報:自然科學(xué)版,2009,31(6):56-58.
YANGHai-tao,ZHUShi-jun,YANGAi-guo,etal.STransformTime-varyingFilteringintheDenoisingApplicationProcess[J].JournalofSouthwestPetroleumUniversity:Science&TechnologyEdition,2009,31(6):56-58.(inChinese)
[14]段俊,張白林,潘樹林.廣義S變換去除面波技術(shù)及應(yīng)用[J].物探化探計算技術(shù),2011, 33(3):280-285.
DUANJun,ZHANGBai-lin,PANShu-lin.TheTechnologyandApplicationtoRemoveSurfaceWavewiththeGeneralizedSTransform[J].ComputingTechniquesofGeophysicalandGeochiemical,2011,33(3):280-285.(inChinese)
[15]陳昕,李新安,劉威序,等.基于S變換的初至拾取方法及應(yīng)用[J].中國煤炭地質(zhì),2009,21:59-63.
CHNEXi,LIXin-an,LIUWei-xu,etal.MethodologyandApplicationofFirstBreakPickingBasedonS-transform[J].CoalGeologyofChina,2009,21:59-63.(inChinese)
[16]李雪英,侯相輝.基于廣義S變換的疊前高頻噪聲壓制[J].石油地球物理勘探,2011,46(4):545-549.
LIXue-ying,HOUXiang-hui.ThePrestackNoiseSuppressionBasedontheGeneralizedStransform[J].OilGeophysicalProspecting,2011,46(4):545-549.(inChinese)
[17]高靜懷,滿蔚仕,陳樹民.廣義S變換域有色噪聲與信號識別方法[J].地球物理學(xué)報,2004,47:869-875.
GAOJing-huai,MANYu-shi,CHENShu-min.RecognitionofSignalsfromColoredNoiseBackgroundinGeneralizedSTransformationDomain[J].ChineseJournalofGeophysics,2004,47:869-875.(inChinese)
[18]李雪英,孫丹,井涌泉,等.利用廣義S變換進(jìn)行等效Q值掃描分析[J].地球物理學(xué)進(jìn)展,2009(5):1696-1702.
LIXue-ying,SUNDan,JINGYong-quan,etal.ScanningAnalysisofEffectiveQUsingGeneralizedSTransform[J].ProgressinGeophysics,2009(5):1696-1702.(inChinese)
[19]劉喜武,劉洪,李幼銘,等.基于廣義S變換研究地震地層特征[J].地球物理學(xué)進(jìn)展,2006(21):440-451.
LIUXi-wu,LIUHong,LIYou-ming,etal.StudyonCharacteristicsofSeismicStratigraphybyGeneralizedS-transform[J].ProgressinGeophysics,2006(21):440-451.(inChinese)
[20]陳學(xué)華,賀振華,黃德濟(jì),等.時頻域油氣儲層低頻陰影檢測[J].地球物理學(xué)報,2009(52):215-220.
CHENXue-hua,HEZheng-hua,HUANGDe-ji,etal.LowFrequencyShadowDetectionofGasReservoirsinTime-frequencyDomain[J].ChinaJournalofGeophysics,2009(52):215-220.(inChinese)
[21]張固瀾,賀振華,張彥斌,等.基于改進(jìn)的廣義S變換的低頻陰影檢測[J].地球物理學(xué)進(jìn)展,2010,25(6):2040-2044.
ZHANGGu-lan,HEZhen-hua,ZHANGYan-bin,etal.DetectionofLow-frequencyShadowBasedonImprovedGeneralizedS-transform[J].ProgressinGeophysics,2010,25(6):2040-2044.(inChinese)
[22]劉蘭鋒,劉全新,雍學(xué)善,等.基于廣義S變換的低頻瞬時能量譜油氣檢測技術(shù)[J].天然氣地球科學(xué),2005,16(2):238-241.
LIULan-feng,LIUQuan-xin,YONGXue-shan,etal.TechnologyofLowFrequencyInstantaneousEnergySpectrumofGeneralizedSTransformBasedonOilandGasDetection[J].NaturalGasGeoscience,2005,16(2):238-241.(inChinese)
[23]盧育霞,石玉成,萬秀紅.近地表速度結(jié)構(gòu)對場地強(qiáng)地震動特征的影響[J].地震工程學(xué)報,2014,36(4):813-819.
LUYu-xia,SHIYu-cheng,WANXiu-hong,etal.InfluenceofNear-surfaceVelocityStructureonSiteCharacteristicsofStrongGroundMotion[J].ChinaEarthquakeEngineeringJournal,2014,36(4):813-819.(inChinese)
[24]WeaverRL,LobkisOI.UltrasonicwithoutaSource:ThermalFluctuationCorrelationsatMHzFrequencies[J].PhysRevLett,2001,87(13):134301.
[25]WeaverRL,LobkisOI.OntheEmergenceoftheGreen’sFunctionintheCorrelationsofaDiffuseField[J].JAcoustSocAm,2001,110:3011-3017.
[26]DerodeA,LaroseE,TanterM,etal.RecoveringtheGreen’sFunctionfromField-fieldCorrelationsinanOpenScatteringmedium[J].JAcoustSocAm,2003,113:2973-2976.
[27]SniederR.ExtractingtheGreen'sFunctionfromtheCorrelationofCodaWaves:ADerivationBasedonStationaryPhase[J].PhysRevE,2004:69.
[28]WapenaarK.RetrievingtheElastodynamicGreen'sFunctionofanArbitraryInhomogeneousMediumbyCrossCorrelation[J].PhysicalReviewLetters,2004,93:254301.