面向基礎設施的三維探地雷達屬性成像及信息可視化研究

朱家松, 雷占占, 羅享寰

(1.深圳大學土木與交通工程學院, 廣東 深圳 518060; 2. 深圳大學城市智慧交通與安全運維研究院,廣東 深圳 518060; 3. 廣東省城市空間信息工程重點實驗室, 廣東 深圳 518060; 4. 深圳市地鐵集團有限公司, 廣東 深圳 518026)

0 引言

探地雷達(GPR)是一種廣泛應用于城市基礎設施內部測量和成像的無損檢測方法,具有分辨率高、非破壞性和非接觸式連續測量等優點[1]。GPR的成像方式分為一維(A-scan)、二維(B-scan)和三維(C-scan)3種。利用B-scan對城市地下基礎設施測量時,需檢查一系列相鄰的剖面才可以確定地下目標的位置和大小。而C-scan可以從俯視圖的角度,以水平切片的方式來呈現測量結果,具有直觀且全覆蓋的優點,可表現目標物體的空間位置及幾何特性。因此,基于C-scan圖像的目標提取研究逐漸成為GPR領域的熱點[2-3]。

由于三維C-scan圖像的空間分辨率可以達到毫米級,這導致包含大粒徑顆粒的非均勻介質,如砂土中的孤石、混凝土大孔隙或地下植被根系等非目標反射物會在C-scan上表現為雪花狀噪聲,從而影響目標物體的準確提取。此外,目前GPR C-scan所展現的測量值僅有雷達波的反射能量,而雷達波包含的相位、色散、波速等信息都未得到充分利用,因此,解譯時可依據的信息非常有限。然而,城市基礎設施多是鋼筋混凝土結構,基礎設施檢測中不同構筑物往往會產生相似的反射強度,難以將其區分,進一步導致C-scan的多解性問題嚴重。目前,針對GPR C-scan圖像質量低、信息少的問題,國內外已有不少研究。例如: Luo等[4]量化了C-scan空間分辨率與對象參數之間的關系,建立了面向目標對象的標準化C-scan成像流程。近年來,為了從地球物理數據中提取更多的信息來解釋和探測目標,以提高數據解譯的準確度和效率,研究人員將探測信號的屬性提取與分析技術引用到GPR數據解釋過程中[5-11],其中基于小波變換和基于復信號技術的瞬時屬性提取方法在應用中取得了良好的效果。Moeller等[12]證明了利用頻率建立C-scan的可行性,然而對于多種屬性融合的C-scan成像效果尚未有量化的比較分析。

C-scan圖像質量通常由空間分辨率(像素尺寸)和顏色(像素值)2部分組成。目前,國內外已有許多學者關注C-scan的分辨率,即由步長、天線頻率和采樣密度等測量參數決定的測量信息密度分布,但是關于GPR C-scan中像素值(信息可視化映射)的研究卻鮮有報道。C-scan的解譯是以人類視覺感知為基準,并沒有關于映射的客觀標準規則,解譯結果常帶有人為主觀誤差,甚至引發漏判或誤判。一般來說,C-scan的像素值由測量信息值映射到色彩區間得到,因而測量信息值的范圍和信息映射方法是決定C-scan色彩的2個主要變量[13-17]。Lu等[18]指出,一個合理的色彩是符合人類視覺感知的。因此,優化信息映射的主要目的是提高對比度和消除噪聲的影響,突出目標物體的特征,從而提高解譯的準確性。

綜上,目前三維GPR C-scan圖像主要存在可直接讀取的信息維度太低以及信息可視化方式主觀2個主要問題。因此,本文的研究內容由2部分組成: 1)提出一種融合波屬性信息的GPR C-scan圖像成像方法,豐富C-scan圖像信息; 2)定量分析C-scan不同測量信息可視化的成像效果,建立適合城市典型應用場景的優化成像策略,最終提高C-scan的解讀性以及GPR測量的解譯精度。

1 屬性C-scan成像及色彩優化方法

1.1 屬性C-scan成像

為實現GPR屬性提取與屬性C-scan成像,本研究基于python開發了GPR成像程序,主要功能包括基于B-scan的屬性提取與分析,提取屬性為瞬時振幅與瞬時頻率、GPR測量值的空間編碼以及屬性C-scan圖像的生成。GPR成像及解譯程序的用戶界面如圖1所示。

(a) 信號坐標編碼

(b) 波屬性提取

1.1.1 基于復信號技術的屬性提取方法

GPR屬性提取過程是通過數學統計與計算處理GPR數據的過程。GPR的屬性存在許多不同的分類方法,可以從計算的角度把GPR屬性提取劃分為2種類別: 1)依據單道波進行計算的GPR屬性,比如頻率、振幅、相位屬性; 2)依據多道波進行計算的GPR屬性,比如相干性和相似性。在這些理論的基礎上,本文采用單道波的復信號分析方法,提取GPR回波在目標物體深度的三瞬屬性,即電磁波信號的瞬時相位、瞬時振幅和瞬時頻率。其具體方法為: 通過采集GPR的實信號進行希爾伯特變換,將各個點的實信號通過卷積運算得到其實信號的虛部,再將兩者構成復信號并以三角函數的形式表達,然后從已分解的三角函數形式中提取各組信號的瞬時反射強度、瞬時頻率和瞬時相位[19]。由于瞬時相位主要由應用在GPR二維剖面上的同相軸連續性進行度量,因此,本文采用瞬時頻率與瞬時振幅作為三維圖像的研究屬性。

1.1.2 空間插值成像

提取GPR的回波屬性后,將其正交均勻排布,然后數值標準化以消除物理量綱的影響,并通過反距離加權法進行空間插值,以填充測量值采樣點間的空隙,形成全覆蓋的、高分辨率的C-scan圖像。

其中,反距離加權法的假設條件是事物間距越接近,則事物彼此越相似。其假定每個測點都對預測位置有局部影響,且離預測位置越近的測點局部影響越大。在使用反距離加權法插值時,測點的分配權重函數隨距離的增大而減小,測點與預測位置間的距離越近,則其被分配到的權重越大[20]。

由于GPR在試驗時采集的數據類型為剖面B-scan,所以,需要進行水平方向與垂直方向的重采樣,以模擬現實中的地表三維GPR圖像。根據Luo等[4]的三維C-scan成像流程,以目標物體直徑為C-scan厚度,對屬性的剖面圖進行垂直方向的重采樣,每個劃分深度為1組平面C-scan數據集,得到空間分辨率一致且臨近像素間無沖突的各屬性C-scan圖像。

1.2 信息映射優化

在形成全覆蓋的C-scan后,需要對C-scan中的測量信息值進行顏色賦值,以具現化測量信息。因C-scan的數據由三維空間信息(坐標編碼)和測量值2部分組成,而測量值僅為單一波段的相對反射強度信息,單色的灰度圖更符合真實情況。把測量值投影到灰度空間有2個步驟,分別為有效測量值閾值和映射函數。

1.2.1 有效測量值閾值

由于海量非目標雜波會在C-scan圖像上表現為雪花狀噪聲,而細碎的噪聲往往集中在測量值分布的高低兩端,定義灰度圖像的上下邊界可有效排除非目標回波[21]。將目標物體和周圍材料產生的GPR信號的測量組成為分布直方圖來描述其振幅分布,將各介質的測量值歸一為0～100,如圖2所示。

本研究采用標準差量化描述有效測量值的區間。首先,從直方圖中確定平均值;然后,根據標準差尋找與平均值有一定距離的上限和下限。其中,定義期望值為μ=0、標準差σ=0.5條件下的正態分布,則當σ=0.5,1,1.5,2,2.5,3,3.5,4時,其灰度取值范圍為全測量值的68.27%～99.99%。上下限之間的測量值為有效測量值,而高于上限和低于下限的測量值則將被去除。

1.2.2 映射函數

確定灰度邊界后,將有效測量值映射成灰度值,如下式所示。

S=Tr。

式中:S為轉化后的像素值;T為映射函數;r為映射前的測量值。

在8 bit圖像中,灰度圖像的像素值為0～255,如何分配灰度是信息映射的關鍵。

目前,GPR成像領域廣泛應用線性變換[13],因為它涉及較少的人為操作且提供了最完整的反射強度值。然而,針對不同的測量環境以及特定的測量物體,簡單的線性轉換未必是最合適的方式。例如: 當成像目標是一個不均勻的界面時,使用線性變換會生成一個噪聲較多且很模糊的C-scan圖像,而模糊的C-scan圖像難以描繪測量區域的細節信息。特別是當相鄰區域具有相似的介電特性時,低對比度的顏色會導致圖像結構信息的損失。除了線性映射以外,廣泛使用的4種不同的灰度映射函數分別為對數函數、平方根函數、平方函數和指數函數。其中: 對數函數和平方根函數關注較暗區域,可將像素值過大的區域做適當灰度減弱;而平方函數和指數函數關注較亮區域,即將像素值過小的區域做適當灰度增大。本文采用定量分析的方法詳細評估每種函數在基礎設施成像中的適用度。

2 基于實驗室環境的成像試驗

2.1 試驗數據

以典型的GPR應用案例——混凝土墻體檢測為例,混凝土墻GPR測線設計圖如圖3所示。

圖3 混凝土墻GPR測線設計圖

在混凝土墻上的正交網格進行GPR檢測,并收集墻體的GPR回波數據。GPR測線由17條平行于X軸和16條平行于Y軸的測線組成,測線間隔為100 mm。為了更容易地識別管道的影響,在設計網格時將網格橫軸縱軸交點中心設置在L形PVC管道的滲流點處。數據采集使用天線中心頻率為2 GHz的GPR設備,設置采集樣本的時間窗為6 s,道間距為3 mm,每個信號采集512個樣本。

在進行三維C-scan成像前,利用Reflexw軟件對采集的原始GPR回波數據進行基本的信號處理,其步驟包括去直流漂移、增益調整、零時校正、帶通濾波以及滑動平均。

信號處理前后的GPR B-scan樣本如圖4所示。雖然墻體內鋼筋與PVC管道呈雙曲線狀反射,且清晰可見,但兩者的圖譜特征太過相似,即使兩者引發的反射信號相位相反,但仍難以區分。

(a) 信號處理前

(b) 信號處理后

2.2 試驗結果

2.2.1 屬性靈敏度分析

在利用波屬性成像前,需確認雷達波屬性是否對不同介質響應不一,且需明確振幅或頻率對介質的靈敏性,以區分屬性在融合成像上的權重。在試驗區中選擇20個點位,涵蓋混凝土、鋼筋以及PVC管道,并提取這些點位的GPR回波信號。對已經過信號處理的回波信號進行復信號分析,最終得到關于混凝土、鋼筋以及PVC管道的三瞬屬性。從中選取瞬時振幅與瞬時頻率,并提取相應的屬性數據,將所提取的瞬時屬性的變化趨勢繪制成如圖5所示的曲線。

(a) 瞬時振幅

(b) 瞬時頻率

在圖5(a)中,鋼筋與PVC管道的曲線起伏變化較大,混凝土的曲線相對平穩,而且鋼筋在每一探測位置的GPR回波峰值振幅均比PVC管道大,因此,可以利用該振幅區分鋼筋與PVC管道。在圖5(b)中,鋼筋曲線的中心頻率變化較大,而PVC管道與混凝土的頻率相對平穩,且每個點位鋼筋的中心頻率都比管道的頻率大。由頻率的靈敏度可見,鋼筋與PVC管道的響應頻率曲線存在明顯不同,可依據該屬性區別鋼筋與PVC管道。

2.2.2 屬性C-scan成像

采用全測量值線性變換的鋼筋混凝土墻體C-scan與理想C-scan對比,結果如圖6所示。

(a) 瞬時振幅

(b) 瞬時頻率

(c) 精調最優

首先,采用自建程序確定試驗區內鋼筋網和PVC管道的深度,以15～30 mm為1個切片深度,根據采樣頻率確定C-scan的合適深度為150～300 mm。然后,根據當前業界普遍采用的全測量值線性投影的方法形成屬性C-scan原始圖像,如圖6(a)所示。靈敏度分析表明,PVC管道與混凝土的介電常數相接近,其振幅反射灰度值也相似,因此,振幅灰度圖無法判斷PVC管道的位置,但頻率灰度圖可以看清PVC管道的部分形態。基于對試驗場景的了解,人為精調得到接近真實情況的理想C-scan(如圖6(c)所示),可清晰界定同深度的鋼筋及PVC管道。

2.2.3 成像優化結果評價與試驗結果討論

為量化評價不同的測量值可視化方式對于C-scan圖像的優化效果,本文選取峰值信噪比(peak signal-to-noise ratio, PSNR)[22]與結構相似性(structural similarity, SSIM)[23]作為成像優化指標。PSNR越大,表示圖像質量越高;而SSIM的取值范圍在0～1,當SSIM為1時,表示試驗C-scan與理想C-scan的結構完全一樣。因為SSIM的計算不受色彩飽和度的影響,僅依據圖像中物體幾何形狀的完整度,適合用于本研究的目標物體區分評價。

對于測量信息可視化的2個變量——測量值閾值和映射函數,本文采用控制變量法研究其合適取值。首先,映射函數固定為常用的線性映射不變,測量信息為振幅,改變測量值閾值σ,形成不同測量值閾值的振幅C-scan圖像,與理想C-scan對比計算SSIM(如圖7所示)。比較SSIM數值后得知,閾值σ=1.0時得到最高SSIM,為0.90。然后,固定映射函數為線性映射不變,測量信息為頻率,改變測量值閾值σ,同樣當閾值σ=1.0時得到最高SSIM,為0.85。

同理可得,本文分別以振幅和頻率為對象,研究映射函數的選取。固定閾值σ=1.0,改變不同映射函數,形成映射函數不同的C-scan圖像,與理想C-scan對比計算SSIM。映射函數為平方函數時,SSIM最高,即成像最接近理想C-scan。試驗結果表明,可視化設置為有效閾值σ=1.0和平方映射函數時,振幅和頻率圖像的PSNR和SSIM均達到最高值,分別為32.74/0.90和29.96/0.85。

為探索融合屬性是否可提高C-scan圖像的可解讀性,將振幅和屬性的圖像按一定權重進行疊加,并根據上述試驗結果,測試測量值閾值屬性和映射函數的成像效果。結果表明,當有效測量值閾值σ=1.0,且灰度映射函數為平方函數時,融合屬性的SSIM最高,為0.90,高于振幅和頻率的最大值,且變化趨勢與振幅和頻率相符。融合屬性C-scan的PSNR高于單屬性頻率圖像,與單屬性振幅圖像接近(見圖8)。

(a) 瞬時振幅優化

(b) 瞬時頻率優化

(c) 屬性融合優化

2.2.3.1 屬性區別信息映射

上述試驗結果表明在選定同樣的有效測量值閾值(相同的σ值)與灰度映射方式時,瞬時振幅屬性C-scan的SSIM與PSNR明顯高于瞬時頻率屬性C-scan,而且多屬性融合后的SSIM結果明顯高于單屬性圖像,因為同時具備振幅和頻率的圖像所能讀取的信息更加豐富,所以與真實值更加接近。

2.2.3.2 有效測量值閾值σ

在屬性與灰度映射方式相同時,C-scan的成像效果僅由有效測量值閾值決定。上述試驗結果表明,噪聲主要集中在分布直方圖的兩端,應適當截選。分析上述結果可知,當σ=1.0時,即包含75.72%原始數據范圍時,優化效果最好。此時由于除掉了大部分的非目標回波,使得目標對象更加突出,因此,PSNR與SSIM均達到了最大值。當有效測量值范圍逐步增大時(σ增大),優化效果逐漸減弱。當σ大于1.5時,試驗結果趨于穩定。

2.2.3.3 信息映射方式

在屬性與有效測量值閾值固定時,C-scan的成像效果僅由灰度映射函數決定。試驗結果表明,在該組數據中平方函數和線性函數明顯優于其他信息映射函數,主要是由于鋼筋和混凝土的介電常數差異大, GPR在界面的反射強度大。平方函數與指數函數能較好地提高部分低反射區域的強度,使得目標對象更加清晰。同時,由于介質較為均勻,去噪需求較低,而線性函數一般較符合現實情況,可保留較多的測量值。對數函數與平方根函數不適用于低反射強度區域,因此,無論是PSNR還是SSIM均低于常規圖像數值。

3 基于道路現場的成像試驗

3.1 試驗數據背景

道路GPR C-scan成像試驗以香港石門一段長為20 m的道路為測量對象。該道路模型分為2個部分,左半部分為標準的瀝青道路結構,右半部分為標準的鋼筋混凝土道路結構,其面層包含0.1 m×0.2 m網狀鋼筋,路基下為不均勻砂土。道路結構中埋設2種地下管線設施,分別為PVC管道與球墨鑄鐵管道。石門道路地下結構及GPR測線圖如圖9所示,正交網格為GPR測線路徑。GPR數據由IDS雙頻系統(200/600 MHz)按照網格測線采集,并經過必須的信號處理。

S為管道直徑; D為管道埋深。

3.2 屬性成像與優化結果

3.2.1 屬性融合成像結果

經相關預處理并通過上述屬性提取與三維成像流程后,完成該案例的屬性C-scan圖像。

由于該試驗場地預埋材料復雜,砂土、混凝土與PVC管道均具有接近的介電常數,且3種材質的反射振幅值接近,導致振幅灰度圖像較為模糊,無法有效識別對應位置的介質。經對比,實地試驗中有效測量值閾值為1.0,映射函數為平方函數時,C-scan成像的SSIM最高,即效果最好。其優化后的C-scan圖像如圖10所示,2根管線均清晰可見,其PSNR與SSIM均為最高,分別為30.78和0.87。

(a) 瞬時振幅

(b) 瞬時頻率

(c) 屬性融合

3.2.2 屬性成像及映射優化的參數評價

石門道路試驗場地變參數成像效果如圖11所示。

上述試驗結果表明,在選定同樣的有效測量值閾值(σ值)與灰度映射方式時,瞬時振幅C-scan的SSIM與PSNR明顯低于瞬時頻率C-scan,且多屬性融合后的圖像SSIM結果明顯高于單屬性圖像。這主要是因為同時具備振幅和頻率的圖像所能讀取的信息更加豐富,能表現出更多的測量區域介質特征。

在確定相同的屬性與灰度映射方式時,當測量值閾值σ=1.0,即包含75.72%原始數據范圍時,成像效果最好。試驗的地下環境較復雜,噪聲較多,需要更窄的有效測量值,此時由于除掉了大部分的背景噪聲,使得目標檢測物體更加突出,因此,PSNR與SSIM達到最大值。當所選取范圍逐步增大時,優化效果逐漸減弱,試驗結果變化趨于穩定。

在確定相同的屬性與有效測量值閾值時,實地試驗中平方函數和線性函數明顯優于其他灰度映射函數,主要是由于鋼筋混凝土、管道和砂土的GPR反射強度差異不顯著,而平方函數與指數函數能較好地提高部分低反射強度區域,使得優化效果更好。綜上所述,城市基礎設施的結構相對復雜,介質種類較多,但單一介質內空間變異性不大,平方函數的映射方式可強化局部特征。

4 結論與建議

本研究提出一種基于波屬性的GPR C-scan成像方法,通過組合振幅和頻率提高C-scan信息維度,研究標準化的C-scan像素值映射方法,建立適合城市典型應用場景的面向對象C-scan成像策略。

通過開發三維GPR C-scan成像及解譯程序,在實驗室和實地道路構建典型GPR應用場景,定量分析三維GPR測量值可視化的2個主要參數,即測量值閾值與灰度映射函數的取值對C-scan成像的影響,并以PSNR與SSIM作為圖像質量的評價指標。混凝土墻體和實地道路的試驗結果一致,對于鋼筋混凝土結構,適當去除高頻低頻噪聲,保留70%有效測量值,信息映射方式為平方函數時,組合屬性的三維GPR圖像成像效果最好,可有效突出目標物體。對于環境不太復雜的城市基礎設施,需要較窄的有效測量值和更強化的灰度映射方式。因此,面向不同的測量場景和需求,建議建立規范的測量值可視化流程,可在保證C-scan目標物體區分精度的同時,減少人為操作誤差,進一步提高GPR測量的精度和準確度。