(桂林電子科技大學(xué)信息與通信學(xué)院,廣西桂林541004)
在過去數(shù)十年中,隨著應(yīng)用于國防、救災(zāi)、醫(yī)療等領(lǐng)域的設(shè)備與相關(guān)技術(shù)的需求不斷增強,實用性探測雷達成像技術(shù)得到一定的發(fā)展,一些針對回波數(shù)據(jù)逆問題求解的成像算法也應(yīng)運而生。其中一類利用時域數(shù)據(jù)的重構(gòu)算法,例如時延求和、時間反轉(zhuǎn)鏡等算法,由于具有較為直觀的可視性、原理簡單等特點,因此具有一定的研究價值,但是此類方法通常存在重復(fù)運算較多、分辨率較差、無法反映目標外形、重構(gòu)結(jié)果存在多重解、數(shù)據(jù)容易受到復(fù)雜場景結(jié)構(gòu)布局和雜波干擾等弊端。而在另一類利用頻域數(shù)據(jù)的算法中,一種由Colton等提出的頻域線性采樣法(LSM)[1]得到了廣泛的關(guān)注,該方法通過構(gòu)造一個可以區(qū)分采樣點是否屬于散射體內(nèi)部的指示函數(shù),對回波數(shù)據(jù)進行逆問題求解,該方法可以在先驗信息較少的情況下對散射目標進行成像重構(gòu),因此具有很高的實用價值。隨后一系列基于該算法的改進方法得到了研究和驗證。例如,Nakamura等利用LSM求解熱導(dǎo)體邊界逆問題[2];Audibert等提出一種利用LSM結(jié)合部分確知散射體外形特征的數(shù)學(xué)模型求解逆問題的方法[3];Bazán等提出一種改進的最大間距準則并結(jié)合LSM方法求解逆問題[4];Erramshetty等提出一種利用LSM求解混合邊界特性散射體外形的方案[5]。盡管頻域線性采樣法有較好的邊界重構(gòu)效果,但相較于時域數(shù)據(jù),頻域數(shù)據(jù)運算通常只能針對單一頻點或多頻點分別運算后疊加,而待測場景整體信息分布在各個頻點中,單獨某個頻點數(shù)據(jù)包含的信息分量較少,因此需要大量觀測角度下的數(shù)據(jù)作為支撐。
針對傳統(tǒng)時域數(shù)據(jù)的特點和頻域方法的優(yōu)勢,已經(jīng)出現(xiàn)了一些時域線性采樣法(TD-LSM)求解逆問題的研究。例如,Aramini等提出一種利用單一功能化函數(shù)的傅里葉變換替代采樣點的離散遠場檢驗函數(shù)[6]。但是此方法本質(zhì)是將時域信號轉(zhuǎn)換為頻域信號,利用主頻點及其附近頻點的頻域數(shù)據(jù)進行單頻點或多頻點數(shù)據(jù)運算后疊加。Carayol等提出一種直接利用時域數(shù)據(jù)的TD-LSM方法[7]。Guo等給出了TD-LSM和多頻點LSM方法在小孔徑雷達、有限觀測角度和不同數(shù)據(jù)支撐條件下的性能比較[8]。Khaji等提出一種在已知彈性邊界數(shù)據(jù)條件下利用TD-LSM探測隱蔽空洞的方案[9]。Haddar等將時域線性采樣擴展到求解滿足羅賓紐曼邊界條件的散射體問題上[10]。Chen等針對對稱外形散射體提出一種半平面TD-LSM方法[11]。盡管TD-LSM數(shù)據(jù)中包含了較多的分量信息,且一定程度上避免了需要大量觀測角度下的數(shù)據(jù)支撐。但由于TD-LSM需要將原始數(shù)據(jù)變換為離散卷積形式,導(dǎo)致運算消耗嚴重增加,且計算檢驗函數(shù)過程中的輻射時延沒有特定的選取標準。
為了克服以上傳統(tǒng)TD-LSM的缺陷,本文主要通過分析時域數(shù)據(jù)波場傳播特性,利用信號互相關(guān)特性,首先提出并驗證了最佳時延值選取準則,并進一步在基于最佳時延選取準則的基礎(chǔ)上,提出一種基于分割思想的TD-LSM改進算法,該方法先利用各個發(fā)射陣元的零時信號與接收信號進行互相關(guān)運算,預(yù)先獲得各道接收信號的時延值,再利用這些時延值求出該發(fā)射陣元到達散射體的平均時延,作為等式右側(cè)最佳時延值,以此替代經(jīng)驗選取時延方法。隨后,在每個采樣點運算過程中根據(jù)該點到各個陣元走時范圍,將等式左右兩端的數(shù)據(jù)矩陣和計算點源輻射函數(shù)分割,截取其中有效數(shù)據(jù)后,再利用TD-LSM進行運算,以此提高運算效率。
圖1為TD-LSM應(yīng)用場景示意圖,其中Ω∈?d,d=2或d=3表示封閉待測區(qū)域,區(qū)域中存在位置固定的散射體D,其邊界?D連續(xù),散射體外部記為電磁波在背景介質(zhì)中速度穩(wěn)定為c0,待測區(qū)域被位于曲線內(nèi)的N個收發(fā)共置陣元測量,采用信息分量較為豐富的時域雷克波信號作為信源,其表達式如下:

式中,ζ=2π2f2,χ=1/f,f表示主頻,x,y∈?d,t∈[0,T]表示時間信號的長度,間隔為Δt,時間序列量化個數(shù)為N=T/Δt,x=[x1,x2]∈Ω表示采樣點在待測區(qū)域中的坐標,y=[y1,y2]∈Γi表示入射陣元位置,根據(jù)時間信號的有序性,當(dāng)t<0時,ui(t;x,y)≡0。根據(jù)惠更斯原理,對于信源照射到的待測區(qū)域中任意采樣點x,若x存在散射系數(shù),則可以將x看作是一個新的子波源,若τ表示子波源相對發(fā)射信號的時延,則反射波可以改寫為

由此得到電磁波在自由空間傳播過程中的波動方程、初始條件和邊界條件:

式中,τ為信號時延,c0為自由空間波速,δ為狄拉克函數(shù),ρ(t)為單位階躍函數(shù),ui(·;x,y)為雷克波信號源激勵,且在遠場模式下,待測區(qū)域內(nèi)的總場u(t;x,y)、入射場ui(t;x,y)以及散射場us(t;x,y)滿足下式:

通過求解式(3)可以得到散射場us(·;x,y)的基本解:


圖1 場景示意圖
因此,通過正演過程求解得到的us(t;x,y)即為通過數(shù)學(xué)表達獲取的仿真實驗數(shù)據(jù),在實測實驗中us(t;x,y)即為雷達天線陣列接收到的數(shù)據(jù)。利用TD-LSM求解逆問題,其本質(zhì)就是在已知陣元位置的條件下,利用在不同入射角y n下的采集回波數(shù)據(jù)us(x,t;y n),對每個采樣點求解式(6),從而在位置目標物理類型和幾何形狀未知的條件下,獲取未知目標(若存在)的位置與邊界外形信息:

由于時域數(shù)據(jù)具有因果時變特性,因此需要將左端數(shù)據(jù)矩陣FTD修改為離散卷積形式,若原始場景中包含M個收發(fā)共置陣元,每個陣元接收時間為t∈[0,T]內(nèi)的數(shù)據(jù),數(shù)據(jù)時間間隔為Δt,單個陣元的數(shù)據(jù)總長度為N=T/Δt,則得到FTD離散化形式:

再將離散形式數(shù)據(jù)矩陣重組為離散卷積形式后得到

式中,y i表示第i個發(fā)射天線,x r表示第r個接收天線。若z為采樣點坐標,τ為選取的時延,則利用所有入射方向y n∈Γi的接收數(shù)據(jù)對所有采樣點z∈?構(gòu)建TD-LSM核心方程組:

式中,ψ∞(x-z,t-τ)為遠場檢驗函數(shù),其數(shù)學(xué)模型與信源發(fā)射信號近似,通過求解式(9),并根據(jù)如下判斷條件可以得到關(guān)于采樣點z的相應(yīng)結(jié)論:
1)若z∈D,且εr>0,則存在函數(shù)


然而,由于求解式(9)為病態(tài)逆問題,因此需要引入Tikhonov正則化方法,利用F′TD的SVD分解參數(shù){Utd,Std,Vtd}求出式(9)的可靠解g(z,τ),并求出其L2范數(shù)G(z,τ)作為像素點的指示器函數(shù):

至此,直接利用時域數(shù)據(jù)對待測區(qū)域進行重構(gòu)的原始TD-LSM主要步驟如下:
1)將式(9)中原始數(shù)據(jù)矩陣變換為離散卷積形式F′TD,并對其進行SVD分解得到{Utd,Std,Vtd};
2)將待重構(gòu)場景進行離散化處理,選取某個特定時刻作為等式右端時延值,從第一個采樣點開始選取,計算該時延下式(9)中的檢驗函數(shù),并進行離散卷積化;
3)根據(jù)式(12)求出該采樣點指示器G(z,τ),若G(z,τ)→0則z∈D,若G(z,τ)→ ∞ 則z∈
4)計算下一個采樣點z n(n=1,2,3,…,N),重復(fù)步驟2)~3),直到遍歷成像區(qū)所有采樣點。
盡管TD-LSM克服了頻域線性采樣法需要大量數(shù)據(jù)點支撐的缺陷,同時避免了對時域數(shù)據(jù)進行傅里葉變換的操作,但傳統(tǒng)TD-LSM同樣存在缺陷,時延τ的選取沒有理論化標準,通常依靠經(jīng)驗選取,因此,需要尋找一種可以替代經(jīng)驗選取時延值的方法;其次,由于等式左端數(shù)據(jù)矩陣變換為離散卷積形式,導(dǎo)致數(shù)據(jù)量異常增大,SVD運算異常緩慢。
針對時延τ的選取沒有理論化標準,通常依靠經(jīng)驗選取的弊端,首先從吉洪諾夫原理入手,若獲得原始數(shù)據(jù)矩陣F′的SVD分解參數(shù){Utd,Std,Vtd}后,則原始數(shù)據(jù)可以表示為

式中,Utd={u1,u2,u3,…,u n},Vtd={v1;v2;v3;…;v n},Std=diag{s1,s2,s3,…,s n}。3個分解參數(shù)包含了整個原始矩陣的全部信息,其信息的重要程度按照奇異值矩陣Std依次遞減。根據(jù)SVD分解原理,若存在一個矩陣與原始數(shù)據(jù)矩陣相似,則可以通過左右奇異矢量與相似矩陣還原原始分解矩陣。而時域線性采樣法作用機理是在遍歷采樣點過程中選取一個特定σ時刻,計算該時刻下從每個采樣點接收到的時域波形矩陣,再利用吉洪諾夫正則化驗證每個采樣點是否存在散射系數(shù)。若檢驗函數(shù)與原始數(shù)據(jù)矩陣相似,則根據(jù)式(9)、式(12)求解得到

式中,f(x-z,t-τ)表示原始數(shù)據(jù)矩陣分解出的奇異值參數(shù)。因此,當(dāng)某個采樣點計算出的檢測函數(shù)與接收到信號相關(guān)程度最大時,可以通過式(12)判斷該點存在散射系數(shù)。因此在TD-LSM中,選取有效的檢驗函數(shù)時延σ是一個重要環(huán)節(jié),σ越接近原始數(shù)據(jù)接收回波中的時延值τ,其檢驗函數(shù)和接收數(shù)據(jù)函數(shù)相關(guān)度越高,其獲得的結(jié)果越好。由于墻后目標通常未知,因此在缺少先驗信息的條件下,無法有效地獲得各個接收陣元的τ值,導(dǎo)致經(jīng)驗選取σ值的過程中容易與τ值產(chǎn)生較大誤差。針對這一缺陷,接下來提出一種利用互相關(guān)原理尋找最佳檢驗函數(shù)時延σ值的方法,若u0(t)表示發(fā)射陣元零時信號,u i(t-τ)表示回波接收信號,τ為接收信號相對零時信號的時延,根據(jù)數(shù)據(jù)互相關(guān)原理得到

因為u0(t)與u i(t-τ)為相同信號,最佳時延準則是利用發(fā)射陣元零時信號與接收數(shù)據(jù)矩陣互相關(guān),得到接收回波的時延時刻,使得在τ時刻達到最大值1,隨后選取τ作為檢驗函數(shù)σ值。
通過正演實驗驗證最佳時延選取準則獲取的σ值,首先分析在TD-LSM測量過程中的信號傳播路徑,如圖2所示。從圖2可以看出,在路徑傳播的過程中,存在散射體上的某點,與位置固定的兩個陣元組成完整的信號路徑,若將該點A看作子波源,可以等效獲得A到達各個陣元的信號集合,對此信號集合補償由陣元D到達A點的距離時延,便可使該信號達到與陣列中部分陣元的接收信號完全相關(guān)狀態(tài)。

圖2 信號傳播路徑示意圖
根據(jù)3.1節(jié)仿真實驗場景相同配置,仿真反射波到達各個陣元的正演過程,選取金屬散射體前表面上正中心采樣點,坐標(0.55,0.4),并在此添加一個發(fā)射陣元,發(fā)射與仿真實驗相同的時域雷克信號波,其他陣元作為接收端,獲取此場景下的驗證接收數(shù)據(jù),如圖3所示。其中,圖3(a)中補償1曲線表示對驗證數(shù)據(jù)補償?shù)?1號陣元到達點(0.55,0.4)的時延后,與該組數(shù)據(jù)的互相關(guān)系數(shù)曲線;補償2曲線為補償?shù)?4號陣元到達該點的時延后,驗證接收數(shù)據(jù)與該組互相關(guān)系數(shù)曲線;未補償曲線表示不補償14號陣元到達該點的距離時延,驗證接收數(shù)據(jù)與陣列陣元的互相關(guān)系數(shù)。圖3(b)中3條曲線分別表示檢測函數(shù)在0.25 ns,0.75 ns和1.25 ns時延下,與第11組接收數(shù)據(jù)的互相關(guān)系數(shù)曲線,這組驗證實驗表明,0.75 ns為較為合適的時延值選擇。圖3(c)表示檢測函數(shù)在0.75 ns時延下與第5組、第11組、第16組接收數(shù)據(jù)的互相關(guān)系數(shù)曲線,這組驗證表明該點信號較為接近第5組陣元的回波散射點。通過正演實驗可以發(fā)現(xiàn),互相關(guān)函數(shù)可以尋找兩列信號最相似的時刻,利用這個特性,可以提高時延選取效率,驗證了最佳時延選取準則的有效性。
由于TD-LSM需要對采集數(shù)據(jù)進行離散卷積變化,導(dǎo)致參與運算的數(shù)據(jù)維數(shù)嚴重膨脹,針對這一缺陷,進一步提出一種基于數(shù)據(jù)分割思想的TDLSM改進。從圖4可以看出,TD-LSM成像過程中以采樣點作為循環(huán)變量。利用這一特性,在遍歷成像區(qū)域采樣點的過程中,首先計算各個天線到達該點的時延范圍τ1~τn,找出時延范圍中的最大值τmax和最小值τmin,隨后從離散卷積變換后的數(shù)據(jù)中截取位于[τmin,τmax]子區(qū)間內(nèi)的數(shù)據(jù)作為有效數(shù)據(jù),再對有效數(shù)據(jù)進行離散卷積化變形和運算,以此達到數(shù)據(jù)壓縮提高運算效率的目的。至此,基于分割思想且符合最佳時延選取準則的TD-LSM改進方法步驟為:
1)利用各個發(fā)射陣元的零時信號與各組接收信號進行互相關(guān)運算獲取接收時延的集合,計算時延集合的均值作為檢驗函數(shù)的時延值;
2)將待重構(gòu)場景進行離散化處理,從第一個采樣點依次遍歷選取,計算每個采樣點到達各個陣元的時間范圍作為該采樣點的有效數(shù)據(jù)區(qū)間,從原始數(shù)據(jù)中剔除落于有效區(qū)間外的冗余數(shù)據(jù),隨后將其離散化重組并進行SVD分解獲得相應(yīng)的分解參數(shù){U′td,S′td,V′td};
3)利用步驟1)獲取的時延值,計算該時延下采樣點到各個陣元的檢驗函數(shù),并重組為離散卷積形式;
5)對下一個采樣點z n(n=1,2,3,…,N)重復(fù)步驟2)~4),直至遍歷所有采樣點。
通過仿真和實測數(shù)據(jù)對本文提出的算法進行驗證,并采用如下幾個評判標準:利用平均運行時間(APRT)判定算法執(zhí)行時間,APRT越小代表算法運行時間越短;利用均方根誤差(RMSE)判定成像效果與原散射體位置偏離程度,RMSE值越小成像結(jié)果越接近原散射體位置;利用一階差分均值法(MFOD)判定輪廓特性,MFOD值越小檢測區(qū)域越平滑;利用圖像熵(ENT)判定成像復(fù)雜度,ENT值越小整幅圖雜波越少,其中RMSE和MFOD為無綱量參數(shù)。
首先仿真實驗信源采用主頻為1 GHz的雷克波,整體場景大小為邊長1.05 m的正方形,數(shù)據(jù)采集裝置采用20個收發(fā)共置陣元組成的線陣,陣元間隔0.05 m,距離墻體0.05 m,從x軸0.05 m處開始向x軸正方向延伸。墻體厚度0.1 m,墻體前表面距離邊界0.1 m,白色為墻體位置,黑色為目標金屬矩形物體,三角形表示各個陣元,左下角坐標位于(0.4,0.4),右上角坐標位于(0.7,0.45),如圖5(a)所示。
從圖5(b)~圖5(d)可以看出,TD-LSM選取時延τ值無論過大還是過小都會影響最終成像結(jié)果的精確度。根據(jù)本文提出的最佳時延選取準則,最佳的時延σ值范圍在真實目標到各個陣元的平均時延上下浮動,因此在這一準則下,最佳時延范圍應(yīng)在整個金屬前表面上各個采樣點到各個陣元的平均時延,直接根據(jù)仿真場景配置計算可得7.564 4 ns。仿真結(jié)果驗證了最佳時延選取準則以及選取時延的改進方法。


圖5 經(jīng)驗選取時延與本文提出的改進方案對比
表1給出了圖5中4組成像結(jié)果的APRT,RMSE,MFOD和ENT的比較參數(shù)。從表中可以看出,盡管第1組結(jié)果中的MFOD值最小,說明其成像結(jié)果輪廓最好,但是RMSE值和ENT值較大,說明其結(jié)果與散射體真實位置和形狀偏差較大且圖像中雜波較多;第2組結(jié)果中,由于選取了接近最佳選取準則計算得到的時延值,因此,其RMSE值相比第1組性能提高約30%,ENT值相比第1組性能提高約17%,結(jié)果說明其結(jié)果較為準確,可信度較高且雜波較少,由于成像結(jié)果較為接近散射體位置和形狀,因此導(dǎo)致MFOD值相對略有增加;第3組結(jié)果中,由于選取時延值再次違背最佳時延準則,盡管其ENT值為4組數(shù)據(jù)中最小,但其RMSE值和MFOD值均為4組結(jié)果中最差,其結(jié)果十分不可靠;最后一組結(jié)果中,相對于第2組數(shù)據(jù),盡管犧牲了小部分MFOD值和ENT值為代價,換取獲得RMSE值最優(yōu)和運算效率顯著提升,表明本文提出的算法不僅可靠性最佳且運算效率最大。

表1 幾種時延取值下運算效果對比
實測實驗選取法國菲涅爾實驗室提供的散射場數(shù)據(jù)來對本文提出的算法進行驗證。實驗?zāi)繕艘?guī)格為80 mm×50 mm,一側(cè)缺口為40 mm的金屬立方體,位于測試區(qū)域中心,發(fā)射天線位于半徑為0.72 m的圓陣上,入射角度范圍0°~350°(間隔10°),接收天線位于半徑為0.76 m的圓陣上,接收角度范圍60°~300°(間隔10°),頻率范圍為1~8 GHz,步進1 GHz。
從圖6(a)~圖6(d)的頻域成像結(jié)果可以看出,成像結(jié)果雜波和虛像較多,而低頻部分結(jié)果邊界連續(xù)性和平滑性較差。相較于單頻或多頻提取方法的結(jié)果,圖6(e)和圖6(f)中的傳統(tǒng)TD-LSM和改進TD-LSM直觀上可以較為良好的兼顧二者。


圖6 幾種不同算法的成像結(jié)果對比
表2給出了4組方法的APRT,RMSE,MFOD和ENT的比較結(jié)果。第1組結(jié)果中,盡管ENT值最佳,但是MFOD值較差;第2組結(jié)果中,相較第1組結(jié)果,僅僅提升了MFOD值,RMSE值和ENT值均下降,造成這一結(jié)果的主因,是該場景中U型金屬缺口部分僅有部分正對陣元可以提供真實的信息情況,而其他部分的陣元無法透過金屬外殼獲取U型缺口內(nèi)部信息;第3組結(jié)果中,運行時間大幅增加,但是RMSE值、MFOD值和ENT值均提升有限;而本文提出的改進算法在盡可能保證重構(gòu)質(zhì)量的情況下,最大幅度地減少運算時間。

表2 實測實驗算法對比
通過對傳統(tǒng)TD-LSM的數(shù)值公式分析,獲得線性采樣類方法的核心步驟為利用正則化方法求解病態(tài)解逆問題。隨后,分析了TD-LSM算法特點和數(shù)據(jù)結(jié)構(gòu)特點,找到了檢驗函數(shù)時延值的選取是影響成像效果的重要參數(shù)。針對檢驗函數(shù)時延值的選取沒有標準化的流程和方法,結(jié)合信號互相關(guān)原理,提出一種TD-LSM最佳時延選取準則,并通過正演實驗進行分析驗證。之后進一步通過分析信源在整個待測區(qū)域的走時過程,進一步提出一種基于分割思想的TD-LSM改進,利用采樣點到各個陣元的距離范圍,有效地剔除原始數(shù)據(jù)矩陣中的冗余信息,大大提高了運算效率。
[1]COLTON D,MONK P.Linear Sampling Method[M].US:Society for Industrial&Applied Mathematics,1999.
[2]NAKAMURA G,WANG H.Linear Sampling Method for the Heat Equation with Inclusions[J].Inverse Problems,2013,29(29):104015.
[3]AUDIBERT L,HADDAR H.A Generalized Formulation of the Linear Sampling Method with Exact Characterization of Targets in Terms of Farfield Measurements[J].Inverse Problems,2014,30(3):035011.
[4]BAZáN F S V,FRANCISCO J B,LEEM K H,et al.Using the Linear Sampling Method and an Improved Maximum Product Criterion for the Solution of the Electromagnetic Inverse Medium Problem[J].Journal of Computational and Applied Mathematics,2015,273(C):61-75.
[5]ERRAMSHETTY M,BHATTACHARYA A.Shape Reconstruction of Mixed Boundary Objects by Linear Sampling Method[J].IEEE Trans on Antennas and Propagation,2015,63(7):3077-3086.
[6]ARAMINI R,BRIGNONE M,PIANA M.The Linear Sampling Method Without Sampling[J].Inverse Problems,2006,22(6):2237.
[7]CARAYOL Q,COLLINO F.A Sampling Method for Inverse Scattering in the Time Domain[J].ESAIM:Mathematical Modelling and Numerical Analysis,2004,38(2):371-394.
[8]GUO Y,MONK P,COLTON D.The Linear Sampling Method for Sparse Small Aperture Data[J].Applicable Analysis,2015,95(8):1599-1615.
[9]KHAJI N,MANSHADI S H D.Time Domain Linear Sampling Method for Qualitative Identification of Buried Cavities from Elastodynamic over-Determined Boundary Data[J].Computers&Structures,2015,153(2):36-48.
[10]HADDAR H.An Improved Time Domain Linear Sampling Method for Robin and Neumann Obstacles[J].Applicable Analysis,2014,93(2):369-390.
[11]CHEN B,MA F,GUO Y.Time Domain Scattering and Inverse Scattering Problems in a Locally Perturbed Half-Plane[J].Applicable Analysis,2017,96(8):1303-1325.