周麗軍

(山西交通科學(xué)研究院集團(tuán)有限公司 智能裝備研究院，太原 030032)

目前對(duì)公路下隱藏病害的無(wú)損探測(cè)主要依賴于探地雷達(dá)(GPR)技術(shù)，即向地下空間發(fā)射高頻電磁波，利用地下物體與環(huán)境之間的介電差異對(duì)接收回波進(jìn)行處理，提取地下介質(zhì)分布情況[1-2]。應(yīng)用最廣泛的GPR探測(cè)技術(shù)主要是基于B-Scan掃描，通常得到的掃描圖像是基于同一場(chǎng)景的多方位多角度測(cè)量，并根據(jù)周圍介電環(huán)境重構(gòu)目標(biāo)像從而得到目標(biāo)的相對(duì)位置，但是并不能通過(guò)信號(hào)幅值或雙曲線峰值區(qū)分不同材質(zhì)的目標(biāo)[3]。而且僅通過(guò)B-scan圖像很難對(duì)目標(biāo)進(jìn)行精確的解析，識(shí)別難度大，主觀因素高[4]。

近年來(lái)，有學(xué)者將人工智能技術(shù)用在探地雷達(dá)公路病害檢測(cè)中。文獻(xiàn)[5]建立人工神經(jīng)網(wǎng)絡(luò)多層模型，用于道路病害的自動(dòng)診斷。文獻(xiàn)[6]利用深度字典學(xué)習(xí)方法處理探地雷達(dá)回波信號(hào)，區(qū)分多個(gè)不同形狀的掩埋目標(biāo)。文獻(xiàn)[7]將向量量化神經(jīng)網(wǎng)絡(luò)方法用于鐵路路基的病害識(shí)別中。文獻(xiàn)[8]利用聚類的方法進(jìn)行雷達(dá)圖像異常邊緣的自動(dòng)識(shí)別。這些基于B-scan圖像的方法還主要集中在物理模型中，多數(shù)是識(shí)別雷達(dá)圖像的雙曲線特征，對(duì)目標(biāo)的介電屬性識(shí)別方面的進(jìn)展還十分有限[9]。文獻(xiàn)[10]通過(guò)不同核匹配追蹤算法能有效識(shí)別水害和空洞，可識(shí)別城市道路土基病害。

隨著時(shí)域脈沖雷達(dá)的發(fā)展，基于雷達(dá)回波瞬態(tài)響應(yīng)的目標(biāo)識(shí)別方法得到了較大的發(fā)展，相比于成像技術(shù)，其優(yōu)勢(shì)在于僅需要目標(biāo)時(shí)域測(cè)量的單個(gè)瞬態(tài)響應(yīng)。根據(jù)SEM技術(shù)，目標(biāo)依賴于其瞬態(tài)響應(yīng)回波晚時(shí)部分的復(fù)自然諧振(CNR)，即極點(diǎn)特征，其分布規(guī)律只與目標(biāo)屬性有關(guān)，而與入射波的入射方式等無(wú)關(guān)，此性能能夠很好的用于目標(biāo)識(shí)別領(lǐng)域[11]。如文獻(xiàn)[12]利用極點(diǎn)差異識(shí)別不同形狀或不同方位的目標(biāo)。近年來(lái)，這些技術(shù)也已經(jīng)應(yīng)用于GPR探測(cè)的地下目標(biāo)識(shí)別。利用瞬態(tài)諧振現(xiàn)象估計(jì)與識(shí)別地下目標(biāo)，往往需要得到目標(biāo)的沖激響應(yīng)。文獻(xiàn)[13]從時(shí)域積分方程的角度研究了目標(biāo)的脈沖響應(yīng)以及其所包含的內(nèi)部諧振電流，并指出目標(biāo)的瞬態(tài)散射特性可以完全由其脈沖響應(yīng)表示。文獻(xiàn)[11]利用解卷積方法獲取目標(biāo)的沖激響應(yīng)，并對(duì)沖激響應(yīng)進(jìn)行奇異值分解，獲取目標(biāo)晚時(shí)響應(yīng)的CNR，研究半空間介電目標(biāo)的諧振行為。文獻(xiàn)[14]與文獻(xiàn)[15]基于瞬態(tài)響應(yīng)的內(nèi)、外部諧振分別研究了CNR與目標(biāo)幾何屬性的關(guān)系。但是上述文獻(xiàn)尚缺乏諧振與目標(biāo)介電特征的深入研究。

針對(duì)上述問(wèn)題提出一直基于極點(diǎn)特征聚類的地質(zhì)雷達(dá)公路隱藏裂縫自動(dòng)識(shí)別方法。通過(guò)對(duì)回波信號(hào)進(jìn)行奇異值分解，提取體現(xiàn)目標(biāo)信號(hào)的大奇異值，重構(gòu)回波信號(hào)，由此不僅能夠?qū)夭ㄐ盘?hào)進(jìn)行降維，還去除了噪聲干擾。再利用解卷積過(guò)程計(jì)算目標(biāo)的沖激響應(yīng)，利用奇點(diǎn)展開法(SEM)提取其晚時(shí)響應(yīng)部分的極點(diǎn)信息，構(gòu)建極點(diǎn)特征空間，并用模糊C均值聚類方法對(duì)極點(diǎn)進(jìn)行分類，根據(jù)已知環(huán)境樣本對(duì)極點(diǎn)進(jìn)行初始聚類中心優(yōu)化并取其均值，再利用此均值初始聚類中心對(duì)待識(shí)別目標(biāo)回波的極點(diǎn)特征空間進(jìn)行聚類，輸出分類結(jié)果。

1 算法描述

1.1 奇異值分解的降維方法

由于采集的數(shù)據(jù)量大，并且伴有很多噪聲，對(duì)回波數(shù)據(jù)進(jìn)行奇異值分解，獲得前M個(gè)主奇異值，用此重構(gòu)原始信號(hào)。方法如下：

對(duì)噪聲數(shù)據(jù)y(kΔt)，其中Δt是采樣間隔，k為采樣點(diǎn)數(shù)，構(gòu)造Hankel矩陣Y：

Y=

(1)

式中:N為數(shù)據(jù)采樣點(diǎn)數(shù)，L為函數(shù)束參數(shù)。

對(duì)Y進(jìn)行奇異值分解以獲得特征矢量及特征值：

Y=U∑VH

(2)

式中：矩陣U、V分別為Y的左奇異矩陣和右奇異矩陣，由矩陣YYH、YHY的特征矢量組成，上標(biāo)H表示共軛轉(zhuǎn)置。對(duì)角矩陣∑由Y的特征值組成:

(3)

由于σc(c=1，2，…，N-L)按照從大到小順序排列，比較每個(gè)σc和最大奇異值的比值直至滿足σc/σmax?10-p(p代表精度，默認(rèn)值為采樣數(shù)據(jù)的小數(shù)位數(shù))，即可確定M的取值。取∑的前M列∑′，重新構(gòu)造數(shù)據(jù)，顯然重構(gòu)數(shù)據(jù)維數(shù)遠(yuǎn)低于原始數(shù)據(jù)。并且由此得到的數(shù)據(jù)去除了噪聲干擾，突出了目標(biāo)信號(hào)信息。

1.2 目標(biāo)沖激響應(yīng)極點(diǎn)提取與極點(diǎn)特征空間構(gòu)成

入射波從目標(biāo)體表面返回到接收機(jī)而不再作用到目標(biāo)上時(shí)，目標(biāo)體上殘留的被激發(fā)(或誘導(dǎo))的電磁波會(huì)產(chǎn)生與目標(biāo)介電性能相關(guān)的諧振。理論上，此諧振電磁波在方位上獨(dú)立，僅依賴于目標(biāo)的形狀與介電特性，被稱為晚時(shí)響應(yīng)。為了得到目標(biāo)的諧振信息，需要將目標(biāo)的瞬態(tài)響應(yīng)從系統(tǒng)響應(yīng)中提取出來(lái)。為了減少數(shù)據(jù)量，去除噪聲，采用各回波信號(hào)經(jīng)過(guò)奇異值分解后的重構(gòu)信號(hào)。

首先測(cè)量背景信號(hào)，即沒(méi)有目標(biāo)時(shí)的回波響信號(hào)yb(t)，再測(cè)量有目標(biāo)時(shí)的目標(biāo)回波響應(yīng)yt(t)，設(shè)目標(biāo)的沖激響應(yīng)為h(t)，則卷積過(guò)程可以寫成：

yb(t)*h(t)=yt(t)

(4)

于是可以通過(guò)解卷積過(guò)程得到目標(biāo)的沖激響應(yīng)h(t).

此時(shí)目標(biāo)的沖激響應(yīng)h(t)可以寫成如下形式：

h(t)=hop(t)+u(t-t0)∑Riesit

(5)

其中hop(t)表示沖激響應(yīng)的早時(shí)響應(yīng)，上式右邊第二部分則是沖激響應(yīng)的晚時(shí)響應(yīng)，與目標(biāo)的固有屬性有關(guān)，其中t0為晚時(shí)響應(yīng)的開始時(shí)刻，其計(jì)算方法可參考文獻(xiàn)[16]。式中Ri表示第i個(gè)諧振態(tài)的復(fù)幅值，si=αi+j2πfi表示第i個(gè)諧振態(tài)的極點(diǎn)，其中αi與fi分別為衰減因子和諧振頻率。

根據(jù)獲得的衰減因子αi與諧振頻率fi構(gòu)成極點(diǎn)特征空間如下：

[A，F(xiàn)]={(α1，f1)，(α2，f2)，…，(αM，fM)}

(6)

其中M為重構(gòu)數(shù)據(jù)獲得的可用極點(diǎn)個(gè)數(shù)。將前M個(gè)諧振態(tài)的極點(diǎn)分布在式(6)的二維極點(diǎn)特征空間，能得到不同目標(biāo)沖激響應(yīng)的特征分布，從而通過(guò)極點(diǎn)特征進(jìn)一步識(shí)別不同目標(biāo)。

1.3 模糊C均值分類算法

模糊C均值(FCM)聚類算法是一種無(wú)監(jiān)督方法，它是通過(guò)對(duì)目標(biāo)函數(shù)進(jìn)行迭代優(yōu)化獲得各個(gè)樣本點(diǎn)對(duì)所屬的聚類中心的隸屬度，從而對(duì)樣本進(jìn)行劃分分類[17]。

(7)

(8)

由于模糊C均值聚類算法依賴于聚類中心的初始化值，本文根據(jù)多組已知的背景、空氣裂縫、充水裂縫的極點(diǎn)特征空間分布調(diào)整優(yōu)化初始聚類中心，使其趨于穩(wěn)定。

1.4 算法流程

第一步：降維。為了減少龐大的數(shù)據(jù)量，利用奇異值分解方法對(duì)對(duì)探地雷達(dá)的背景回波、空氣裂縫回波、充水裂縫回波信號(hào)進(jìn)行數(shù)據(jù)降維，同時(shí)也能對(duì)采集信號(hào)進(jìn)行降噪。

第二步：構(gòu)建特征空間。對(duì)回波信號(hào)提取沖激響應(yīng)，計(jì)算極點(diǎn)特征，并構(gòu)建以衰減因子與諧振頻率為極點(diǎn)的特征空間模型，分別計(jì)算出背景、空氣裂縫、沖水裂縫的極點(diǎn)特征空間。

第三步：對(duì)背景、空氣裂縫、沖水裂縫的極點(diǎn)特征空間中的衰減因子-諧振頻率對(duì)進(jìn)行模糊C均值聚類，根據(jù)已知的極點(diǎn)特征空間與聚類結(jié)果調(diào)整聚類算法的初始聚類中心。

第四步：對(duì)投影到極點(diǎn)特征空間的待識(shí)別信號(hào)極點(diǎn)根據(jù)調(diào)整優(yōu)化后的初始聚類中心進(jìn)行聚類，輸出分類結(jié)果。

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

本文使用GSSI公司的SIR-20系列天線對(duì)地下空氣裂縫，充水裂縫進(jìn)行實(shí)驗(yàn)分析，為了模擬空氣裂縫與充水裂縫，分別取空PVC管與充滿水的PVC管進(jìn)行試驗(yàn)，管長(zhǎng)1 m，直徑4 cm.天線的中心頻率為1 GHz，維度為49.5 cm×21 cm×55.6 cm.由于瞬態(tài)響應(yīng)是基于平面波信號(hào)源，為了將發(fā)射信號(hào)近似成平面波，預(yù)先對(duì)天線放置高度進(jìn)行經(jīng)驗(yàn)值估算，如圖1所示，在天線架高77.9 cm時(shí)，可以近似形成面積范圍139 cm×111 cm的平面波區(qū)域，即目標(biāo)在此區(qū)域內(nèi)其入射波可近似為平面波，后向接收回波為瞬態(tài)響應(yīng)回波。實(shí)測(cè)場(chǎng)景如圖2所示。

圖1 地質(zhì)雷達(dá)入射波近似平面波示意圖Fig.1 Gram of approximated plane wave from GPR

圖2 實(shí)測(cè)場(chǎng)景Fig.2 Real scene

分別獲取30組背景雷達(dá)回波信號(hào)、50組空氣裂縫雷達(dá)回波信號(hào)、50組充水裂縫雷達(dá)回波信號(hào)作為極點(diǎn)特征空間訓(xùn)練樣本，再各取5組空氣裂縫與充水裂縫回波作為測(cè)試樣本。

首先對(duì)回波數(shù)據(jù)進(jìn)行降維，既減少了數(shù)據(jù)量，又進(jìn)一步去除了噪聲干擾。在奇異值分解的降維方法中，選取合適的特征值重構(gòu)原始信號(hào)，這里精度取0.01，使得重構(gòu)信號(hào)既能完整的表征目標(biāo)信息，又能去除冗余與噪聲，50組空氣裂縫雷達(dá)回波樣本與50組充水裂縫雷達(dá)回波樣本的主要特征值個(gè)數(shù)及重構(gòu)信號(hào)與原始信號(hào)的誤差如圖3所示。

圖3 重構(gòu)信號(hào)與原始信號(hào)的誤差Fig.3 Error between reconstructed signal and original signal

從圖中可以看出，空氣裂縫與充水裂縫的重構(gòu)信號(hào)與原始信號(hào)誤差都在0.01以內(nèi)，在此誤差范圍內(nèi)，充水裂縫重構(gòu)信號(hào)所需的主特征值數(shù)量要大于空氣裂縫重構(gòu)信號(hào)的主特征值數(shù)量。而相比于分解得到的1 024個(gè)特征值，已經(jīng)極大減少了計(jì)算維數(shù)。

從圖4可以看出，比較空氣裂縫與充水裂縫的沖激響應(yīng)函數(shù)，在早時(shí)響應(yīng)(<3.2 ns)兩者的沖激響應(yīng)函數(shù)相似，而在5.38 ns處，充水裂縫的沖激響應(yīng)回波明顯平緩于空氣裂縫，即在此處充水裂縫回波體現(xiàn)的頻率低于空氣裂縫。

圖4 空氣裂縫與充水裂縫的沖激響應(yīng)函數(shù)Fig.4 Impulse response functions of air crack and water-filled crack

根據(jù)SEM算法對(duì)沖激響應(yīng)計(jì)算復(fù)自然諧振信息，即極點(diǎn)特征，得到極點(diǎn)的衰減因子與諧振頻率，如圖5所示。

圖5 探地雷達(dá)回波的極點(diǎn)特征分布Fig.5 Pole feature distribution of GPR echoes

從圖中可以看出，對(duì)沙地背景進(jìn)行探測(cè)，獲得的極點(diǎn)如“□”所示，主要分布在0.2 GHz，0.8～1 GHz，1.3 GHz附近，并且衰減因子的絕對(duì)值比較小，靠近0.空氣裂縫的極點(diǎn)分布如“○”所示，主要分布在1 GHz以上，即高頻分量較多。而充水裂縫的極點(diǎn)分布如“*”所示，其中一部分與空氣裂縫中的高頻分量類似，而另一部分則出現(xiàn)在0.3～1 GHz之間，相比于背景，出現(xiàn)了較多的高頻成分，相比于空氣裂縫，出現(xiàn)了較多的低頻成分。這一部分極點(diǎn)正是區(qū)分空氣裂縫與充水裂縫的關(guān)鍵特征信息。

對(duì)上述極點(diǎn)特征空間進(jìn)行模糊C均值聚類，得到的聚類結(jié)果如圖6所示。

從圖6中可以看出，利用模糊C均值聚類方法將極點(diǎn)特征空間分成4類，同時(shí)根據(jù)背景、空氣裂縫、充水裂縫的已知經(jīng)驗(yàn)信息調(diào)整初始聚類中心，每組樣本調(diào)整優(yōu)化后的初始聚類中心與其均值如圖箭頭處所示，對(duì)測(cè)試樣本的極點(diǎn)特征利用圖6的初始聚類中心均值作為其初始聚類中心，獲得的分類結(jié)果可以按照如下規(guī)則進(jìn)行判斷：

圖6 極點(diǎn)特征聚類Fig.6 Clustering of pole feature

對(duì)于極點(diǎn)p，如果滿足：

p∈①，p∈②，且p?③，則屬于空氣裂縫；

p∈③，則屬于充水裂縫。

分別取5組空氣裂縫與充水裂縫雷達(dá)回波信號(hào)的測(cè)試樣本，按同上方法計(jì)算其極點(diǎn)，得到的極點(diǎn)在極點(diǎn)特征空間中的分布如圖7所示。

圖7 測(cè)試樣本的極點(diǎn)在極點(diǎn)特征空間中的分布Fig.7 Distribution of poles from test samples in pole feature space

圖7中“+”的極點(diǎn)主要分布在①與②區(qū)域，并且沒(méi)有極點(diǎn)分布在③區(qū)域，根據(jù)上述規(guī)則可判定為空氣裂縫，而“△”的極點(diǎn)大多分布在①與③區(qū)域，根據(jù)上述規(guī)則，凡是出現(xiàn)在③區(qū)域的極點(diǎn)，可判定為充水裂縫的極點(diǎn)。結(jié)合給定的樣本參數(shù)，此判定均正確。

3 小結(jié)

針對(duì)地質(zhì)雷達(dá)無(wú)損檢測(cè)中回波數(shù)據(jù)量大，人工解讀困難等問(wèn)題，本文提出基于極點(diǎn)特征聚類的公路隱藏裂縫自動(dòng)識(shí)別方法。首先對(duì)回波信號(hào)進(jìn)行奇異值分解，提取體現(xiàn)目標(biāo)信號(hào)的大奇異值，重構(gòu)回波信號(hào)，由此不僅對(duì)回波信號(hào)進(jìn)行降維，還去除了噪聲干擾。在精度范圍允許值內(nèi)，分析了不同介電目標(biāo)所需的大奇異值數(shù)量。再通過(guò)解卷積過(guò)程獲得目標(biāo)的沖激響應(yīng)，利用奇點(diǎn)展開法提取其晚時(shí)響應(yīng)部分的極點(diǎn)信息，構(gòu)建極點(diǎn)特征空間，并對(duì)極點(diǎn)特征進(jìn)行模糊C均值聚類，為了降低初始聚類中心對(duì)聚類結(jié)果的影響，對(duì)已知環(huán)境下的學(xué)習(xí)樣本進(jìn)行聚類，并根據(jù)聚類結(jié)果調(diào)整優(yōu)化初始聚類中心，獲取調(diào)整后的均值，將待識(shí)別目標(biāo)回波的極點(diǎn)投影到極點(diǎn)特征空間按照均值初始化聚類中心進(jìn)行分類。通過(guò)100組學(xué)習(xí)樣本與10組測(cè)試樣本的實(shí)驗(yàn)，驗(yàn)證了算法的有效性。