LabVIEW時頻分析在地空瞬變電磁信號處理中的應(yīng)用

摘 要： 地空瞬變電磁信號中含有大量干擾噪聲，尤其是晚期信號干擾影響較大，若不對其進(jìn)行濾波處理，將會影響瞬變信號的分析及反演成像的質(zhì)量。傳統(tǒng)的傅里葉變換不能分辨正常信號中夾帶的瞬態(tài)反?，F(xiàn)象，也不能展示其成分，無法滿足非平穩(wěn)信號的分析要求；而瞬變電磁信號是一種動態(tài)、瞬態(tài)的信號。因此采用短時傅里葉變換（STFT）進(jìn)行數(shù)據(jù)分析，并開發(fā)了一套基于LabVIEW的預(yù)處理軟件，用于現(xiàn)場快速分析地空瞬變電磁信號。實(shí)驗(yàn)結(jié)果表明，預(yù)處理軟件對地空瞬變電磁的分析和濾波處理是可行有效的，便于及時評價濾波效果，為進(jìn)一步確定信號濾波方法提供了依據(jù)。

關(guān)鍵詞： LabVIEW； 地空瞬變電磁信號； 短時傅里葉變換； 濾波處理

中圖分類號： TN919?34； TP391.76 文獻(xiàn)標(biāo)識碼： A 文章編號： 1004?373X（2015）24?0107?04

Application of LabVIEW time?frequency analysis in ground?to?air transient electromagnetic single processing

XIONG Jing， GAO Song， ZHANG Liang， REN Peng， ZHANG Zhenhai

（College of Information Science and Technology， Chengdu University of Technology， Chengdu 610059， China）

Abstract： There are a lot of interference noises in the ground?to?air transient electromagnetic signal， especially the late signal is greatly influenced by the interference. If the filtering processing isn’t adopted for the interference， it may influence on the analysis of transient signal and the quality of inversion imaging. The traditional Fourier transform can’t distinguish the transient abnormal phenomenon in normal signal and show its composition， so it can’t satisfy the analysis requirement of the non?stationary signal. Since the transient electromagnetic signal is a dynamic and transient one， the short time Fourier transform （STFT） is adopted to analyze the data. A set preprocessing software based on LabVIEW was developed to rapidly analyze the ground?to?air transient electromagnetic signal on site. The experimental results show that the preprocessing software is feasible and effective for the analysis of ground?to?air transient electromagnetism and filtering processing， and is convenient for evaluating the filtering effect in time. It provides the basis for further determining the signal filtering method.

Keywords： LabVIEW； ground?to?air transient electromagnetic signal； STFT； filtering processing

0 引 言

地空瞬變電磁法[1?2]以其受地形影響小、探測深度大、信噪比更高、空間分辨率更好等優(yōu)點(diǎn)成為近幾年來瞬變電磁研究領(lǐng)域的重點(diǎn)研究對象。而實(shí)測的電磁信號具有如下特點(diǎn)[3]：動態(tài)范圍大；中晚期信號幅值較小，幾乎淹沒在噪聲中；線圈運(yùn)動噪聲幅值較大，分布范圍廣；工頻干擾和飛機(jī)的電磁噪聲也是主要噪聲來源。若不對這些干擾進(jìn)行剔除，將直接影響后期反演成像的質(zhì)量。為了選擇合適的濾波方法，需對信號進(jìn)行分析處理，而傳統(tǒng)的傅里葉變換要么全時域，要么全頻域，適合分析平穩(wěn)信號；針對非平穩(wěn)信號，采用將時域和頻域聯(lián)合分析的時頻分析法，它能分辨正常信號中夾帶的瞬態(tài)反?，F(xiàn)象并展示其成分，比較適合動態(tài)、瞬態(tài)的信號處理，而在瞬變電磁法中得到的就是動態(tài)、瞬態(tài)的信號[4?5]。所以本文提出了在LabVIEW平臺上采用時頻分析法對地空瞬變電磁信號進(jìn)行分析，不僅為濾波方法的選擇提供了依據(jù)，而且對濾波效果進(jìn)行了及時評估。

1 短時傅里葉變換

1.1 工作原理

目前應(yīng)用較多的三種時頻分析方法為：短時傅里葉變換（STFT）、Wigner?Ville分布（WVD）和小波變換。WVD處理多頻率分量的信號時存在嚴(yán)重的交叉干擾現(xiàn)象，影響頻譜的物理解釋；小波變換中小波基的不同選擇將會有不同的效果，甚至差別很大。所以本文使用STFT對瞬變信號進(jìn)行分析。

STFT的思想是[6?7]：假定某個時間窗內(nèi)的信號是平穩(wěn)的，然后通過傅里葉變換確定時間窗內(nèi)信號的頻率特性，再移動窗函數(shù)分析下一個時間窗內(nèi)的信號，最終獲得整個信號的時頻特性。其定義為：

[STFTxt， f=-∞+∞xt′g*（t′-t）?e-j2πft′dt′]

式中：上標(biāo)“*”表示復(fù)共軛。當(dāng)窗函數(shù)g（t）=1時，STFT即為傳統(tǒng)的傅里葉變換。圖1（b）和圖1（c）對應(yīng)的頻率分別為100 Hz，200 Hz，300 Hz，采樣率為1 000，采樣數(shù)為300的正弦信號STFT時頻圖與FFT幅頻特性圖，從時頻圖中可明確判斷時間、頻率、幅度三者的關(guān)系，而FFT變換無法將時域的瞬變信息在頻域中表現(xiàn)出來。

在Matlab中可使用spectrogram函數(shù)對信號進(jìn)行STFT變換，圖1（d）為上述正弦仿真信號的STFT時頻圖，由圖可得信號的時頻特性與LabVIEW所得結(jié)果相同，但是Matlab

時頻分析需編寫相應(yīng)的代碼，而且需人工對指定區(qū)域的信號進(jìn)行截取，不能得到信號整體的時頻特性。所以本文使用LabVIEW對指定長度的信號進(jìn)行周期性STFT變換，從而分析信號的整體時頻特性。

圖1 正弦仿真信號時頻圖

1.2 窗函數(shù)的選擇

進(jìn)行STFT時，關(guān)鍵要考慮窗函數(shù)對STFT變換的影響。不同的窗函數(shù)會影響信號的主頻分布，選擇窗函數(shù)時遵循的準(zhǔn)則是：使窗函數(shù)的主瓣包含更多的能量，相應(yīng)旁瓣沒有多少能量；選擇窗函數(shù)時，不但需要分析窗函數(shù)的不同性質(zhì)，還需根據(jù)信號的不同性質(zhì)選擇合適的窗，如：矩形窗主瓣集中，旁瓣較高，有負(fù)旁瓣，對于頻率分辨率要求很高，譜估計(jì)精度要求不是很高的信號，處理時可選用矩形窗；高斯窗譜的主瓣較寬，故而頻率分辨率低，適合用來截?cái)嘁恍┓侵芷谛盘?；如果要分析窄帶信號且具有較強(qiáng)的干擾噪聲時，則應(yīng)選用旁瓣幅度較小的窗函數(shù)，如Hamming窗等[8?9]。圖2所示為分別使用矩形窗、Hanning窗、Hamming窗和Blackman對正弦信號進(jìn)行分析的時頻圖，由圖2可知，不同窗的時間和頻率分辨率不同。結(jié)合瞬變電磁信號頻譜分量較多，頻譜表現(xiàn)復(fù)雜，而且更多關(guān)注頻率點(diǎn)進(jìn)行濾波的目的，所以選用Hanning窗處理數(shù)據(jù)。

圖2 四種常見窗的比較

由STFT的不確定性原理[10?11]（?t·?f≥[14π]，?t，?f分別為時間和頻率分辨率）可知，STFT的時頻表示在時間分辨率和頻率分辨率之間是矛盾的，若要提高時間分辨率，則應(yīng)用短的時間窗函數(shù)；若要提高頻率分辨率，則應(yīng)選用長的時間窗函數(shù)。所以在實(shí)際分析過程中，需對時間分辨率與頻率分辨率進(jìn)行折中，以求兩者都相對較高。圖3為使用[18，14，12]數(shù)據(jù)長度的Hanning窗對正弦信號進(jìn)行STFT的時頻分析圖，[18]數(shù)據(jù)長度的時間分辨率高，但頻率分辨率低，而[12]數(shù)據(jù)長度則相反。在時間分辨率與頻譜分辨率間進(jìn)行折中，選擇窗長為[14]數(shù)據(jù)長度對瞬變信號進(jìn)行時頻分析。

2 數(shù)據(jù)的分析處理

對瞬變電磁信號進(jìn)行分析處理首先需讀取信號，并對其做時頻分析，得到主要頻段范圍；然后對主要頻段范圍的噪聲進(jìn)行濾波處理，可通過“雙通道譜測量VI”計(jì)算濾波前后的信號響應(yīng)和相干情況，選擇合適的濾波器；最后與濾波前后的譜分析圖進(jìn)行比較，分析濾波效果。本文需要的數(shù)據(jù)分析處理主要包括：圖4所示的讀取數(shù)據(jù)、STFT時頻分析、濾波處理、顯示及存儲。

圖3 不同長度Hanning窗時頻分析比較

圖4 數(shù)據(jù)分析處理結(jié)構(gòu)圖

2.1 讀取數(shù)據(jù)

數(shù)據(jù)采集中存儲格式的選擇需與數(shù)據(jù)格式和采樣率等因素相匹配。TDMS文件是NI公司最新推出的數(shù)據(jù)管理系統(tǒng)。使用這種格式的數(shù)據(jù)可以在測試結(jié)果中存儲大量的數(shù)據(jù)，并可以將測量數(shù)據(jù)導(dǎo)入Excel中，所以本文也使用TDMS格式存儲采集數(shù)據(jù)，并使用讀取測量文件VI讀取已存儲的瞬變信號，方便后期分析處理。

2.2 STFT時頻分析

讀取得到的瞬變信號需進(jìn)行STFT時頻分析得到噪聲的主要頻段范圍和出現(xiàn)的時間，為后期的濾波和分析異常體的地點(diǎn)提供依據(jù)。通過調(diào)用STFT時頻圖.VI對信號進(jìn)行STFT時頻分析，根據(jù)前面的分析，選擇窗長為[14]頻率區(qū)間的Hanning窗進(jìn)行分析，并根據(jù)實(shí)際采樣率設(shè)置STFT時頻圖的橫、縱坐標(biāo)顯示范圍。本文通過調(diào)用時頻圖的色碼表屬性設(shè)置幅值強(qiáng)度表示的顏色，為了使顯示美麗直觀，本文將強(qiáng)度設(shè)置為類似彩虹的顏色，具體實(shí)現(xiàn)如圖5所示。本文還調(diào)用了自功率.VI分析指定頻率區(qū)間內(nèi)信號幅值和頻率的相互關(guān)系。

2.3 濾波處理

使用STFT時頻分析得到信號噪聲的主要頻段范圍后需要進(jìn)行濾波處理，本文主要對指定頻段范圍的信號進(jìn)行陷波處理，所以本文使用條件結(jié)構(gòu)控制帶阻濾波器的啟動與關(guān)閉。選擇濾波器時可通過使用雙通道譜測量VI將原始信號與濾波后信號分別輸入A和B通道中，通過對各通道中的信號有序?qū)Φ姆治觯?jì)算濾波信號與未濾波信號的信號響應(yīng)和相干情況，從而得到適合的濾波器。本文使用10階Butterworth帶阻濾波器進(jìn)行濾波處理，為后期硬件設(shè)計(jì)時濾波器的選擇提供了依據(jù)，具有方便、快捷、節(jié)約資源等優(yōu)點(diǎn)。

圖5 STFT時頻分析

2.4 顯示與存儲

顯示與存儲部分主要是對原始信號與濾波信號進(jìn)行顯示比較，并可通過條件結(jié)構(gòu)的調(diào)用選擇是否對濾波后的信號進(jìn)行存儲，方便后期做進(jìn)一步分析處理，程序會在原信號文件的路徑下自動創(chuàng)建一個TDMS文件保存濾波后的信號。

3 實(shí)驗(yàn)測試

2015年1月15日在成都市高威公園，采用成都理工大學(xué)自主研發(fā)的微功率發(fā)射儀以20 A，25 Hz的電流進(jìn)行發(fā)射，空中通過無人機(jī)搭載接收系統(tǒng)進(jìn)行電磁信號采集分析實(shí)驗(yàn)。接收系統(tǒng)采用24位同步采集模塊NI 9239，以25 kHz采樣率，全波段連續(xù)采集，并通過cDAQ?9191無線機(jī)箱以WiFi形式進(jìn)行地面實(shí)時監(jiān)控。當(dāng)線圈距離地面20 m，放大100倍時采集到的信號如圖6（a）所示。使用[14]頻率區(qū)間的Hanning窗對指定長度的信號進(jìn)行周期性STFT變換得到的時頻圖，如圖6（b）所示，對應(yīng)的頻率區(qū)間為25 000×0.03=750。由圖6可知，信號的主要能量集中在100～200 Hz內(nèi)，使用10階的Butterworth帶阻濾波器進(jìn)行處理后的時域信號和時頻圖如圖6（c）所示，由圖可知100～200 Hz信號的幅值趨于0，表明濾波已達(dá)到指定要求，濾波后信號如圖6（d）所示。

4 結(jié) 語

實(shí)驗(yàn)結(jié)果表明，時頻分析克服了傳統(tǒng)快速傅里葉變換分析瞬變信號的不足，準(zhǔn)確知道各頻率及幅值隨時間的變化，而LabVIEW相較于Matlab可實(shí)時對采集的信號進(jìn)行分析，也可對保存的數(shù)據(jù)的指定長度進(jìn)行周期性分析。而本文設(shè)計(jì)的基于LabVIEW的地空瞬變電磁信號STFT分析處理軟件綜合了兩者的優(yōu)勢，既能有效確定干擾信號，為濾波處理提供依據(jù)，又大大減少了開發(fā)時間，具有方便、快捷、直觀顯示等優(yōu)點(diǎn)，適合需要快速評估實(shí)驗(yàn)成效的野外作業(yè)。

圖6 濾波前后信號比較

參考文獻(xiàn)

[1] 毛立峰，王緒本.磁性源半航空瞬變電磁矢量有限元法三維正演[C]//中國地球科學(xué)聯(lián)合學(xué)術(shù)年會：電磁地球物理學(xué)研究應(yīng)用及其新進(jìn)展論文集.北京：中國地球物理學(xué)會，2014：177?179.

[2] 李肅義，林君，陽貴紅，等.電性源時域地空電磁數(shù)據(jù)小波去噪方法研究[J].地球物理學(xué)報(bào)，2013，56（9）：3145?3152.

[3] 陽貴紅.時域電性源地空電磁探測數(shù)據(jù)預(yù)處理研究[D].長春：吉林大學(xué)，2012.

[4] 冉云.瞬變電磁測深數(shù)據(jù)的小波分析處理與電性層劃分[D].太原：太原理工大學(xué)，2010.

[5] 何臘梅，王宇航.STFT與FIR在航空瞬變電磁數(shù)據(jù)處理中的應(yīng)用[J].工程地球物理學(xué)報(bào)，2013，10（1）：15?20.

[6] 王曉華，彭程，梁靜.STFT時頻分析算法研究及其在列車震動分析中的應(yīng)用[J].噪聲與振動控制，2010（1）：65?68.

[7] 遲華山，王紅星，郭奇，等.短時傅里葉變換在線性調(diào)頻信號時頻濾波中的應(yīng)用[J].電訊技術(shù)，2012，52（2）：155?159.

[8] 毛青春，徐分亮.窗函數(shù)及其應(yīng)用[J].中國水運(yùn)，2007，7（2）：230?232.

[9] 肖雪亮.窗函數(shù)應(yīng)用中的幾個問題[J].動 態(tài)分析與測試技術(shù)，1989（2）：10?19.

[10] 呂東偉.吊艙式時間域直升機(jī)航空電磁數(shù)據(jù)處理方法研究[D].成都：成都理工大學(xué)，2011.

[11] 何臘梅.航空瞬變電磁數(shù)據(jù)時頻分析及去噪方法研究[D].成都：成都理工大學(xué)，2013.