甘 露 ，周 龍 ，尤新革 ，2

（1.武漢工業(yè)學(xué)院 湖北 武漢 430023；2.華中科技大學(xué) 湖北 武漢 430074）

探地雷達(dá)的信號(hào)處理結(jié)果決定著其在應(yīng)用領(lǐng)域中的有效性，其處理過程主要是通過分析探地雷達(dá)接收到的有效反射電磁波的特征來推斷地下介質(zhì)的空間分布狀態(tài)。其中關(guān)于回波信號(hào)特征提取是最重要、最關(guān)鍵也最困難的問題之一，它直接關(guān)系到對(duì)地下目標(biāo)探測(cè)的準(zhǔn)確性和可靠性。傳統(tǒng)探地雷達(dá)資料常用的分析方法是Fourier分析[1-2]，然而，在傳統(tǒng)的探地雷達(dá)回波信號(hào)的表達(dá)中，一般都認(rèn)為回波信號(hào)是平穩(wěn)的。但是，在實(shí)際應(yīng)用中，大地的電磁特性呈現(xiàn)隨機(jī)性，導(dǎo)致接收的回波信號(hào)具有不同的幅度和相位，因此探地雷達(dá)的回波信號(hào)具有非平穩(wěn)的特征。在回波數(shù)據(jù)中，收發(fā)天線直接耦合波、地面反射波與它們之間多次波疊加，引起的反射波稱統(tǒng)稱為直達(dá)波。直達(dá)波的能量很大，某些情況下甚至能夠掩蓋目標(biāo)。

希爾伯特一黃變換（Hilbert-Huang Transform，HHT）是一種新的非平穩(wěn)信號(hào)的處理技術(shù)，該方法由經(jīng)驗(yàn)?zāi)B(tài)分解（Empirical Mode Decomposition，EMD）與Hilbert譜分析（Hilbert Spectral Analysis，HSA）兩部分組成：任意的非平穩(wěn)信號(hào)首先經(jīng)過EMD方法處理后被分解為若干個(gè)本征模態(tài)函數(shù)（Intrinsic Mode Function，IMF），然后對(duì)每個(gè)IMF分量進(jìn)行Hilbert譜分析得到相應(yīng)分量的Hilbert譜，最后匯總所有IMF分量Hilbert譜就得到了原始非平穩(wěn)信號(hào)的Hilbert譜[3-5]。按照這種方法分解所得到的IMF分量也具備明確的物理意義。這些IMF刻畫了信號(hào)在每一個(gè)局部的振蕩結(jié)構(gòu)或頻率結(jié)構(gòu)。為此，本文以探地雷達(dá)實(shí)測(cè)數(shù)據(jù)位研究對(duì)象，引入HHT將各種頻率成分以IMF的形式從原始信號(hào)中分離出來，對(duì)探地雷達(dá)回波信號(hào)進(jìn)行特征提取，并給出探地雷達(dá)信號(hào)EMD分解算法。

1 HHT方法

Hilbert-Huang變換按照如下的2個(gè)步驟來分析數(shù)據(jù)。1）用經(jīng)驗(yàn)?zāi)Ｊ椒纸鈱?duì)數(shù)據(jù)作預(yù)處理，通過分解我們可以得到一組“本征模式函數(shù)”，也即我們可以用得到的“基”展開數(shù)據(jù)。2）將分解得到的“本征模式函數(shù)”作希爾伯特變換并構(gòu)造能量-時(shí)間-頻率分布，即希爾伯特譜，它將保存事件的時(shí)間局部性。換句話說，我們需要的是瞬時(shí)的頻率和能量，而不是傅立葉譜分析中的全局頻率和能量。經(jīng)驗(yàn)?zāi)Ｊ椒纸猓‥MD）的本質(zhì)是根據(jù)數(shù)據(jù)的特征時(shí)間尺度來經(jīng)驗(yàn)地識(shí)別固有振蕩模式，然后據(jù)此來分解數(shù)據(jù)。不同于信號(hào)的Fourier分解和小波分解.EMD沒有確定的基，它的“基”是根據(jù)信號(hào)而自適應(yīng)產(chǎn)生的，這使得它不僅具有很高的分解效率，同時(shí)也具有良好的時(shí)頻局部性。EMD分解的具體步驟如下：

設(shè) x（t）為原始數(shù)據(jù)。

步驟 1 初始化:令 r0（t）=x（t），i=1

步驟 2 提取第i個(gè)IMF：

1）初始化:h0（t）=xi-1（t），j=1

2）求 hj-1（t）的局部極大值和局部極小值；

3）以局部極大值作為節(jié)點(diǎn)作三次樣條插值形成hj-1（t）的上包絡(luò) u pper（t），類似的形成下包絡(luò) l ower（t）；

5）令 hj（t）=hj-1（t）-mj-1（t）；

6）若hj（t）的穿越零點(diǎn)的個(gè)數(shù)與極值點(diǎn)的個(gè)數(shù)最多相差為 1 ，且

小于某個(gè)預(yù)先給定的常數(shù)，則 I MFi（t）=hj（t），否則轉(zhuǎn)至2），并令 j =j+1；

步驟 3 令 ri（t）=ri-1（t）-IMFi（t）；

步驟 4 若ri（t）至少有2個(gè)極值點(diǎn)，轉(zhuǎn)至步驟2并令i=i+1；否則，ri（t）分解結(jié)束，為余量。這樣，原信號(hào) X（t）可表示為

2 系統(tǒng)硬件設(shè)計(jì)

雷達(dá)系統(tǒng)硬件主要由控制顯示單元、雷達(dá)主機(jī)和雷達(dá)天線及發(fā)射接收系統(tǒng)3部分組成，系統(tǒng)硬件結(jié)構(gòu)圖如圖1所示[6]。系統(tǒng)采用“PC＋DSP＋FPGA”三層控制方案。底層的控制器由以Altera公司的EP1K30TC144-3為核心的高速控制邏輯電路構(gòu)成，該芯片具有1728個(gè)邏輯單元和171個(gè)I/O引腳，可以滿足雷達(dá)時(shí)序系統(tǒng)的存儲(chǔ)容量需要和接口需要。FPGA不僅可以實(shí)現(xiàn)對(duì)時(shí)變?cè)鲆娣糯笃饕约澳?shù)轉(zhuǎn)換器的同步控制，還同時(shí)生成發(fā)射和取樣觸發(fā)脈沖源。該脈沖源經(jīng)放大整形后送往天線系統(tǒng)；中層控制由AD公司的ADSP21065L芯片完成，該芯片采用超級(jí)哈佛（Super Harvard）結(jié)構(gòu)，其峰值運(yùn)算速度可以達(dá)到198MFLOPS，它還提供了多種外部接口，具有10個(gè)DMA通道，可以為雷達(dá)數(shù)據(jù)的采集傳送提供快速通道。DSP不僅可實(shí)現(xiàn)與PC機(jī)的通信，接收并執(zhí)行由PC機(jī)發(fā)送的控制命令，實(shí)施對(duì)底層控制器的控制，同時(shí)還可將采集數(shù)據(jù)讀出并迅速上傳至PC機(jī)進(jìn)行處理；上層控制由PC機(jī)完成對(duì)探地雷達(dá)功能和參數(shù)的設(shè)置，并向DSP發(fā)送控制命令，同時(shí)還實(shí)時(shí)地接收DSP上傳的回波信號(hào)數(shù)據(jù)，進(jìn)行基于HHT的時(shí)頻分析并顯示圖象及波形。

3 基于HHT的探地雷達(dá)信號(hào)特征提取

根據(jù)上述分析，使用HHT方法進(jìn)行探地雷達(dá)信號(hào)特征提取[7]的過程如下：

圖1 雷達(dá)硬件結(jié)構(gòu)圖Fig.1 Structure diagram of the radar hardware system

1）移動(dòng)雷達(dá)對(duì)地下埋藏物的狀態(tài)進(jìn)行數(shù)據(jù)采集，得到探地雷達(dá)回波信號(hào)的時(shí)域波形。

2）1）中所得的時(shí)域波形為一維信號(hào)，信號(hào)中包含直達(dá)波及目標(biāo)信號(hào)。設(shè)接收到信號(hào)表示為d（x，t），其中x表示測(cè)點(diǎn)，t表示時(shí)間，則任意測(cè)點(diǎn)x處探地雷達(dá)信號(hào)的EMD分解可表示為：

式中，k表示第級(jí)IMF分量。在對(duì)實(shí)際探地雷達(dá)信號(hào)進(jìn)行EMD分解時(shí)，當(dāng)殘留分量剖面的能量遠(yuǎn)遠(yuǎn)小于原始剖面能量時(shí)，EMD分解就可以停止，則確定了其分解級(jí)數(shù)。

3）通過EMD分解，得到具有窄帶特性的內(nèi)蘊(yùn)模式函數(shù)剖面imfk（x，t），進(jìn)一步求取其瞬時(shí)頻率作為探地雷達(dá)回波數(shù)據(jù)特征。

4 實(shí)驗(yàn)結(jié)果與分析

文中選取探地雷達(dá)實(shí)測(cè)數(shù)據(jù)中的一組數(shù)據(jù)進(jìn)行測(cè)試。實(shí)驗(yàn)場(chǎng)地選在戶外的一段水泥路面上，已知該路段下方有一水管道。將探地雷達(dá)工作參數(shù)的設(shè)置為：發(fā)射脈沖重復(fù)頻率為300 kHz，掃描速率為 256 scans/Sec，采樣率為 512 samp/Scan。采集到的回波信號(hào)如圖2所示。

圖2 探地雷達(dá)回波信號(hào)Fig.2 GPR echo signal

由圖2中可以看出在100點(diǎn)左右與400點(diǎn)左右均存在著明顯的幅值越變，而且幅值變化較劇烈。而位于前端的即為直達(dá)波信號(hào)，后端的為埋藏物信號(hào)。顯然，當(dāng)直達(dá)波信號(hào)較強(qiáng)烈時(shí)，容易引起誤判。對(duì)圖2所示信號(hào)進(jìn)行EMD分解，即可得到其IMF分量，如圖3所示。

圖3 回波信號(hào)的IMF分量Fig.3 IMF component of echo signal

由圖3可見，對(duì)圖2中的回波信號(hào)用EMD方法分解可得到3個(gè)IMF（imf1，imf2，imf3）分量和一個(gè)殘余分量（res），分解較為徹底，imf1～imf3頻率由高到低。然而，無法分辨出直達(dá)波與目標(biāo)信號(hào)。進(jìn)一步求取圖3中IMF分量的瞬時(shí)頻率，如圖4所示。

圖4 回波信號(hào)IMF分量的瞬時(shí)頻率Fig.4 EMD result of echo signal

由圖4可見，imf1的瞬時(shí)頻率表征了直達(dá)波與埋藏物信號(hào)處存在明顯變化。而imf2的瞬時(shí)頻率僅表征了埋藏物所處位置的頻率變化，說明在雷達(dá)經(jīng)過埋藏物上方時(shí)的回波信號(hào)頻率為明顯的高頻分量，體現(xiàn)了沖擊頻率。imf3和imf4的瞬時(shí)頻率分量則不包含過多信息。實(shí)驗(yàn)結(jié)果表明：IMF的瞬時(shí)頻率能有效的提取出信號(hào)的主要特征信息，較好消除直達(dá)波的影響。

5 結(jié) 論

HHT作為一種新的信號(hào)處理方法，在非平穩(wěn)非線性信號(hào)的分析上有著獨(dú)特的優(yōu)勢(shì)，本文將HHT方法應(yīng)用于探地雷達(dá)信號(hào)特征提取，并通過對(duì)實(shí)測(cè)數(shù)據(jù)的分析，得到各階固有模態(tài)函數(shù)信號(hào)。分析結(jié)果表明，采用HHT方法得到的IMF分量的瞬時(shí)頻率可較好的消除直達(dá)波影響，提取出地下埋藏物的特征，可供信號(hào)的進(jìn)一步資料解釋。

[1]Maida A，Pennock S，Shepherd P.Improving ground penetrating radar signal analysis through FFT superimposition[J].Antennas and Propagation Society International Symposium，2005 IEEE，2005（4）:118.

[2]SONG Jia-yu，LIU Qing-Huo，Torrione P，et al.Two-dimensional and three-dimensional NUFFT migration method for landmine detection using ground-penetrating Radar[J].Geoscience and Remote Sensing，IEEE Transactions on，2006（6）:1462-1466.

[3]Huang N E，Shen Z，Long S R.A new view of nonlinear water waves：the Hilbert spectrum[J].Annual Review of Fluid Mechanics，1999（31）:417-457.

[4]Montesinos M E，Munoz-Cobo J L，Perez C.Hilbert-Huang analysis of BWR neutron detector signals：application to DR calculation and to corrupted signal analysis[J].Annals of Nuclear Energy，2003（30）：715-727.

[5]Phillips S C，Gledhill R J，Essex J W，et al.Application of the Hilbert-Huang transform to the analysis of molecular dynamics simulations[J].Journal of Physical Chemistry Ser A，2003（24）:4869-4876.

[6]GAN Lu，GAN Liang-cai，TIAN Mao.Step-system of high resolution ground penetrating radar[J].Chinese Journal of Radio Science，2008（3）:555-559.

[7]孫百紅，宋少偉.基于小波分析的發(fā)動(dòng)機(jī)轉(zhuǎn)動(dòng)慣量測(cè)量信號(hào)特征提取[J].火箭推進(jìn)，2010（06）:52-55.SUN Bai-hong，SONG Shao-wei.Feature extraction for measurement signal of engine rotary inertia moment based on wavelet analysis[J].Journal of Rocket Propulsion，2010（06）:52-55.