吴荣新,张斯薇,吴海波
(安徽理工大学地球与环境学院,安徽 淮南 232001)
探地雷达是利用高频脉冲电磁波在地下介质中的传播特性来探测地下目标体分布形态与特征的一种地球物理方法。探地雷达具有快速、无损、高效的特点,现广泛应用于路面、桥梁和隧道等工程结构体检测以及地下管道、溶洞等异常体探测,可为工程评价及地质解释提供可靠的资料。但实测雷达数据不仅有反射和绕射等有效波,还混有干扰噪声,给目标体的识别造成困难。
电磁波与地震波具有一定的相似性,所以雷达数据处理常借鉴地震数据处理方法,如静校正、反褶积、偏移和滤波等。滤波方法多样,如自适应中值滤波在压制混叠噪声的同时可保持有效信号;二阶拉普拉斯边缘检测算子能突出探地雷达有效异常并提高信号识别能力。多尺度-多方向的邻域滤波器能有效去除数据中的杂波;奇异值分解可消除地面穿透雷达信号的随机噪声和直达波;而结合奇异值分解的Shearlet变换能保留相对完整、明显的同相轴;地震数据处理的f-k滤波结合扇形滤波器使用可滤除探地雷达噪声;雷达信号包含人为设施引起的杂乱反射波及尖峰或突变状非平稳成分。雷达信号是多散射中心、具有局部细微特征的合成信号,小波去噪技术能够聚焦到信号细节对雷达信号处理有效。
为此,本文针对探地雷达数据中的椒盐噪声,尝试采用小波变换的方法进行滤波处理,对探地雷达记录进行模拟和加噪处理,并使用基于不同小波基的离散小波变换对雷达信号进行滤波去噪。通过计算滤波后峰值信噪比,并对比分析其波形特征,选出适合于探地雷达数据椒盐噪声去噪的小波基,并应用于实测数据的滤波处理。
a
和平移参数b
的离散小波变换,表示为(1)
其中小波函数ψ
(t
)定义为平方可积函数,需满足条件(2)
通常情况下,经小波变换分解后的数据信号由两部分组成,一部分是有效信号对应的低频部分,也称为近似部分,另一部分为噪声对应的高频部分,也称细节部分。因此,运用小波变换对数据进行去噪时需选取合适的阈值λ
,用来控制去噪精度,常用的阈值有硬阈值与软阈值函数。本文采用硬阈值小波去噪,滤波处理的思路如图1所示。图1 波变换去噪流程
首先输入需处理的数据信号,再选择合适的小波基对其进行离散小波变换,并保留由此得到的近似部分。然后对细节部分做判断处理,若细节部分的系数小于λ
,判断结果为“是”,则直接对其归零;若判断结果为“否”,即系数大于等于λ
,则重复上一步。再次对其进行小波变换,保留近似部分,判断细节部分的系数,直到细节部分的系数小于λ
。最后对所有保留下来的近似部分进行重构,即得到去噪后的信号。在本文的研究中,主要分析小波基的选取对探地雷达椒盐噪声的滤除效果,通过对比分析不同小波基的小波变换对单道波形的重构效果,并结合其峰值信噪比的提升程度,选择合适的小波基应用于正演记录,进而推广到探地雷达数据椒盐噪声的去噪处理。
采用不同信噪比的正演模拟记录来测试不同小波基的滤波效果。设置尺寸为2.5m×1m的矩形模型,模型中含两个水平位置相同但深度不同的圆形空气管。
模型上层为湿砂岩, 其介电常数为6F/m, 电导率为0.005S/m, 相对磁导率为1, 其余物理参数均设置为0, 厚度0.5m。下层为混凝土,介电常数为20F/m,电导率为0.1S/m,相对磁导率为1,其余均设置为0,厚度0.5m。空气管水平位置均为1.15m, 埋深分别为0.25m与0.75m, 半径均为0.1m。
正演模拟参数设置如下:发射天线的起始位置坐标为(0.087 5m,0.152 5m),沿水平方向的移动步长为0.02m,频率为900MHz。发射信号为雷克子波,延迟时间为0。网格步长dx
、dy
均为0.002 5m,时窗为24ns,收发距为0.025m。得到的探地雷达响应记录如图2所示。如图2所示,时间深度1ns处为直达波。埋深0.25m处的空气管上下管壁的反射波振幅强;埋深
图2 正演记录
0.75m的空气管反射弧不明显,上下管壁反射波难以区分。7ns附近的水平向反射同相轴为湿砂岩与混凝土层的分界面。
提取正演记录中的第60道信号,增益处理后其波形如图3(a)所示。对增益后的记录添加高斯白噪,设置其信噪比分别为1、2、5、10,得到不同信噪比的探地雷达记录如图3(b)~图3(e)所示。
图3单道记录图
选用haar小波、紧支撑正交小波、双正交样条小波、Coiflet小波和近似对称的紧支集正交小波对不同信噪比的单道记录波形进行去噪处理,均分解重构5次,并计算去噪后的峰值信噪比,不同小波基处理后的信噪比值如表1所示。表1中可看出db1与haar小波去噪后的峰值信噪比值一致,因为当小波基阶数为1时,紧支撑正交小波相当于haar小波,即在实际应用中也可用db1小波代替harr小波。信噪比为1、2、5、10的记录波形在滤波后信噪比值提升最大的小波基分别为sym8、db10、db10、bior3.7(见图4)。
表1 不同小波基去噪后的峰值信噪比值
图4 不同小波基去噪结果对比
图4为信噪比为2时,正演记录采用不同小波基的去噪结果。黑色曲线为去噪后的信号,红色曲线为去噪信号与原始记录的残差。采用三种小波基去噪后,单道记录包含的椒盐噪声均被滤除,且整体波形相似,但局部存在差异。采用sym8小波基去噪后,浅层空气管反射波前的波形极值与增益信号存在时差。采用bior3.7小波基的去噪结果中2~3ns处波形过于圆滑,细节丢失。采用db10与bior3.7小波基去噪后,7.5~9ns处波形有抖动现象。结合去噪后峰值信噪比值的变化,综合认为db10小波的去噪效果较好。为验证效果,用db10小波处理不同信噪比记录,结果如图5所示。
图5 不同信噪比记录的去噪结果
图5中锯齿状的高斯白噪被滤除,信噪比为1的去噪图0~1ns波形与尾部均有抖动。其余信噪比的去噪图0~2ns波形较平缓,无突出的波峰波谷,与信号一致。4.5~6ns的波形平缓无极值,而增益信号有5处极值,这部分波形的重构效果较差。7ns附近与增益图一致,均有极小值,深层空气管的反射重构效果很好。
不含噪声的二维正演记录如图2所示,对二维正演记录添加高斯白噪声后得到信噪比为0.1的记录如图6所示。由于记录的信噪比值较低,有效的反射信号被大量椒盐噪声掩盖,反射弧和界面的反射同相轴难以准确识别,尤其是图6中红色箭头所示位置。
对含噪声记录使用db10小波基进行滤波处理,去噪结果与残差分别如图7(a)和图7(b)所示。去噪后记录的椒盐噪声明显减少,不同深度的空气管反射弧成像清晰,易于识别。不同深度的界面反射同相轴清晰可识别。图7(a)中箭头所示的位置,反射信号信噪比明显增强,细节得以体现。
图6 含噪声正演记录
(a)去噪结果 (b)残差图7 含噪声正演记录的小波变换去噪效果
实测记录为埋深接近、水平排列的三根管线探测剖面。对实测记录进行静校正和能量增益后,其记录如图8所示。图8中实测剖面100ns后广泛分布 “椒盐式”噪声,部分有效波被掩盖,造成管线位置和深度的识别困难。
选用db10小波基对实测记录进行去噪分析,分解重构均为2次,小波变换后的去噪结果与残差如图9所示。对比分析去噪及残差图看出原始记录100ns后“椒盐式”分布的噪声经小波变换后基本被滤除,探地雷达记录剖面有效反射波清晰并凸显。水平位置为10m、50m、 80m处的三处管线反射弧清晰,部分界面的反射波同相轴清晰可见,管线的反射弧和界面的反射同相轴均易于识别。
图8 实测探地雷达记录
(a)去噪后 (b)残差图9 实测记录去噪效果
(1)为分析不同小波基在一维探地雷达正演记录中对椒盐噪声的去噪效果,通过对比去噪信号与原始波形的差异,并计算去噪前后的峰值信噪比,认为db10紧支撑正交小波的去噪效果最好,滤波前后的峰值信噪比由2提高到17.366 3。信噪比为0.1的加噪正演记录经小波变换去噪后,能去除记录中含有的大量高斯白噪,突出管线的反射弧和界面的反射同相轴。
(2)利用db10小波基进行实测记录的滤波处理,发现小波变换能有效去除记录中“椒盐式”分布的噪声,突出管道的反射弧和界面的反射同相轴。
因此,小波变换能用于探地雷达椒盐噪声的滤除,但重点在于选取合适的小波基,才能保证去噪效果,服务于后续的处理以及解释。