APP下载

连续小波变换在探地雷达信号分析中的应用研究

2012-01-11朱德兵郭政学王树敏

物探化探计算技术 2012年5期
关键词:探地小波剖面

黄 敏,朱德兵,郭政学,王树敏

(1.河南省煤田地质局 物探测量队,河南 郑州 450009;2.中南大学 地球科学与信息物理学院,湖南 长沙 410083)

连续小波变换在探地雷达信号分析中的应用研究

黄 敏1,朱德兵2,郭政学1,王树敏1

(1.河南省煤田地质局 物探测量队,河南 郑州 450009;2.中南大学 地球科学与信息物理学院,湖南 长沙 410083)

广泛应用于工程勘探的探地雷达,由于其发射的高频电磁波在介质中迅速衰减、受环境噪声大等因素影响,探测剖面难以对相对较深的异常体有较好的响应。复信号分析提取三瞬参数可多角度分析信号,并能突出弱反射信息,但易受噪声干扰,且常规的滤波方法不能有效地去除这一干扰。将具有时~频双重局限性的小波变换引入到复信号分析中,解决了复信号分析易受噪声干扰的问题,恢复了其对异常提取的能力。将此算法用于实际雷达探测数据处理中,取得了良好的应用效果。

连续小波变换;复信号分析;探地雷达;异常提取

0 前言

探地雷达具有轻便灵活,适应性强,工作效率高,无损、分辨率高等诸多优点,在浅层地球物理勘探中有着广泛的应用[1]。探地雷达数据剖面直观,在信噪比较高、探测效果较好的情况下,能直接根据原始剖面进行解释、解译。但在大多数情况下,环境噪声大,浅层干扰严重,地质结构复杂或者所需探测深度较大,剖面往往不能直观解释。并且由于雷达发射的高频电磁波在介质中传播时衰减迅速,造成利用所接收的弱反射电磁波信号解释困难。所以怎样从复杂的原始数据剖面中提取弱反射信号,进行图像解译和异常分析对雷达资料的解释尤为重要。许多学者进行过相关的研究,柳刚等人[2~5]研究过小波变换在信号去噪和弱信号提取上的应用,但所用的都是离散小波变换,分解和重构较为粗糙。张英德等人[6~8]采用 Hilbert-Huang变换对信号进行分析和处理,但其存在瞬时频率波动,严重时还会出现附加的虚假成份,对信号的最高分析频率只能达到采样频率的1/4以下,信号的分析受到了限制[9]。

1 小波复信号分析弱信号提取的理论基础

1.1 相位、频率分析法的优势

球面电磁波在均匀介质中传播可用式(1)表示[10]。

其中 E0为天线的初始辐射电场强度;φ0为初始相位;波数k=α+jβ;α为相位常数;β为衰减系数。

此时电磁场的相位表示为式(2)。

φ1=ω0t+φ0+αr (2)

当地下存在反射系数为R的反射体时,其反射的电磁波表示为式(3)。

这时反射电磁波的相位为式(4)。

对比式(2)和式(4)可以得出这样的结论:

(1)当存在异常体反射时,其电磁波的相位会引入由异常体产生的相移arg[Rδ(t-t1)]。

(2)角频率为相位的导数,同样由于异常体的引入,会产生arg′[Rδ(t-t1)]/2π的频移。

由于信号的相位、频率均与振幅无关,因此对信号的相位分析、频率分析可以同等大小来反映异常体反射信号,能更好地从弱反射信号中提取出相应的地质信息,这是振幅分析所不能比拟的。

1.2 探地雷达数据的复信号分析

由于相位、频率分析法在弱信号提取上具有优势,提取信号中的相位、频率信息能更好地对雷达记录做出正确的解译。复信号分析能分离信号中的振幅、相位信息。

探地雷达系统中的反射信号可由式(5)表示[11]:

a(t)为信号的包络,其频带宽度Δω远小于ω0,一般我们认为雷达反射信号f(t)为窄带信号。f(t)的希尔伯特变换)如式(6)表达:

其中 φ(t)=ω0t+θ(t)。

从式(7)中可以得到分析信号的三瞬参数:瞬时振幅a(t),瞬时相位φ(t)及其导数瞬时角频率ω(t),分别如式(8)、式(9)、式(10)所示。

在许多文献[6、12、13]中,都将瞬时相位表示为:,其主值范围为(-π/2,π/2)。实际情况下f(t)取值可以为零,但在计算和成图时,要求f(t)≠0,这就造成相位-π/2和π/2的缺失;在过f(t)取值为零的点时,出现相位的跳变,并在求导计算瞬时频率时,出现尖峰脉冲。因此用该式表示相位所计算出的瞬时相位剖面和雷达记录剖面,不能保持良好的对应关系,且使瞬时频率剖面分辨率极低。由于反正弦和反余弦函数不存在这样的问题,用其表达瞬时相位优于反正切。用反余弦函数表示相位时,为使其取值和成图时与原始雷达记录剖面保持一致,将其相移π/2,使其取值区间为 [-π/2,π/2],其表达式为:

图1(a)和(b)分别展示的是200MHz天线正对铝板,前后移动范围为1m~23m时,对采集的数据采用反正切和反正弦表达式所计算出的瞬时相位剖面图。可以看出,采用反正切计算出来的瞬时相位剖面,不如反正弦计算出来的对弱反射信息的反映效果好。

下页图2中的(a)和(b)分别展示了对用反正切和反正弦表达式求取瞬时相位求导而计算出的单道瞬时频率的对比图。同样可看出,由反正切计算的瞬时相位求导得到的瞬时频率,由于相位的突变而出现大量的尖峰脉冲,这使得利用瞬时频率剖面对探测数据进行解释变得极为困难。

图1 计算瞬时相位剖面Fig.1 Arctangent(a)and arcsine(b)calculate the instantaneous phase profile

在对信号进行复信号分析过程中,信号中的噪声也被放大,使得多参数信息受噪声干扰严重,这给信号多参数分析带来了困难,甚至对数据的解释和判读还不如原始信号。复信号分析方法的这种缺陷限制了其应用范围。图2 计算单道瞬时频率Fig.2 Arctangent(a)and arcsine(b)calculate the instantaneous frequency of single-channel

1.3 探地雷达数据的小波复信号分析

1992年 Mallat提出李氏(Lipschitz)常数,可以用来衡量信号的奇异性[15],并证明通过某点小波变换的极大模可以求取此点的李氏常数。若某点x0的李氏常数为a,当|x-x0|<δ时,δ为任意的正实数,函数f(x)的小波变换满足|Wsf(x)|≤ksa,其中k为正实数,s为小波变换的尺度。信号的李氏常数取值区间为[0,1],因而信号的小波变换的极大模将随尺度s的增加而增大;相反的,噪声的李氏常数小于零,噪声的小波变换的极大模将随尺度s的增加而减小。现证明随机噪声的李氏常数a为负数。

设n1(x)为随机噪声,其方差为σ2,则:

其中 ψs=ψs(x-v)ψs(x-u)。

那么 E(|Wsn1(x)|2)= ‖ψ‖2σ2/s,其中E(·)表示求其期望,‖·‖表示求其L2(R)的范数。可看出噪声在小波算子作用下,其模平方的平均值随尺度s的增加而以的速度下降。因此,噪声将随着小波变换的尺度的增大而得到很好的压制。

利用有用信号和噪声小波变换的极大模随小波变换尺度的变化的规律,不仅可以区别有用信号和噪声,还可以有效地压制噪声。小波复信号分析就是将小波分析引入到复信号分析中,能更好地对信号进行有效的多参量分析,从多个角度分析雷达剖面使得解译更准确。

图3展示的是对一探测数据进行瞬时频率分析,不但受到了高频的噪声干扰,还存在一定的尖峰脉冲。对其进行连续小波正变换,其时间~尺度域幅值图如图4所示,重构尺度11~32,得到处理后的瞬时频率曲线如图5所示(见下页),可看出噪声和尖峰得到很好的压制。因此对复信号分析的剖面做小波处理能大大增强其可视效果。

1.4 连续小波变换及反变换的程序实现

式(13)给出了连续小波变换(CWT)定义:

连续小波变换重构公式为:

a≠0,a、b∈R。由函数ψ(t)经过尺度伸缩和平移得到。

待处理的信号进入计算机前须离散化,连续小波变换在计算机上程序实现过程为:

(1)首先对离散信号xn做离散傅里叶变换,如式(15)所示。

(2)再对小波基进行傅氏变换后就得到信号xn连续小波变换,如式(16)所示。

其中 s为尺度;dt为采样间隔,ωk如式(17)。

(3)最后进行连续小波变换重构,其重构公式为:

图5 单道信号小波瞬时频率曲线Fig.5 Single-channel signal wavelet instantaneous frequency curve

2 实测数据的小波复信号处理

下面是一段铁路实测数据,由于噪声的干扰且深部反射信号微弱,深部异常体难以识别。

其原始数据剖面和去除水平干扰后的剖面分别如图6、图7所示。

对去水平干扰后的原始数据进行复信号分析得到多参数剖面,如图8所示(见下页)。可以明显看出,虽然复信号分析剖面中瞬时相位剖面比原始数据剖面对深部弱信号具有更好的反映,但复信号分析受到严重的噪声干扰。

经小波复信号分析,得到小波复信号分析剖面,如图9所示(见下页)。可以看出,复信号分析中的噪声得到了很好的压制,剖面更加清晰,尤其是小波瞬时相位、小波瞬时频率能更好地提取深部介质弱反射信号。

在原始数据中,很难辨认绕射弧,在小波瞬时相位、小波瞬时频率剖面图中得以清晰的展现。小波瞬时振幅剖面也能较好地展示道床的起伏形态。因此采用小波复分析提取深部的弱信号的方法切实可行。

3 结论

在当前技术条件下,地质雷达的发射功率有限,因而接收的深部反射信号甚是微弱,其应用的勘探深度受到限制,弱信号的提取问题一直备受关注。小波复信号分析综合了复信号分析对弱信号的提取能力和连续小波变换对噪声的压制能力,能更好地对弱信号进行识别,提高地质雷达的勘探能力和勘探效果。

[1] 李大心.探地雷达方法与应用[M].北京:地质出版社,1994.

[2] 柳刚,李术才,薛翊国,等.基于小波变换的雷达低信噪比信号处理技术及应用研究[J].工程勘察,2009(9):85.

[3] 詹毅,梁昌洪,方广有.探地雷达回波信号处理中小波变换域滤波方法的研究[J].西安电子科技大学学报,1999,26(3):305.

[4] 师学明,张剑,刘梦花,等.小波变换在地质雷达干扰波压制中的应用[J].工程地球物理学报,2008,5(3):279.

[5] 许建文,汪坚,吴江宁,等.小波分析在地质雷达数据去噪中的运用[J].重庆交通学院学报,2007,26(1):108.

[6] 张英德,刘江平,刘良琼.Hilbert变换在地质雷达数据处理中的应用[J].工程地球物理学报,2004,1(4):349.

[7] 汤井田,化希瑞,曹哲民,等.Hilbert-Huang变换与大地电磁噪声压制[J].地球物理学报,2008,51(2):603.

[8] 皮红梅,刘财,王典.利用 Hilbert—Huang变换提取地震信号瞬时参数[J].石油地球物理勘探,2007,42(4):418.

[9] 黄海,陈祥献.Hilbert—Huang变换应用中的预处理方法研究[J].浙江大学学报,2007,41(3):431.

[10]韦宏鹄,杨顺安.探地雷达波的相位参数及其应用[J].地质科技情报,1999,18(1):102.

[11]肖兵,鲍光淑,赵秋梅,等.探地雷达复信号分析及改进[J].中南工业大学学报,1997,28(1):1.

[12]张恩厚,陈前新.瞬时振幅、瞬时相位、瞬时频率在地质雷达工程中的应用[J].安徽地质,2004,14(1):58.

[13]谢雄耀,万明浩.复信号分析技术在地质雷达信号处理中的应用[J].物探化探计算技术,2000,22(2):108.

[14]张志勇,胡敬浓,张丽娇,等.探地雷达复信号分析的几点讨论[J].物探化探计算技术,2006,28(2):146.

[15] MALLAT S,HWANG W L.Singularity detection and processing with wavelets[J].IEEE Trans Inform Theory,1992,38:617.

[16]杨新安.地质雷达检测铁路路基新技术[J].中国铁路,2004:41.

[17]DAVIS JL,ANNAN AP.Ground penetrating radar for high resolution mapping of soil and rock stratigraphy[J].Geophys Prospect,1989,37:531.

[18]MASER KR.Condition assessment of transportation infrastructure using ground penetrating radar[J].J Infrastruct Syst,1996(2):94.

[19]廖立坚,杨新安,杜攀峰.铁路路基雷达探测数据的处理[J].中国铁道科学,2008,29(3):18.

[20]杨峰,彭苏萍,张全升,等.地质雷达精细处理技术方法的应用研究[J].工程勘察,2007(4):70.

[21]黄吉林,黄忠来,张建中,等.探地雷达小目标信号时频特征分析[J].厦门大学学报,2011,50(1):17.

[22]JACK R,JACKSON P.Imaging attributes of railway track formation and ballast using ground probing radar[J].NDT&E International,32-1999:457.

[23]华剑飞.物探技术在路基病害检测中的应用[J].铁道标准设计,2002(7):56.

[24]DANIELS DJ.Surface Penetrating Radar electronics Communication Engineering Journal[J].IEEE,1996,8(4):121.

P 631.4

A

10.3969/j.issn.1001-1749.2012.05.17

1001—1749(2012)05—0593—06

2011-12-19 改回日期:2012-04-18

黄敏(1985-),男,回族,湖南常德人,硕士,主要从事电磁法勘探工作。

猜你喜欢

探地小波剖面
ATC系统处理FF-ICE四维剖面的分析
探地雷达法检测路面板脱空病害的研究
基于多小波变换和奇异值分解的声发射信号降噪方法
构造Daubechies小波的一些注记
基于超表面的探地雷达增强探测研究
全极化探地雷达系统
基于MATLAB的小波降噪研究
基于探地雷达法的地下管线探测频谱分析
基于改进的G-SVS LMS 与冗余提升小波的滚动轴承故障诊断
复杂多约束条件通航飞行垂直剖面规划方法