S變換在地質雷達隧道襯砌檢測數據處理中的應用研究

朱自強, 王玉琴, 魯光銀

(中南大學 地球科學與信息物理學院，長沙 410083)

S變換在地質雷達隧道襯砌檢測數據處理中的應用研究

朱自強, 王玉琴, 魯光銀

(中南大學 地球科學與信息物理學院，長沙 410083)

通過將S變換應用于地質雷達探測隧道襯砌的數據處理中，解決了由于襯砌表層密集鋼筋網的影響導致雷達數據剖面受干擾嚴重、深部信號弱、病害異常難以識別的問題。進行了地質雷達二維有限差分數值模擬，研究了淺部密集鋼筋網對深部鋼拱架和脫空區(qū)造成的影響情況，并對原始剖面進行偏移成像壓制了繞射波和多次波，提高了鋼筋網下方弱反射信號的成像效果；最后通過對模擬數據和實際資料的S變換時頻分析，準確清晰地識別出了深部鋼筋和脫空區(qū)的特征。通過利用S變換時頻處理方法，為受干擾嚴重的隧道襯砌檢測數據的處理提供了一種方案。

S變換； 地質雷達； 偏移成像； 隧道襯砌檢測

0 引言

隧道襯砌檢測是隧道質量安全檢測最重要的環(huán)節(jié)，檢測結果影響著整個隧道的施工和運營。地質雷達技術作為一種無損方式的淺層地球物理勘探方法[1-2]，以其分辨率高、圖像直觀、工作效率高等技術優(yōu)勢，使地質雷達成為隧道襯砌檢測的常規(guī)儀器。

關于地質雷達的數值模擬研究已經比較成熟[3-4]，在襯砌數值模擬上，肖都[5]模擬了不同雷達天線探測襯砌的效果研究；周斌[6]對雷達襯砌檢測的數值模擬和工程實例進行了研究；舒志樂等[7]、黃玲等[8]研究了隧道襯砌病害的探地雷達三維正演模擬及工程驗證。但是在實際工程中，施工環(huán)境復雜，干擾嚴重，實際數據與數值模擬還是有一定的差別，所以要注重在復雜環(huán)境下的資料處理，數據處理的準確性直接決定了整個檢測結果。由于地質雷達受鋼筋的影響，整個雷達剖面中存在較多干擾信號，影響正常的數據解釋，因此如何在鋼筋影響下準確進行異常識別就成為了關鍵。張家松等[9]研究了希爾伯特變換的原理并將其應用在隧道含水裂隙識別處理中，但是對含水裂隙通道的具體尺寸大小難以斷定；朱自強等[10]研究了Curvelet變換處理隧道裂隙水GPR數據，但是實際資料處理比較困難，對于強干擾下目標的識別還是不夠清晰與準確直觀。

Stockwell等[11]提出了 S 變換時頻分析方法，該方法的時頻窗口具有可調性質，頻率的倒數決定了窗口的尺度大小，適合用于非平穩(wěn)信號的分析與處理。S 變換集中了短時傅里葉變換和小波變換的優(yōu)點[12]，其時窗寬度隨頻率呈反向變化，即在低頻段的時窗較寬, 從而獲得較高的頻率分辨率；而高頻段的時窗較窄，故可獲得很高的時間分辨[13]。S變換方法克服了其他濾波方法的缺陷，和小波變換和短時傅立葉變換相比，S變換具有較好的無損可逆性和時頻分辨率[14-15]，在地震和地質雷達上具有廣泛的研究價值，孫長蘋[16]將S 變換用于地震信號薄層預測處理；張先武等[17]研究了將廣義S變換應用在地質雷達的層位識別中；黃航如[18]對S變換在探地雷達數據處理中的應用進行了研究。

筆者提出將F-K偏移[19]與S變換結合運用于密集鋼筋干擾下的地質雷達襯砌檢測數據處理中，更加清晰地識別出了深部鋼筋位置和異常信息，準確地判斷出缺陷的范圍。

1 F-K偏移成像與S變換原理

1.1 F-K偏移成像原理

設雷達記錄剖面為e(x,z=0,t),其二維傅立葉變換為

E(kx,z=0,ω)

(1)

其中：

E(kx,z=0,ω)=?E(x,z=0,t)·e-i(kx+ωt)dxdt

在頻率-波數域作波場外推，得到深度Z處的波場為式(2)。

E(kx,z,ω)=E(kx,z=0,ω)e-ikzz

(2)

由于地質雷達的電磁波傳播是雙程走時，而波場偏移成像是基于單程傳播，所以電磁波的傳播速度為實際傳播速度的一半。通過對E(kx,z,ω)的傅立葉反變換得到式(3)。

(3)

(4)

1.2 S變換原理

S變換是類似于Gabor變換[20]、短時傅立葉變換、連續(xù)小波變換[21]的一種方法，它的推導方法有很多種，可以由短時傅立葉變換和小波變換退推導而來，Stockwell使用的是由小波變換推導而來，其過程如下：

信號h(t)的小波變換定義為式(5)。

(5)

將S變換定義為特殊基本小波與相位因子的乘積：

S(τ,f)=ei2πfτW(τ,d)

(6)

而基本小波定義為式(7)。

(7)

其中擴張系數d是頻率f的倒數。

基本小波式(7)不是嚴格意義上的小波，明確寫出S變換的公式為式(8)。

(8)

(9)

2 數值模擬及S變換處理

2.1 模型一

利用有限差分方法模擬了地質雷達探測隧道襯砌結構的剖面，模型為雙層鋼筋，淺部鋼筋密度大，模擬密集鋼筋網，深部鋼排列間距稍大。天線中心頻率為1 GHz，背景介質的介電常數為6，電導率為0.01 S/m。鋼筋的參數比較大，選擇適中的參數為介電常數為100，電導率為1 S/m。模型大小為長2 m,深1 m，在距模型深0.25 m處有9根間距為0.2 m、半徑為0.02 m的鋼筋；在深0.55 m處有5根間距0.4 m、半徑0.03 m的鋼筋。模型示意圖和雷達剖面圖見圖1和圖2。

圖1 模型一Fig.1 Model 1

圖2 模型一雷達原始剖面圖Fig.2 Radar original profile model 1

圖3 模型一偏移處理后的剖面圖Fig.3 Profile after migrationto model 1

由圖2可知，淺部密集鋼筋網的雙曲線反射強烈，并影響了深部鋼筋的反射，使得深部鋼筋的雙曲線反射弱，從原始圖中難以識別出鋼筋的數量與位置信息，需要進一步的處理。

將圖2中的雷達原始圖進行F-K偏移處理，得到如圖3所示的剖面圖，從圖3中可以看出，剖面圖中的雙曲線反射得到了很好地校正，消除了多次波和繞射波的干擾。淺部鋼筋的數量和位置都能清晰地辨別出，但是深部的鋼筋信息還是很難確定。

不同介質有不同的結構特征，內部反射波的高、低頻率特征明顯不同，這可以作為區(qū)分不同物質界面的依據，根據反射回波的頻譜信息進行識別目標，S變換能同時反映信號的頻率與時間信息，這是它作為信號處理技術的優(yōu)勢。我們提取剖面中水平位置0.2 m處的信號(圖4)，在旅行時5 ns和10 ns處具有強振幅且此兩處的頻率不同，根據這個特點，我們對此處信號進行S轉換，得到該處的時頻圖(圖5)。

圖4 圖3中水平位置0.2 m處的wiggle圖Fig.4 The wiggle diagram of Figure 3 in the horizontal position 0.2 m

圖5 圖3中水平位置0.2 m處信號的S變換時頻圖Fig.5 S transform time-frequency map of Figure 3 in the horizontal position 0.2 m

圖6 圖3中水平位置0.4 m處信號的S變換時頻圖Fig.6 S transform time-frequency map of Figure 3 in the horizontal position 0.4 m

根據圖5，在旅行時5 ns和10 ns處有頻率為1 GHz的異常信號，由此可分辨出兩個不同頻率的信號，分析是淺部鋼筋網和深部鋼筋的信號。對比模型，我們根據時頻圖分辨出了淺部密集鋼筋和深部鋼筋的位置信息，與實際位置有很好地吻合。對比在深部沒有鋼筋的位置，比如在水平位置0.4 m處的時頻圖(圖6)。由圖6得知，此處深部不存在鋼筋，所以在10 ns處并不見其他頻率的信號，符合實際。為了更進一部說明整個剖面的情況，我們對剖面進行S變換頻率切片處理(圖7)。從圖7中，我們能更加清晰地識別出剖面中的鋼筋信息。對這個模型的處理說明了，S變換能解決密集鋼筋網干擾下，深部鋼筋的識別問題。

圖7 圖3中S變換頻率切片圖Fig.7 S transform frequency slices map of Figure 3

圖8 模型二Fig.8 Model 2

2.2 模型二

本次模擬地質雷達對襯砌后的空洞缺陷的檢測，天線裝置和介質參數與模型一相同，只是在深0.70 m～-1.00 m范圍內設置了一個巖石層，且初支與巖石層的接觸處存在不平整面，在巖石層中設置了一個空洞，模型示意圖和剖面圖見圖7和圖8。

圖9 模型二雷達原始剖面圖Fig.9 Rdar original profile model 2

圖10 模型二偏移處理后的剖面圖Fig.10 Model 2 profile after migration

圖11 圖10中水平位置0.8 m處信號的S變換時頻圖Fig.11 S transform time-frequency map of Figure10 in the horizontal position 0.8 m

圖8中的原始雷達剖面圖受多次波和繞射波干擾嚴重，深部鋼筋和空洞的信號基本在干擾下難以識別，我們對其進行偏移歸位處理，得到偏移后的剖面圖(圖9)。圖9中剖面圖能得到淺部鋼筋的信息，但是仍然受強淺部信號和不平整接觸面的影響,導致深部的鋼筋和空洞信號很弱，對于深部的信號難以識別，所以需要進行處理，我們選擇對相應位置的信號進行S變換。水平位置0.8 m、0.9 m、1.0 m三處的S變換如圖10～圖12所示。

圖12 圖10中水平位置0.9 m處信號的S變換時頻圖Fig.12 S transform time-frequency map of Figure 10 in the horizontal position 0.9 m

圖13 圖10中水平位置1.0 m處信號的S變換時頻圖Fig.13 S transform time-frequency map of Figure 10 in the horizontal position 1.0 m

圖10顯示了在5 ns處的有異常區(qū)域，且頻率大小范圍為0.5 GHz～1.2 GHz，可知此處只是淺部鋼筋網的異常，結合圖8可見時頻圖劃定的鋼筋異常與實際十分符合；圖11顯示了在5 ns和18 ns兩處有異常區(qū)域，且頻率大小范圍為1 GHz左右，此處是淺部鋼筋網和深部空洞的異常，結合圖8可見時頻圖劃定的異常與實際符合；圖12顯示了在5 ns、13 ns、19 ns三處有異常區(qū)域，且頻率大小范圍為0.5 GHz～1.5 GHz，可知此處是淺部鋼筋網、深部鋼筋、深部空洞的異常，結合圖8可見時頻圖劃定的異常與實際吻合。為了更進一部說明整個剖面的情況，我們對剖面進行S變換頻率切片處理，見圖14。從圖14中，我們能更加清晰地識別出剖面中的鋼筋和空洞信息。

圖14 圖10中S變換頻率切片圖Fig.14 S transform frequency slices map of Figure10

圖16 實測數據雷達偏移處理后的剖面圖Fig.16 Actual radar profile after migration

2.3 實例

為了驗證S變換在實際資料處理中的效果，對實際含兩層鋼筋的混凝土進行探測，工作雷達天線頻率是1 GHz，圖15是探測到的雷達原始剖面圖，圖16是偏移處理后的剖面圖。

由圖16仍然很難識別出深部鋼筋的信息，于是我們將其進行S變換頻率切片處理，得到圖17的切片圖。從圖17中，我們能更加清晰地識別出淺部和深部的鋼筋信號。

圖17 實際數據S變換頻率切片圖Fig.17 S transform frequency slices map of actual data

在S變換圖中，通過對每個位置的數據道進行時頻分析，我們能清楚地得到該位置上的異常變化，通過時頻圖中的時間信息我們可以得到異常的深部位置及深部延伸，通過頻率信息，我們可以得到異常的特性。上述兩個模型及實例，說明了S變換在解決此類干擾下識別異常的可行性。

3 結論

筆者利用S變換對模擬隧道襯砌檢測數據和實例進行了處理，表明了S變換在異常識別上的優(yōu)點，能在表層密集鋼筋網的強干擾下識別出深部的異常信號。通過模型模擬與實例研究得出結論：利用不同異常體反射系數的頻譜特征，S 變換可以在得到反射回波的時間和頻率的信息，據此可以有效地提高信號的時頻分辨率，而且還可以根據實際信號具體情況，選擇合適的基本小波進行計算，靈活度高、實用性強。基于 S 變換的單頻曲線或單頻剖面，可以真實、細致地刻畫出深部異常的分布，使其具有很高的分辨率，為雷達資料高分辨處理研究工作，為下一步雷達資料處理和解釋工作提供借鑒。

[1] 李大心.探地雷達方法與應用[M].北京:地質出版社,1994.

LI D X. The methods and application of GPR[M].Beijing: Geological Publishing House,1994. (In Chinese)

[2] DAVID J DANIELS. Ground-penetrating radar. The Institution of Electrical Engineers[M]. London, United Kingdom,2004.

[3] 戴前偉,馮德山.時域有限差分法在地質雷達二維正演模擬中的應用[J]．地球物理學進展,2004,19(4):898-902.

DAI Q W, FENG D S. The apply of finite difference time domain method in ground penetrating radar two-dimension forward simulation[J]. Progress in geophysics, 2004,19(4):898-902. (In Chinese)

[4] 劉四新,曾昭發(fā).頻散介質中地質雷達波傳播的數值模擬[J].地球物理學報,2007,50(1):320-326.

LIU S X, ZENG Z F. Numerical simulation for ground penetrating radar wave propagation in the dispersive medium[J].Geophys,2007,50(1):320-326. (In Chinese)

[5] 肖都,李文杰,郭鵬.基于GPRMax的隧道襯砌檢測數值模擬及應用[J].物探與化探,2015,39(4):855-859.

XIAO D, LI W J, GUO P. The tunnel lining detection numerical simulation based on GPRMax and its application[J].Geophysical and Geochemical Exploration,2015,39(4):855-859. (In Chinese)

[6] 趙峰,周斌,武永勝.探地雷達在隧道襯砌空洞檢測的正演模擬應用研究[J].鐵道建筑,2012(8):99-103.

ZHAO F, ZHOU B, WU Y S.Study on application of forward modelling for inspecting void in tunnel lining by ground probing radar[j].Railway Engineering,2012(8):99-103. (In Chinese)

[7] 舒志樂,劉新榮.隧道襯砌病害探地雷達三維正演模擬及工程驗證[J].中國鐵道科學，2013,34(4):46-53.

SHU Z L, LIU S X.GPR three-dimensional forward modeling of defects in tunnel lining and engineering verification[J].China Railway Science,2013,34(4):46-53. (In Chinese)

[8] 黃玲,曾昭發(fā).鋼筋混凝土缺陷的探地雷達檢測模擬與成像效果[J].物探與化探，2007,31(2):181-185.

HUANG L, ZENG Z F .Modeling and imaging of GPR detection of reinforced concrete defects[J]. Geophysical and Geochemical Exploration,2007,31(2):181-185. (In Chinese)

[9] 張家松,朱自強,魯光銀,等.偏移成像與希爾伯特變換在地質雷達探測隧道裂隙流水通道中的應用[J].地球物理學進展,2012,27(6):2722-2728.

ZHANG J S, ZHU Z Q, LU G Y, et al. Application of migration and Hilbert transform on exploration of karst fractured ground water channel in the tunnel with GPR[J]. Progress in geophys,2012,27(6):2722-2728. (In Chinese)

[10] 朱自強,朱賀,魯光銀,等.基于Curvelet變換的隧道裂隙水GPR數據處理研究[J].物探化探計算技術,2014,36(5):571-576.

ZHU Z Q, ZHU H, LU G Y, et al. Processing of GPR data in tunnel fissure water based on Curvelet transform[J].Computing Techniques for Geophysical and Geochemical Exploration,2014,36(5):571-576. (In Chinese)

[11] STOCKWELL R G, MANSINHA L, LOWE R P.Localization of the Complex Spectrum: TheSTransform[J].IEEE Transactions on signal Processing,1996,44(4):998-1001.

[12] 趙淑紅,朱光明.S變換時頻濾波去噪方法[J].石油地球物理勘探,2007,42(4):402-406.

ZHAO S H, ZHU G M. Time-frequency filtering to denoise by S transform[J]. Oil Geophysical Prospecting,2007,42(4):402-406. (In Chinese)

[13] 操峰,陳淑珍,肖博勛.一種改進的STFT在PGPR回波信號處理中的應用[J].計算機應用研究,2005,22(6):127-129.

CAO F, CHEN S F, XIAO B X. Analysis of echo signal by an improved STFT algorithm for phased-array ground penetrating radar[J].Application Research of Computers, 2005,22(6):127-129. (In Chinese)

[14] 王曉凱,高靜懷,陳文超,等.基于高維連續(xù)小波變換的地震資料不連續(xù)性檢測方法研究[J].地球物理學報,2016,59(9):3394-3407.

WANG X K, GAO J H, CHEN W C, et al.Detecting method of seismic discontinuities based on high dimensional continuous wavelet transform[J].Geophys,2016,59(9):3394-3407. (In Chinese)

[15] 孫月.基于廣義S 變換和閾值函數的地震信號去噪研究[D].吉林：吉林大學,2012:13-26.

SUN Y. Seismic signal denoising research based on generalized S transformation and threshold function[D].Jilin:Jilin University, 2012:12-26. (In Chinese)

[16] 孫長蘋.S變換地震信號薄層預測處理技術[D].大慶：東北石油大學,2011:36-52.

SUN C P. Thin bed prediction processing technique of seisemic signal based on S transform[D].Daqing: Northeast Petroleum University,2011:36-52. (In Chinese)

[17] 張先武,高云澤,方廣有.帶有低通濾波的廣義S變換在探地雷達層位識別中的應用[J].地球物理學報,2013,56(1):309-316.

ZHANG X W, GAO Y Z, FANG G Y. Application of generalized S transform with lowpass filtering to layer recognition of ground penetrating radar[J]. Geophys,2013,56(1):309-316. (In Chinese)

[18] 黃航如.S變換在探地雷達數據處理中的應用研究[D].廣西:桂林理工大學,2012:24-37.

HUANG H R. The study of the application of S transform on ground penetrating radar data progressing[D].Guangxi: Guilin University of Technology,2012:24-37. (In Chinese)

[19] 李學聰,劉伊克,常旭,等.基于F-K偏移和反偏移的地震道插值方法研究[J].地球物理學進展,2009,24(6):2020-2029.

LI X C, LIU Y K, CHANG X, et al. Interpolation for seismic trances based on F-K migration and demigration[J].Progress Geophys,2009,24(6):2020-2029. (In Chinese)

[20] GABOR D. Theory of communication[J].Journal ofthe Institution of Electrical Engineers,1946,93(3):429-457.

[21] GOUPILLAUD P,GROSSMANN A,MORLET J.Cycle-octave and related transforms in seismic analysis[J].Geoexploration,1984,23(1):85-102.

ApplicationofStransformindataprocessingofdetectionoftunnelliningbygroundpenetratingradar

ZHU Ziqiang, WANG Yuqin, LU Guangyin

(School of Geosciences and Info-Physics,Central South University,Changsha 410083,China)

In this paper, the S transform is applied to the data processing of ground penetrating radar to detect the tunnel lining. To solve the problem that the radar data profile is seriously disturbed and the deep signal is weak and the disease abnormity is difficult to be identified due to the influence of the dense steel mesh. We carried out two-dimensional finite-difference numerical simulation of ground penetrating radar, studied the influence of shallow dense steel mesh on the deep steel arch and the void space, used the migration imaging method to deal with the original profile, repression of the diffraction wave and multiple waves, improve the deep weak reflection signal imaging results; Finally, the characteristics of the deep rebar and the void are identified clearly and accurately by the S transform time-frequency analysis. By using the S transform time-frequency processing of the simulation data and the actual data. This paper provides a effective solution to deal with the data of tunnel lining which is seriously disturbed.

S transform; ground penetrating radar; migration imaging; detection of tunnel lining

2016-11-03 改回日期： 2017-01-04

朱自強(1964-)，男，教授，博士生導師，從事地質災害探測與監(jiān)測工作，E-mail：13507319431@139.com。

王玉琴(1991-)，男，碩士，主要從事地質雷達信號處理研究，E-mail：940603995@qq.com。

1001-1749(2017)06-0755-07

P 631.4

A

10.3969/j.issn.1001-1749.2017.06.07