奇異值分解法和模糊C均值法在探地雷達(dá)圖像處理中的應(yīng)用研究

聞鋮， 嚴(yán)玲玲， 牟京亞， 屈乾達(dá)， 熊宇， 陳秋航

(1.國(guó)網(wǎng)湖北省電力有限公司檢修公司， 湖北 武漢 430050;2.三峽大學(xué) 電氣與新能源學(xué)院, 湖北 宜昌 443000)

0 引言

探地雷達(dá)(Ground Penetrating Radar，GPR)作為一種無(wú)損的探測(cè)設(shè)備，目前在高速公路、鐵路路基檢測(cè)，考古發(fā)掘探測(cè)及隧道內(nèi)部探傷等各項(xiàng)工程探測(cè)都使用其作為主要探測(cè)手段[1]。探地雷達(dá)技術(shù)可以用來(lái)埋藏在地下的金屬和非金屬物體，并能夠進(jìn)行二維或三維模式下的定位。探地雷達(dá)與其他常規(guī)的地下探測(cè)方法相比，具有探測(cè)速度快、探測(cè)過(guò)程連續(xù)、分辨率高、操作方便靈活探測(cè)費(fèi)用低等優(yōu)點(diǎn)，故其在工程勘察領(lǐng)域的應(yīng)用日益廣泛[1]。

由于淺部地下介質(zhì)復(fù)雜多變，且大都具有頻散特性，電磁波在其中傳播衰減速度快，散射非常強(qiáng)烈，加之各種人為設(shè)施的雜亂回波、天線的耦合等原因，天然干擾很難與目標(biāo)體的反射波分離，使得在進(jìn)行實(shí)際探測(cè)中，得到的數(shù)據(jù)不可避免地包含了雜波和噪聲[2-3]。所以探地雷達(dá)接收到的數(shù)據(jù)包括目標(biāo)信號(hào)、雜波以及噪聲。在探地雷達(dá)的應(yīng)用過(guò)程中，經(jīng)典的雜波抑制方法是通過(guò)在頻域和時(shí)域中進(jìn)行數(shù)據(jù)濾波來(lái)實(shí)現(xiàn)的。但當(dāng)被探測(cè)物體接近于地表，或與地表具有相同的時(shí)間響應(yīng)的情況下，這種雜波抑制方法便無(wú)法奏效。文獻(xiàn)[4]提出了一種基于頻率域-多普勒域?qū)拵Ю走_(dá)雜波抑制方法，雖然該方法能夠有效減小濾除雜波過(guò)程中的能量損失問(wèn)題，但在探測(cè)目標(biāo)和雜波有較多重疊的區(qū)域時(shí)，該方法便無(wú)法保障目標(biāo)的探測(cè)精度。其他一些雜波抑制方法包括卡爾曼濾波法[5]和基于神經(jīng)網(wǎng)絡(luò)的雜波抑制方案[6]也用于減少探地雷達(dá)成像中雜波的減少。這些方案均是基于目標(biāo)的雙曲線特征，如果目標(biāo)被埋在地面附近，目標(biāo)的雙曲線特征將因?yàn)榈孛娣瓷涠丿B，影響對(duì)探測(cè)目標(biāo)的識(shí)別。

近年來(lái)，相關(guān)學(xué)者對(duì)基于奇異點(diǎn)分析、主成分分析、獨(dú)立成分分析的統(tǒng)計(jì)分析方法在雜波抑制領(lǐng)域的應(yīng)用有了廣泛的研究。這些方法將接收的數(shù)據(jù)分解成以下成分：目標(biāo)、雜波和噪聲。剔除雜波和噪聲后得到的成分即為包含目標(biāo)特征的成分。文獻(xiàn)[7]提出了一種基于主元分析主分量估計(jì)的方法，確定了目標(biāo)特征成分選取的限制條件。

針對(duì)上述問(wèn)題，本文提出了一種基于奇異值分解法(Singular Value Decomposition，SVD)和模糊C均值法(Fuzzy C-Means，F(xiàn)CM)的GPR圖像處理方法。使用改進(jìn)的SVD將接收的數(shù)據(jù)辨別成探測(cè)目標(biāo)、雜波和噪聲的子空間。針對(duì)探測(cè)目標(biāo)、雜波和噪聲信號(hào)存在重疊區(qū)域的情況時(shí)，使用FCM聚類(lèi)算法對(duì)其進(jìn)行分離，最后通過(guò)不同光譜成分的加權(quán)得到所需的目標(biāo)圖像。

1 基于奇異值分解和模糊C均值法的圖像處理方法

1.1 探地雷達(dá)探測(cè)信號(hào)的組成

探地雷達(dá)探測(cè)信號(hào)一般由雜波、探測(cè)目標(biāo)和噪聲信號(hào)組成，運(yùn)用探地雷達(dá)進(jìn)行地下物質(zhì)探測(cè)。其探測(cè)示意圖,如圖1所示。

圖1 探地雷達(dá)探測(cè)示意圖

由此可見(jiàn)在使用探地雷達(dá)進(jìn)行地下物質(zhì)探測(cè)時(shí)，接收的探測(cè)信號(hào)中，地層反射波(雜波)信號(hào)和噪聲信號(hào)是不可避免的。

1.2 奇異值分解(SVD)定義

矩陣的奇異值分解的主要理論定義如下。

1.3 模糊C均值法(FCM)定義

FCM算法是一種基于目標(biāo)函數(shù)的聚類(lèi)算法，通過(guò)目標(biāo)函數(shù)的最小化實(shí)現(xiàn)聚類(lèi)。一般使用X={x1(i,j),x2(i,j),…,xn(i,j)}來(lái)描述圖像，其中n為像素點(diǎn)總數(shù)，x(i,j)為第k個(gè)像素點(diǎn)的特征向量(i≤k≤j)利用聚類(lèi)方式進(jìn)行圖像分割的實(shí)質(zhì)是將像素點(diǎn)集合分為C類(lèi)的問(wèn)題，每個(gè)類(lèi)屬中均含有唯一的聚類(lèi)中心，在不斷的迭代中實(shí)現(xiàn)聚類(lèi)中心的更新，通過(guò)目標(biāo)函數(shù)的最小化實(shí)現(xiàn)聚類(lèi)結(jié)果的優(yōu)化，使用隸屬度矩陣來(lái)描述各個(gè)像素點(diǎn)的類(lèi)屬性質(zhì)，通過(guò)單個(gè)像素點(diǎn)與聚類(lèi)中心的相似程度來(lái)判斷其在不同類(lèi)別中的隸屬程度[8]。

1.4 基于SVD和FCM的圖像優(yōu)化方法

在GPR圖像增強(qiáng)算法中，對(duì)于掃描得到的圖像X(X為N行M列的數(shù)據(jù)矩陣，M代表采樣道數(shù)，N代表每道數(shù)據(jù)的樣點(diǎn)數(shù))，對(duì)其用SVD法將其分解成不同的光譜分量,如式(1)。

X=USVH

(1)

式中，U和V分別為有N行N列、M行M列的酉矩陣。S=diag(s1,s2…,sM)，且s1≥s2≥…≥sM>0，其均為X的奇異值。雜波圖像記為Xc，目標(biāo)圖像記為Xt，噪聲圖像記為Xn。探測(cè)圖像X的分解,如式(2)。

(2)

式中，k1為雜波的奇異值，k2-k1為探測(cè)目標(biāo)以及其他噪聲的奇異值。

文獻(xiàn)[3]提出了兩種提取目標(biāo)光譜圖像的方法，方法一提出，第一類(lèi)光譜分量包含了地層雜波，如式(3)。

(3)

而其他的光譜分量則包含了探測(cè)目標(biāo)和噪聲等，如式(4)。

(4)

方法一的局限性在于，只過(guò)濾掉了由于地表反射等產(chǎn)生的雜波，然而噪聲成分卻沒(méi)有分離，即方法一得到的圖像是包含噪聲的。

基于此，文獻(xiàn)[3]提出了方法二。在方法二中，第一類(lèi)光譜分量包含了雜波，第二類(lèi)光譜分量則包含了目標(biāo)圖像，噪聲則包含在其余的光譜分量中,如式(5)。

(5)

然而這種方法卻無(wú)法對(duì)多個(gè)目標(biāo)進(jìn)行識(shí)別。事實(shí)上，在探地雷達(dá)進(jìn)行探測(cè)的過(guò)程中，即使探測(cè)場(chǎng)景僅有一個(gè)探測(cè)目標(biāo)，該目標(biāo)子空間的維度也可以超過(guò)一維，基于此，本文提出一種改進(jìn)的SVD圖像分解方法,如式(6)。

(6)

式中,X0雜波剔除后的圖像，Y是目標(biāo)圖像，Z則是雜波圖像。矩陣Y的秩RY=k2

(7)

(8)

原始圖像X基于探測(cè)目標(biāo)、噪聲、雜波的子空間,如式(9)。

(9)

k2的值需要通過(guò)計(jì)算來(lái)確定。此時(shí)根據(jù)文獻(xiàn)[9]提出的模糊C均值法來(lái)計(jì)算k2準(zhǔn)確的的值。應(yīng)用模糊C均值法的損失函數(shù)Jc最小求解表達(dá)式,如式(10)。

(10)

(11)

式中，C1、C2分別是目標(biāo)和噪聲類(lèi)別的聚類(lèi)中心。聚類(lèi)中心Cl的計(jì)算方式,如式(12)。

(12)

最后目標(biāo)圖像Xt3由不同光譜成分的加權(quán),如式(13)。

(13)

2 實(shí)驗(yàn)驗(yàn)證

2.1 實(shí)驗(yàn)結(jié)果對(duì)比

本文采用WGPR200探地雷達(dá)系統(tǒng)進(jìn)行現(xiàn)場(chǎng)實(shí)驗(yàn)。選用900 MHz天線，512個(gè)采樣點(diǎn)，步長(zhǎng)為0.015 m。實(shí)驗(yàn)過(guò)程中的計(jì)算參數(shù)的選取參照文獻(xiàn)[10]，其設(shè)定值為：干砂介電常數(shù)為3，空氣介電常數(shù)為1，電導(dǎo)率為1. 0×10-10S/m。探測(cè)區(qū)域選定為一個(gè)2 m×4 m×1 m的布滿(mǎn)干沙的沙坑。將3個(gè)盛滿(mǎn)空氣，其大小分別為0.25 m×0.05 m×2 mm、0.20 m×0.05 m×2 mm、0.15 m×0.05 m×2 mm的小型塑料箱埋藏在沙坑中。

根據(jù)文獻(xiàn)[9]中提出的FCM聚類(lèi)算法，可得奇異值sm的幅值圖和基于FCM聚類(lèi)算法的目標(biāo)和噪聲的隸屬度圖像，如圖2(a)、(b)所示。

探測(cè)所得原始圖像轉(zhuǎn)換成灰度圖的結(jié)果,如圖3所示。

采用文獻(xiàn)[3]中的方法二和本文方法進(jìn)行了對(duì)比，對(duì)比實(shí)驗(yàn)結(jié)果,如圖4、圖5所示。

通過(guò)圖像對(duì)比可以看出，相較于文獻(xiàn)[3]的處理結(jié)果，本文所提出的算法能夠有效地區(qū)分探測(cè)目標(biāo)、雜波以及噪聲的重疊邊界，獲取更清晰準(zhǔn)確的目標(biāo)圖像。

2.2 性能參數(shù)對(duì)比

通過(guò)計(jì)算峰值信噪比(Peak Signal to Noise Ratio，PSNR)和均方誤差(Mean Square Error)來(lái)評(píng)判經(jīng)過(guò)算法處理后的圖像處理效果，PSNR是一種評(píng)價(jià)圖像的客觀標(biāo)準(zhǔn)，為了衡量經(jīng)過(guò)處理后的影像品質(zhì)，我們通常會(huì)參考PSNR值來(lái)衡量圖像處理方法能否令人滿(mǎn)意。MSR表示原圖像和處理圖像之間的均方誤差。PSNR值越高，MSR值越小，則說(shuō)明處理后的圖像和原圖越接近。本文通過(guò)Matlab2016a計(jì)算出文獻(xiàn)[3]中的方法二與本文所提出算法的PSNR值和MSR值,如表2所示。

圖2(a) 奇異值sm幅值圖

圖2(b) 目標(biāo)和噪聲的隸屬度函數(shù)

圖3 原始探測(cè)圖像灰度圖

圖4 文獻(xiàn)[3]中方法二處理結(jié)果

圖5 本文方法處理結(jié)果

表2 PSNR值與MSR值對(duì)比

通過(guò)以上結(jié)果，證明本文所提出的方法MSR值更低，PSNR值更高，較于文獻(xiàn)[3]中的方法，對(duì)圖像的處理效果更好。

3 總結(jié)

本文提出了一種基于奇異值分解法和模糊C均值法的探地雷達(dá)圖像處理方法，尤其對(duì)于探測(cè)目標(biāo)、雜波以及噪聲信號(hào)存在重疊區(qū)域的情況，本文能夠有效地對(duì)雜波和噪聲信號(hào)進(jìn)行分離，并能通過(guò)不同光譜成分的加權(quán)得到所需的目標(biāo)圖像。實(shí)驗(yàn)結(jié)果很好地證明了該方法的有效性。