基于波場(chǎng)互相關(guān)的探地雷達(dá)快速自聚焦成像

楊忠委 郭聰隆 孫浩然 蘭 天 楊小鵬

(1. 北京理工大學(xué)信息與電子學(xué)院, 北京 100081; 2. 北京理工大學(xué)重慶創(chuàng)新中心, 重慶 401120)

1 引言

探地雷達(dá)偏移成像技術(shù)通過(guò)改變雷達(dá)回波中能量的分布情況,可以更為直觀地表現(xiàn)出地下場(chǎng)景及其中目標(biāo)的幾何或物理特征,大幅降低了目標(biāo)檢測(cè)和識(shí)別的難度,便于目標(biāo)解譯工作的開(kāi)展。伴隨著雷達(dá)成像技術(shù)的飛速發(fā)展,國(guó)內(nèi)外學(xué)者提出了許多探地雷達(dá)偏移成像算法,并取得了很好的應(yīng)用效果[1-2]。但這些算法均是建立在介質(zhì)參數(shù)已知的條件下,當(dāng)介質(zhì)參數(shù)未知時(shí),將會(huì)產(chǎn)生成像散焦和深度偏移等問(wèn)題。為解決此類問(wèn)題,國(guó)內(nèi)外學(xué)者研究了介質(zhì)參數(shù)測(cè)量方法[3-5]。Ahmad等基于高階統(tǒng)計(jì)特性提出了一種自聚焦成像算法,當(dāng)介質(zhì)參數(shù)估計(jì)誤差較小且目標(biāo)處于遠(yuǎn)場(chǎng)區(qū)時(shí)能夠取得很好的自聚焦成像效果[6]。崔國(guó)龍等提出了用于時(shí)分MIMO穿墻雷達(dá)的后向投影成像算法,通過(guò)建立介質(zhì)參數(shù)已知條件下的目標(biāo)回波模型,使用費(fèi)馬定理和斯涅耳定律計(jì)算聚焦補(bǔ)償量,然后通過(guò)后向投影算法獲得高質(zhì)量圖像[7]。Dehmollaian與Sarabandi通過(guò)求解非線性優(yōu)化問(wèn)題得到介質(zhì)參數(shù),然后使用合成孔徑成像算法來(lái)實(shí)現(xiàn)目標(biāo)的自聚焦成像[8]。但以上這些算法往往存在計(jì)算復(fù)雜度大且穩(wěn)定性差等問(wèn)題。

基于逆時(shí)偏移成像的思想,本文提出了一種探地雷達(dá)快速自聚焦成像算法。由于逆時(shí)偏移不存在對(duì)波動(dòng)方程的近似,允許電磁波在任意方向進(jìn)行傳播,能夠?qū)θ我忸愋偷碾姶挪▓?chǎng)進(jìn)行成像,因此逆時(shí)偏移成像算法相比其他偏移方法具有無(wú)法比擬的優(yōu)越性[9-11]。傳統(tǒng)時(shí)域逆時(shí)偏移主要分為波場(chǎng)延拓和成像兩部分,主要過(guò)程如下。首先,將激勵(lì)源信號(hào)沿時(shí)間方向延拓,得到正傳波場(chǎng);將天線在接收位置處接收的數(shù)據(jù)進(jìn)行時(shí)間反轉(zhuǎn),然后反傳播到探測(cè)媒質(zhì)中去,得到反傳播場(chǎng)。正傳和反傳波場(chǎng)通常采用時(shí)域有限差分來(lái)進(jìn)行計(jì)算。然后,對(duì)正傳與反傳波場(chǎng)按成像條件提取成像值,得到單對(duì)收發(fā)天線位置(單炮)偏移結(jié)果。最后,將每一炮的偏移結(jié)果疊加起來(lái),即得到完整的偏移圖像?？梢钥闯鰝鹘y(tǒng)時(shí)域逆時(shí)偏移算法需要多次的正向和反向波場(chǎng)外推,然后應(yīng)用成像條件將每個(gè)炮點(diǎn)的成像結(jié)果進(jìn)行累加。因此時(shí)域逆時(shí)偏移需要很大的計(jì)算與存儲(chǔ)資源,導(dǎo)致其并不能在實(shí)際場(chǎng)景中對(duì)大數(shù)據(jù)量進(jìn)行實(shí)時(shí)成像處理,故并未得到廣泛應(yīng)用。

逆時(shí)偏移成像條件有很多種,如激發(fā)時(shí)間成像條件、上下行波振幅比成像條件、互相關(guān)成像條件等[12],其中零延遲互相關(guān)是逆時(shí)偏移最常用的成像條件,可看作收發(fā)波場(chǎng)之間的匹配濾波過(guò)程?；谠摮上駰l件,本文推導(dǎo)了共偏移距探地雷達(dá)數(shù)據(jù)的頻域逆時(shí)偏移表達(dá)式。通過(guò)應(yīng)用分層介質(zhì)格林函數(shù)的平移不變性,避免了多次波場(chǎng)外推過(guò)程,減小了對(duì)計(jì)算和存儲(chǔ)資源的需求。考慮到實(shí)際場(chǎng)景中,數(shù)據(jù)采集位置并不是均勻分布,采用非均勻快速傅里葉變換將空域非均勻采樣數(shù)據(jù)變換到均勻波數(shù)域,實(shí)現(xiàn)了頻率域多個(gè)炮點(diǎn)成像結(jié)果的快速累加。由于所使用的成像條件為零延時(shí)互相關(guān),因此當(dāng)介質(zhì)參數(shù)準(zhǔn)確估計(jì)時(shí),零成像時(shí)刻所得結(jié)果便是最優(yōu)聚焦圖像。然而當(dāng)介質(zhì)參數(shù)未被準(zhǔn)確估計(jì)時(shí),目標(biāo)會(huì)出現(xiàn)距離向偏移、方位向散焦的問(wèn)題,算法聚焦性能嚴(yán)重下降。從逆時(shí)偏移成像算法原理出發(fā),這種現(xiàn)象可以解釋為成像時(shí)刻的不準(zhǔn)確選取。當(dāng)成像時(shí)刻提前或滯后最優(yōu)成像時(shí)刻,所得圖像均會(huì)出現(xiàn)一定程度的散焦。鑒于此,所提算法向逆時(shí)偏移成像表達(dá)式中引入了時(shí)間相位補(bǔ)償因子,通過(guò)計(jì)算偏移圖像的熵值并使其最小化來(lái)確定補(bǔ)償量。因此所提算法不需要已知精確的介質(zhì)先驗(yàn)信息就能實(shí)現(xiàn)良好的自聚焦成像效果。

2 快速自聚焦成像算法

2.1 基于波場(chǎng)互相關(guān)的快速頻域成像

基于零延時(shí)互相關(guān)成像條件的逆時(shí)偏移算法可以表示為:

(1)

其中us(·)與ur(·)分別是發(fā)射和接收天線位置處的發(fā)射信號(hào)前向遞推波場(chǎng)與接收信號(hào)時(shí)間反轉(zhuǎn)后的后向遞推波場(chǎng);Tmax為接收到的A掃描數(shù)據(jù)時(shí)長(zhǎng),r=(x,y,z)為成像網(wǎng)格點(diǎn)坐標(biāo),對(duì)于二維垂直截面和水平切片成像則網(wǎng)格坐標(biāo)分別取r=(x,z),r=(x,y);rs為發(fā)射天線坐標(biāo),rd=(xd,yd,zd=0)是收發(fā)天線之間的坐標(biāo)差值,zd=0是因?yàn)槭瞻l(fā)天線對(duì)僅在xoy平面步進(jìn)移動(dòng)。對(duì)每個(gè)A掃描數(shù)據(jù)進(jìn)行時(shí)域逆時(shí)偏移并進(jìn)行累加便得到最終成像結(jié)果I(r)。式(1)可以寫(xiě)成以下卷積形式:

(2)

根據(jù)卷積定理,函數(shù)卷積的傅里葉變換是函數(shù)傅里葉變換的乘積[13],則式(2)可表示為:

(3)

其中Us(·)與Ur(·)分別為前向與后向時(shí)域波場(chǎng)值ur(·)、us(·)的頻域形式。根據(jù)電磁場(chǎng)理論,頻域波場(chǎng)是分層介質(zhì)格林函數(shù)與相應(yīng)收發(fā)位置處激勵(lì)源頻譜的乘積,即:

(4)

(5)

(6)

考慮到分層介質(zhì)格林函數(shù)的平移不變性,即:

G(r,rs,ω)=G(r-rs,r0,ω)

(7)

其中G(r,r0,ω)為參考坐標(biāo)點(diǎn)r0=(0,0,0)處的頻域分層介質(zhì)格林函數(shù)。

(8)

其中,

M(r,r0,ω)=G(r,r0,ω)G(r+rd,r0,ω)

(9)

在某些遠(yuǎn)場(chǎng)探測(cè)場(chǎng)景中,收發(fā)天線間距很小,式(9)可簡(jiǎn)化為

M(r,r0,ω)=G2(r,r0,ω)

(10)

可以發(fā)現(xiàn)式(8)是對(duì)rs+rd的卷積運(yùn)算,對(duì)其應(yīng)用卷積定理可得:

(11)

從式(11)可以看出,基于波場(chǎng)互相關(guān)的快速頻域成像算法主要分為以下幾步:

(1)根據(jù)成像場(chǎng)景先驗(yàn)信息,計(jì)算數(shù)據(jù)采集平面中心點(diǎn)處分層介質(zhì)格林函數(shù)分布G(·);

2.2 基于最小熵準(zhǔn)則的自聚焦處理

應(yīng)用上述算法的先決條件是必須精確已知介質(zhì)的介電常數(shù)等先驗(yàn)信息,但在實(shí)際應(yīng)用場(chǎng)景中,介質(zhì)參數(shù)并非精確已知,從而導(dǎo)致成像結(jié)果會(huì)出現(xiàn)一定程度的散焦。為準(zhǔn)確獲得介質(zhì)相對(duì)介電常數(shù)εr,可以通過(guò)非線性優(yōu)化問(wèn)題進(jìn)行求解,但優(yōu)化過(guò)程一般比較費(fèi)時(shí)并且不穩(wěn)定[6]。

為解決該問(wèn)題,向式(11)引入一個(gè)時(shí)間相位因子,式(11)可以認(rèn)為是成像時(shí)刻t=0的成像結(jié)果,即:

(12)

當(dāng)介質(zhì)參數(shù)準(zhǔn)確估計(jì)時(shí),則最優(yōu)聚焦時(shí)刻為t=0;如果當(dāng)介質(zhì)參數(shù)估計(jì)誤差較大,則最優(yōu)成像時(shí)刻并非t=0,故將式(12)作如下修正:

(13)

因此當(dāng)介質(zhì)參數(shù)估計(jì)不準(zhǔn)確時(shí),通過(guò)選擇最優(yōu)的成像時(shí)刻t,即可得到目標(biāo)的最優(yōu)聚焦圖像,成像時(shí)刻提前或滯后最優(yōu)成像時(shí)刻均會(huì)導(dǎo)致成像結(jié)果散焦。最優(yōu)成像時(shí)刻通過(guò)最小熵準(zhǔn)則來(lái)確定,即計(jì)算不同成像時(shí)刻的偏移圖像熵值,圖像熵最小的成像時(shí)刻便是最優(yōu)成像時(shí)刻,所對(duì)應(yīng)的圖像會(huì)作為最終自聚焦成像結(jié)果并輸出。圖像熵定義為:

(14)

其中xij為大小為m×n的圖像在像素點(diǎn)(i,j)處的像素值。

3 仿真分析

為論證所提自聚焦成像算法的有效性,本節(jié)將通過(guò)對(duì)場(chǎng)景模擬回波數(shù)據(jù)進(jìn)行自聚焦成像處理來(lái)說(shuō)明。場(chǎng)景回波數(shù)據(jù)采用探地雷達(dá)正演軟件gprMax進(jìn)行模擬[14],激勵(lì)信號(hào)為中心頻率為1.5 GHz的Ricker脈沖。其中地下介質(zhì)的相對(duì)介電常數(shù)為6,相對(duì)磁導(dǎo)率為1,電導(dǎo)率為0.001 S/m,收發(fā)天線以1.6 cm的步進(jìn)間隔在介質(zhì)表面移動(dòng)。

第一個(gè)仿真場(chǎng)景為二維半空間模型,如圖1所示,位于場(chǎng)景中心的目標(biāo)為半徑3.5 cm金屬圓柱體。收發(fā)天線緊貼地表,從x=0.2 m處步進(jìn)移動(dòng)至x=1.8 m。利用基于波場(chǎng)互相關(guān)的成像算法對(duì)仿真回波數(shù)據(jù)進(jìn)行成像處理,圖2為介質(zhì)相對(duì)介電常數(shù)εr分別估計(jì)為4、6、8、10時(shí)的零時(shí)刻成像結(jié)果?？梢钥闯?在零成像時(shí)刻,當(dāng)εr準(zhǔn)確估計(jì)為6時(shí),所得圖像聚焦效果良好;當(dāng)εr估計(jì)值低于或高于其真實(shí)值時(shí),點(diǎn)目標(biāo)位置會(huì)偏離其真實(shí)位置,成像結(jié)果會(huì)出現(xiàn)一定程度散焦,且散焦程度會(huì)隨εr估計(jì)誤差的增大而增強(qiáng)。

圖1 點(diǎn)目標(biāo)二維場(chǎng)景模型Fig.1 The two-dimensional model of point target

圖2 不同介電常數(shù)估計(jì)值零成像時(shí)刻算法 對(duì)點(diǎn)目標(biāo)的成像結(jié)果Fig.2 Imaging results of point target for different estimated relative dielectric constant at focusing time zero

在實(shí)際應(yīng)用中,介質(zhì)的相對(duì)介電常數(shù)并不能準(zhǔn)確估計(jì),故所提算法在零成像時(shí)刻的聚焦圖像必然會(huì)出現(xiàn)散焦現(xiàn)象。為論證所提自聚焦成像算法性能,估計(jì)εr分別為4、6、8、10時(shí)(真實(shí)εr=6),利用所提自聚焦成像算法對(duì)回波數(shù)據(jù)進(jìn)行成像處理。聚焦圖像熵值隨成像時(shí)刻變化曲線如圖3所示,最優(yōu)成像時(shí)刻t所對(duì)應(yīng)的算法自聚焦成像結(jié)果如圖4所示。

圖3 圖像熵隨成像時(shí)刻的變化曲線Fig.3 Image entropy at different focusing time

圖4 不同介電常數(shù)估計(jì)值點(diǎn)目標(biāo)自聚焦成像結(jié)果Fig.4 The autofocusing imaging results of point target for different estimated relative dielectric constant

其中,圖3中的時(shí)間步進(jìn)為0.1 ns,從圖中可知,利用所提自聚焦成像算法,當(dāng)介質(zhì)εr估計(jì)值小于真實(shí)值時(shí),最優(yōu)成像時(shí)刻滯后于0時(shí)刻;當(dāng)介質(zhì)εr估計(jì)值等于真實(shí)值時(shí),最優(yōu)成像時(shí)刻即為0時(shí)刻;當(dāng)介質(zhì)εr估計(jì)值大于真實(shí)值時(shí),最優(yōu)成像時(shí)刻提前于0時(shí)刻。圖4為不同的介質(zhì)εr估計(jì)值,單目標(biāo)場(chǎng)景回波經(jīng)所提算法處理后的自聚焦圖像?？梢园l(fā)現(xiàn),盡管介電常數(shù)存在誤差,所提算法仍能將點(diǎn)目標(biāo)有效聚焦,驗(yàn)證了所提算法的自聚焦成像性能,但介電常數(shù)誤差越大,所得圖像旁瓣電平有所提升,聚焦效果相對(duì)變差。

為進(jìn)一步論證所提算法的自聚焦成像性能,仿真了V形目標(biāo)二維場(chǎng)景回波,其數(shù)據(jù)采集方式與單目標(biāo)場(chǎng)景相同,如圖5所示。圖6是不同的介質(zhì)εr估計(jì)值,V形目標(biāo)仿真回波數(shù)據(jù)經(jīng)所提自聚焦成像算法處理后的圖像。對(duì)比成像結(jié)果可知,在介質(zhì)介電常數(shù)存在不同程度估計(jì)誤差的情況下,所提算法自聚焦處理圖像仍能夠反映目標(biāo)的原始幾何形狀,也進(jìn)一步驗(yàn)證了所提算法的有效性。但是,隨著介電常數(shù)誤差的增大,V形目標(biāo)會(huì)發(fā)生一定程度形變,即張角有所變化,且目標(biāo)能量的聚焦效果也有所降低。

圖5 V形目標(biāo)二維場(chǎng)景模型Fig.5 The two-dimensional model of V-shaped targets

圖6 不同介電常數(shù)估計(jì)值V形目標(biāo)自聚焦成像結(jié)果Fig.6 The autofocusing imaging results of V-shaped targets for different estimated relative dielectric constant

通過(guò)引入時(shí)間相位因子,所提自聚焦處理算法僅能補(bǔ)償由介質(zhì)介電常數(shù)估計(jì)偏差引起的線性相位誤差,而對(duì)于探地雷達(dá)在空間移動(dòng)引起的二次相位誤差并不能完全補(bǔ)償。介質(zhì)介電常數(shù)誤差越大,相位因子的補(bǔ)償效果會(huì)相應(yīng)下降,所提自聚焦處理方法的有效性也會(huì)隨之降低。因此,使用所提算法對(duì)探地雷達(dá)回波數(shù)據(jù)進(jìn)行成像處理時(shí),應(yīng)盡量保證介質(zhì)介電常數(shù)估計(jì)誤差在一定合理范圍內(nèi)。

4 結(jié)論

提出了一種基于波場(chǎng)互相關(guān)的探地雷達(dá)快速自聚焦成像算法。相比于傳統(tǒng)時(shí)域逆時(shí)偏移算法,所提算法通過(guò)引入水平分層介質(zhì)頻域格林函數(shù),減小了對(duì)計(jì)算和存儲(chǔ)資源的需求,通過(guò)快速傅里葉變換大幅提升了算法的計(jì)算效率,因此能夠快速對(duì)地下目標(biāo)進(jìn)行高分辨成像。針對(duì)介質(zhì)參數(shù)未知而導(dǎo)致的圖像散焦問(wèn)題,通過(guò)引入時(shí)間相位因子可以得到不同聚焦時(shí)刻的圖像,然后基于圖像熵最小準(zhǔn)則得到最優(yōu)成像時(shí)刻及其對(duì)應(yīng)的自聚焦圖像。仿真結(jié)果表明當(dāng)介質(zhì)相對(duì)介電常數(shù)存在一定估計(jì)誤差時(shí),所提算法仍能夠得到目標(biāo)的高分辨率自聚焦成像結(jié)果。