遥感科学与应用技术

基于包络检波和STFT谱分析的探地雷达土壤分层信息识别

  • 李俐 1, 3 ,
  • 付雪 2 ,
  • 崔佳 2 ,
  • 张超 1, 3 ,
  • 朱德海 , 1, 3, * ,
  • 吴克宁 4
展开
  • 1. 中国农业大学土地科学与技术学院,北京 100083
  • 2. 中国农业大学信息与电气工程学院,北京 100083
  • 3. 农业农村部农业灾害遥感重点实验室,北京 100083
  • 4. 中国地质大学(北京)土地科学技术学院,北京 100083
朱德海(1962— ),山东单县人,教授,博士生导师,主要从事农业遥感与土地信息技术研究。 E-mail:

李 俐(1976— ),河南南阳人,副教授,主要从事微波农业应用研究。E-mail: lilixch@163.com

收稿日期: 2019-05-30

  要求修回日期: 2019-09-16

  网络出版日期: 2020-04-13

基金资助

中国农业大学基本业务费项目(2019TC117)

版权

版权所有,未经授权,不得转载、摘编本刊文章,不得使用本刊的版式设计。

Soil Layer Identification based on Envelope Detector and STFT Spectrum Analysis of Ground Penetrating Radar Signals

  • LI Li 1, 3 ,
  • FU Xue 2 ,
  • CUI Jia 2 ,
  • ZHANG Chao 1, 3 ,
  • ZHU Dehai , 1, 3, * ,
  • WU Kening 4
Expand
  • 1. College of Land Science and Technology, China Agricultural University, Beijing 100083, China
  • 2. College of Information and Electrical Engineering, China Agricultural University, Beijing 100083, China
  • 3. Key Laboratory of Remote Sensing for Agri-Hazards, Ministry of Agriculture and Rural Affairs, Beijing 100083, China
  • 4. School of Land Science and Technology, China University of Geosciences, Beijing 100083, China
ZHU Dehai, E-mail:

Received date: 2019-05-30

  Request revised date: 2019-09-16

  Online published: 2020-04-13

Supported by

The Fundamental Research Funds for the Central Universities(2019TC117)

Copyright

Copyright reserved © 2020

摘要

土壤分层信息,特别是表土层结构,对土地生产力具有重要影响,是评价土壤质量的一个重要指标。为了快速、准确地获取土壤分层信息,本文利用探地雷达对分层土壤进行了回波信号采集,并分别在时域和频域分析土壤层位置和层厚信息。首先在信号预处理的基础上,借助包络检波方法确定在土壤分层界面在时域上的位置;然后获取电磁波速度,得到土壤分层厚度。考虑到土壤介电常数与电磁波在土壤中传播速度的相关性,采用短时傅里叶变换方法(Short-time Fourier Transform,STFT)获取各土壤层时频域特征值,并利用回归分析建立特征值与介电常数之间的数学关系,实现对各土壤层的介电常数估算,从而计算出电磁波传播速度,进而确定土壤各层厚度。为验证算法的有效性,分别对理想模拟实验环境和农田环境进行了探地雷达实验,结果表明利用包络检波对探地雷达回波信息进行分析,土壤层检出率达到94.5%,借助STFT谱分析进行探地雷达回波速度估计,对于70 cm深度以上土层厚度计算误差大都保持在10%以下,但随着土壤深度的增加,误差变大。总体来说,本方法能有效识别浅层土壤的分层信息,可应用于实际生产中耕层厚度的估测。

本文引用格式

李俐 , 付雪 , 崔佳 , 张超 , 朱德海 , 吴克宁 . 基于包络检波和STFT谱分析的探地雷达土壤分层信息识别[J]. 地球信息科学学报, 2020 , 22(2) : 316 -327 . DOI: 10.12082/dqxxkx.2020.190265

Abstract

Soil stratification information, especially surface soil structure, has important impact on land productivity and is an important index for evaluating soil quality. The present study aimed to obtain the information of soil layers quickly and accurately, for which Ground Penetrating Radar (GPR) technology was used. The echo signals of GPR were processed in both the time and frequency domains. In the time domain, the envelope detection method was used to determine the transience of the echo signals and therefore to get the location of soil layers on the time axis. To get the soil layer location in spatial coordinates, the velocity of electromagnetic wave propagation in soil was needed. Considering the velocity of electromagnetic wave propagation in soil layers varying with the soil dielectric constants, the Short-Time Fourier Transform (STFT) method was applied to the echo signals for dielectric constant analysis in the frequency domain. Soil layers with different dielectric constant exhibited different characteristics in the STFT signals. After clustering analysis of the soil layers, the relationship between STFT characteristic value and dielectric constant in a certain layer was established based on regressive analysis. Then, the velocity of electromagnetic wave propagation in each soil layer was determined using the dielectric constants. After the electromagnetic wave velocity of Ground Penetrating Radar (GPR) was estimated, the location of layers' interface was further determined and then the thickness of each soil layer was computed. To valid the effect of the above-mentioned methods, the echo signals of soil, for both the ideal simulated experimental environment which has obvious layered interface and the farmland environment whose layers have changed naturally, were collected. The experimental results show that, with the envelope detection method, layers not deeper than 70cm in the ideal simulated experimental environment were 100% detected and for both the ideal simulated experimental environment and the farmland environment, the detection rate of ground penetrating radar echo information reached 94.5%. The estimation of ground penetrating radar echo velocity using STFT spectrum analysis shows that the calculation error of soil thickness above 70 cm depth was mostly below 10%. Our findings suggest that the proprosed methodology can effectively identify the stratification information of shallow soil and estimate the thickness of the soil. However, surface vegetation, film mulching, soil voids, soil salinity, moisture heterogeneity, gradual change of soil layers, and soil layer depth will all affect the accuracy of the detection. For example, with the increase of soil depth, the error becomes larger. So, if the data acquisition spot is selected rationally, the proposed methodology can be applied to plough layer thickness detection in practical fileds.

1 引言

土壤是陆地上万千生物赖以生存的基础,是农业发展和人类赖以生存的重要环境。土壤分层信息和土壤质地信息是影响土壤质量的重要指标,特别是以土层厚度为主要内容的土壤分层信息在一定程度上决定了土地适合耕种的程度[1]。因此快速、准确地获知土壤分层,特别是土层厚度信息,是土壤研究的一个重要问题。常规对土壤分层信息的探测主要采用挖剖面或使用麻花钻、手持锥形贯入器等工具来获取土壤特定采样点的分层信息[2]。这种方法直观明了,不需要复杂的数据解译,但普遍存在费时费力、结果依赖于实验人员的个人经验、采样点可能缺乏代表性等缺点。
地球物理学研究和物性探测技术的发展为土壤分层信息获取提供了快速无损探测的思路。研究者尝试了使用地震波[3]和土壤电阻率断层扫描[1]等技术进行土壤分层信息获取,雷达波技术[4,5]也开始应用于浅层的地表信息探测中。国内的研究还主要集中于探地雷达在土壤分层应用中的尝试。在国内,胡振琪等[5]将探地雷达应用于复垦土壤,通过实验检测定性地分析了探地雷达信号用于土壤分层信息获取的可行性;彭亮等[6]应用低频探地雷达大致刻画出了农田土壤剖面层位结构;于秀秀等[7]也通过实地验证对用探地雷达测量土层厚度变化进行了尝试。国外的研究则开始了对算法和模型的尝试。Kataev等[8]利用散射矩阵技术来分析探地雷达回波信号以期实现地下结构层的探测;Ardekani等[9]针对层状土壤结构借助有限差分时域方法建立了探地雷达的B-scan图像模拟模型;Léger等[10]评估了使用探地雷达(Ground Penetrating Radar,GPR)来估计美国阿肯色州巴罗附近冰楔为主的冻土区的解冻层厚度、积雪深度和冰楔特征,并将所得的土层结构信息与电阻率层析成像(Electrical Resistivity Tomography,ERT)得到的结果进行了对比分析。总体来说,探地雷达由于其无损、快速、探测目标范围广、分辨率高、操作方便等特点,在土壤信息识别中的应用研究已受到了研究者的重视[11,12]。然而,就目前国内外的研究情况来看,土壤分层信息的研究还处于初级阶段[13],大都只是验证了探地雷达用来进行土壤分层研究的可行性,特别是对于各层界面清晰、介质存在较大的典型差异时,土层分界面及层厚度较易辨识的研究较多,但对于各层土壤质地差异不明显的大田土层定量化的表达和高精度的结果还需进一步研究。
本文针对理想模拟实验环境和实际农田环境,提出基于包络检波和STFT谱分析的探地雷达土壤分层信息识别方法。首先借鉴调制解调的原理,利用包络检波的方法分析雷达回波信号幅度以及上下包络Hilbert瞬时相位谱的变化,找到土层分界面的时域位置。然后,为了定量确定土层厚度信息,利用短时傅里叶变换跟踪频谱随时间(深度)的变化,分析各层信号的频域特性,通过回归分析建立各层土壤的介电常数与各层信号的频域振幅峰值加权值之间的关系,估测不同层电磁波传播速度,进而计算各层厚度。本研究在时域分析的基础上对回波信号进行谱分析,以期为土壤分层信息提取提供自动化、高精度的方法。

2 实验设计与数据获取

2.1 模拟实验环境搭建及数据获取

为研究探地雷达对不同土层深度和不同土层界面的敏感性,在北京通州中国农业大学实验站搭建模拟实验环境。设计10 m×10 m的区域,等分成5 m×10 m的两块(左边以石子为障碍层,右边以瓷砖为障碍层),挖掘出10条长度为10 m、宽度为1 m、高度以10 cm为步长逐渐递增的阶梯状平面,分别以2.0系列石子和瓷砖为界面在长度10 m的阶梯上左右两边均匀铺设5 m×1 m(长×宽)的障碍层,最后覆土整平。覆土后最深处障碍层位于土壤下1 m位置,障碍层厚度约1 cm,其上下土壤为实验田常规砂质壤土,如图1所示。
图1 北京通州实验站模拟环境搭建及测线分布

Fig. 1 Simulation environment construction and line distribution of the experimental station in Tongzhou, Beijing

考虑到模拟环境的阶梯式土层厚度最深处为1 m,因此选用探测深度为0.8~1.2 m的CAS-S800(工作频率800 MHz)探地雷达在其上方以轮测法沿4条测线进行了数据采集,每条测线长度为10 m左右,采样时窗设置为40 ns,采样频率设置为 26 324 MHz。电磁波在介质中的平均速率根据第一障碍层位置估算。

2.2 实际农田环境研究区概况及数据获取

为进一步验证算法在实际农田中对土壤分层信息的测量准确性,本文选取了位于河北省曲周县的土地综合整治项目区作为实验采样区域。研究区位于114°50'22"E—115°13'27"E,36°35'43"N—36°57'00"N之间,面积约为667 km2,属于太行山山前平原南段、漳河冲积扇下缘,大地貌为山前平原,但县内小地貌繁多,在低平地为主的地形中,还有局部的二坡地、缓岗、决口冲积扇和河间洼地,南北方向大体呈带状分布,而从西到东,岗坡洼相间,呈波状起伏。地貌分布与土属变化关系密切,土壤类型以潮土为主,质地类型主要为轻壤土。土壤在垂直方向上砂壤黏交替沉积呈层状,土壤呈现出一定的分层特征。
在2016年4月开展了外业实验,根据项目区常见的4种土体结构类型(通体粘型、夹粘型、倒蒙金型、底漏型)分别选取了8个采样点进行了探地雷达数据采集、剖面数据采集和土壤样品采样。由于调研时间为2016年4月,项目区大部分农田生长着小麦、玉米等作物,有的样点数据不是很理想,因此剔除干扰、噪声较大的样点,最终选择了裸露地表的样点进行测算和验证,每种土体构型选择一个样本,共计4个采样点。
探地雷达天线中心频率越高,分辨率越高,但探测深度会下降[14];反之亦然,因此实际农田环境中选用探测深度更深、工作频率为400 MHz的GR-IV型便携式探地雷达进行数据采集。采样点位置如图2所示,其中图2(b)底图来自于Google Earth高清影像,每个采样点根据田块大小设置多条10~25 m长度的测线进行往返数据采集,采样时窗设置为25 ns,采样频率设置为100 000 MHz。
图2 位于河北省曲周县的研究区农田环境采样图

Fig. 2 Study area and the location of samples in Quzhou County, Hebei Province

为验证探地雷达结果的准确性,在探地雷达数据获取的同时,进行土壤剖面挖掘和分层土样获取工作。在外业根据土壤颜色、紧实度和触感黏度等特征判断分层质地并测量分层位置的基础上,进行内业实验,对各层土壤样品进行化验,获取土壤质地、土壤容重、含水量、土壤含盐量等指标参数。

3 土壤分层信息识别

3.1 土壤分层信息识别过程

为了快速、准确地获取土壤分层信息,需对探地雷达回波信号进行数据预处理,借助包络检波方法确定土壤时域分层信息,为获取电磁波在不同层的传播信息,采用短时傅里叶变换方法(STFT)计算基于时域分层的各土壤层的频域振幅峰值加权值,由频域振幅峰值加权值计算其各层土壤介电常数,并估计出回波信号的速度,进而确定土壤各层厚度。系统整体流程图如图3所示。
图3 基于包络检波和STFT谱分析的探地雷达土壤分层信息识别方法流程

Fig. 3 Flowchart of soil layer identification based on envelope detector and STFT spectrum analysis of ground penetrating radar signals

由于电磁波在土壤介质中传播将发生衰减,并且由于地面坡度、粗糙度等几何特性的不均一性,以及地表植被散射影响,接收机接收到的回波信号往往具有复杂的散射过程[15],包含各种各样的杂波干扰,使得土层识别十分困难,故首先对雷达回波信号进行增益处理来弥补信号在传播过程中的能量损失,并进行滤波与降噪处理以降低噪声,从而提高回波信号信噪比。本文所涉及的滤波处理过程主要包括:通过均值法进行去除背景的影响;通过带通滤波去除外界射频信号的干扰;通过中值法去除来自地面反射的直达波。然后,分别在时域和频域对回波信号进行处理以获取时域分层信息和估计各层电磁波传播速度:① 采用包络检波的方法进行土层分界面的时域确定,进而确定土层分界面对应传播时间信息;② 采用短时傅里叶变换对信号进行时频域分析,计算基于时域分层的各土壤层的频域振幅峰值加权值,利用建立的频域振幅峰值加权值与土壤介电常数的线性回归关系计算各层的土壤介电常数,从而确定电磁波在各土层的传播速度,最后结合土层界面时间信息确定各土层厚度信息。

3.2 基于包络检波的土壤层次划分

探地雷达发射的电磁波信号遇上障碍物或障碍层会产生回波,因此会改变回波信号的幅度、相位等信息。土壤物理及化学组成差异形成的不同土壤分层,因此也会引起探地雷达回波信号波形幅度和相位的变化。然而,由于干扰噪声的存在,回波信号本身幅度和相位的变化有时不能完全反映不同土层的情况,因此本文利用包络检波的方法从弱回波信号中检测出信号的变化并结合Hilbert变换进行上下包络瞬时相位分析[16,17],进而确定土壤分界面的时域位置。
包络检波是信号处理中一种常用的滤波检波方法,能从被低频信号调制的高频信号中解调出低频信号,对于信噪比较低的信号具有较强的识别能力。本文采用包络检波正是利用了其对信号信噪比要求较低的特点,从高频变化的信号中解调出包含土壤分层信息的低频回波信号。将回波信号表示为 r t ,搜索 r t 信号的所有极大值点和极小值点,然后进行三次样条拟合,获得 r t 的上包络信号 r u t 和下包络信号 r d t 。对上下包络线分别求希尔伯特(Hilbert)变换:
r u ˆ t = 1 π - r u τ t - τ
r d ˆ t = 1 π - r d τ t - τ
式中: τ 为上下包络希尔伯特变换中卷积运算引入的中间积分变量; r u ˆ r d ˆ 是上下包络希尔伯特变换结果。
进而求得上包络和下包络对应的解析信号 R u t R d t
R u t = r u t + i r u ˆ t
R d t = r d t + i r d ˆ t
则上包络信号和下包络信号的瞬时相位 Φ u t Φ d t 为:
Φ u t = arctan r ˆ u t r u t
Φ d t = arctan r ˆ d t r d t
通过识别上下包络信号Hilbert谱瞬时相位跳变,确定土壤分界面的时域位置。
为进一步降低噪声影响,在对土壤分层界面识别时,本文分别对一条测线的多个天线位置的回波信号进行识别,并采用大数判决的方法决定各地块土壤层数和各层对应的电磁波双程时延位置。

3.3 基于短时傅里叶变换的土壤介电常数与电磁波速度估算

土层由于其深度与介电常数不同会引起探地雷达回波信号频率信息的变化,因此分析其频谱信息是利用探地雷达超宽带频谱特性,获取土层介电常数信息的重要方式。短时傅里叶变换[18,19]可以通过获取信号在某一时间 t 的短时傅里叶频谱从而有效跟踪信号频谱特性对时间(深度位置)的变化。
回波信号 r t 的短时傅里叶变换是由以 t 为中心的窗函数对 r t 进行截断处理后再进行傅里叶变换得到的,可以表示为:
R STFT t , ω = - + r τ w t - τ e - jωτ
式中: τ 为短时傅里叶变换中卷积运算引入的中间积分变量; w t 为窗函数,实现对回波信号 r t 的截取,经过加窗截断处理可使非平稳信号在截取窗内近似可看成平稳信号。
短时傅里叶变换中一个重要问题是进行窗函数设计。本文选取汉明窗为窗函数。所得短时傅里叶变换信号可以理解为信号 r t 在分析时刻 t 的局域频谱。根据前面数据采集参数的设计,考虑到每道信号点数,本文窗函数宽度分别选取8/16/32/64进行了性能比较,最终选取了性能最佳的32点为窗宽进行后面的计算。
利用 R STFT t , ω 计算信号频域振幅峰值加权值:
A t = ω R STFT t , ω d ω
分析发现信号频域振幅峰值加权值 A t 与介电常数 ε 具有良好的线性相关性,利用回归分析方法建立 A t 与土壤介电常数 ε 间关系模型。
对于土壤介质来说,电磁波在土壤中传播的波速 v 可由介电常数 ε 确定[20,21]
v = c ε - 1 2
式中: c = 3 × 10 8 m s , 为电磁波在真空中的传播速度。

3.4 土层厚度信息获取

通过时域包络检波方法求得土层厚度对应的电磁波双程时延,通过频域STFT方法获取电磁波在各层土壤中传播速度,则各土层厚度可以表示为:
h = 1 2 v τ
式中: τ 为电磁波双程时延; v 为电磁波在各层土壤中传播速度; h 为各土层厚度。各层土壤厚度均可以根据探地雷达回波时延和由介电常数决定的电磁波传播速度计算获得。

4 结果及分析

4.1 基于包络检波的土壤分层检测结果分析

将基于包络检波的方法用于模拟实验环境中获取的雷达数据,得到如图4所示的土壤层识别结果。为更明确看出每一条回波数据上土壤层数据的表现特性,图中上半部分给出了基于包络检波及Hilbert谱瞬时相位跳变的三维识别结果(红黄色的条纹是根据识别的位置人为添加的障碍层识别结果示意),下半部分分别给出了天线位置在0.65、2.62和6.07 m处回波数据及其上下包络情况的二维曲线。
图4 基于包络检波的模拟实验环境探地雷达土壤分层信息识别结果

Fig. 4 Soil layer identification results based on ground penetrating radar by envelope detection in the simulated experimental environment

图4可以看出,三维土壤分层识别图中准确识别出了前6层台阶的位置,信号随着障碍层深度的增加会逐渐减弱,障碍层深度不大于65 cm时都能提取出障碍层位置,障碍层深度大于65 cm时回波信号极为微弱,因此为了图像显示简便,这里只保留了天线位置0~7 m的测量数据(对应障碍层位于0~70 cm)。从二维结果中可以看出,每个天线位置获取的数据中包含了当地对应台阶的位置,通过包络检波得到的上下包络线已经比回波曲线本身更能直观地显示出回波跳变的位置,但只通过看上下包络线距离的变化有时候还不够直观,特别是当障碍层深度增加,信号变化较小时。因此,本文进一步通过识别其Hilbert谱瞬时相位跳变顺利找到了造成跳变的障碍层位置,结果如图4上半部分三维结果所示。
为进一步验证本文算法在实际农田环境中的应用效果,将基于包络检波土壤分层信息识别技术用于了研究区的4个采样点。实验表明,本算法能在复杂的农田环境中检测出不同土壤层的存在, 图5以样点1和样点3为代表给出了分层界面三维显示。为以更直观的二维图形式展现回波数据及其上下包络线情况,图5中也给出了天线在随机选取的几个位置上样点1和样点3的二维数据及包络线。样点2和样点4的结果类似,限于篇幅,这里不再赘述。
图5 基于包络检波的实际农田环境探地雷达土壤分层信息识别结果

Fig. 5 Soil layer identification results based on ground penetrating radar by envelope detection in the actual farmland environment

与实地剖面结果比较可看出,土壤分层信息基本被检测出来,但由于实际农田环境中复杂的散射环境,在没有明显分层的位置出现了包络信号的突变,存在一定的虚警现象。样点1的探地雷达回波信号在天线位于某些位置时,如在双程时间10.8 ns处,检测出分层,然而实际剖面图中不存在此分层。出现这种虚警,可能是由于土壤中水分或盐分不均匀引起了土壤后向散射的变化。对于这种虚警现象,若其长度小于天线行走距离的1/3,则在后面确定层位置的运算中,都将其作为误差剔除。综合模拟实验环境和实际农田环境中所有测线数据的分层识别结果,实际土壤层检出率达到94.5%。

4.2 基于STFT的土壤介电常数确定

根据实际农田环境外业采集的土壤理化特性计算各层土壤的介电常数 ε ,获取与之对应的探地雷达数据并采用STFT方法计算其各层土壤的频域振幅峰值加权值A,典型样点时频域振幅峰值加权值与介电常数关系如图6所示。采用回归分析建立介电常数 ε 与信号频域振幅峰值加权值A之间的线性关系,可得:
ε = 8.5866 + 2.4373 × 10 - 6 A R 2 = 0.87
计算任意土壤层频域振幅峰值加权值A,并利用式(11)可以获取实际农田环境任意土壤层介电常数。
图6 频域幅度峰值加权值-土壤相对介电常数关系

Fig. 6 Relationship between peak weighted value of frequency amplitude and soil relative dielectric constant

4.3 实际农田环境土壤分层结果分析

对于实际农田环境获取的4个样点的探地雷达数据,分别利用基于包络检波的方法获取土壤分界面的时域位置、利用基于STFT的方法和式(9)-(11)获取电磁波在不同土壤层的传播速度及相应时间的传播深度,得到研究区典型样点土壤分层情况(图7)。图7(a)为4个样点探地雷达数据分层识别结果的三维图像在土壤深度天线位置二维坐标面的投影,图7(b)为对应样点实地剖面及土壤分层人工辨读结果。
图7 河北省曲周县研究区典型样点土壤分层情况

Fig. 7 Soil stratification of typical sample sites in study area in Quzhou County, Hebei Province

可以看出,基于包络检波和STFT谱分析的方法可以实现土壤分层层数和大概位置的准确识别。但由于土壤本身结构的复杂性和渐变性,存在有漏检现象。样点2在实地剖面中在地面下98 cm的位置有分层,但探地雷达数据检测中未检测出该界面,出现漏检。对照土壤取样和实地专家辨读结果发现,此样点土体构型为通体黏型,土壤物理化学特性变化较小,呈渐变性,土壤层次不够分明,并且由于黏土的穿透性较差,回波信号比较弱,因此在50 cm以下土层出现漏检。样点3的回波测定中第3层位置在53 cm处,而根据剖面实地专家辨读第3层位于68 cm左右,产生了较大误差。对比实地剖面照片发现,在50 cm到55 cm的位置土壤有空洞等存在,空洞回波信号强于渐变土层差异造成的回波变化,因此检测的土层位置上移,同时由此带来了速度估计误差及下面各层误差增大的现象。
为定量比较STFT算法获取的土壤介电常数和基于包络检波的土壤分层识别结果的精度,将探地雷达回波数据处理结果与实地野外测量和取样化验结果进行对比,如表1所示。表格中分层“0”对应空气和土地分界面,其双程走时前对应的回波信号在数据预处理时被截除,绝对误差为土壤层厚计算值与实测值之差,相对误差为绝对误差的绝对值与实测值比值的百分比。
表1 研究区各土体构型剖面点分层数据汇总

Tab. 1 Hierarchical data aggregation of soil configuration profile points in the study area

样点 分层 双程走时/ns 介电常数/(F/m) 土壤层厚/cm 层厚误差分析
实测值 计算值 实测值 计算值 绝对误差/cm 相对误差/%
样点1 0 4.65
1 8.32 9.23 9.49 18 17.87 -0.13 0.71
2 15.97 7.38 7.39 44 42.22 -1.78 4.05
3 18.43 10.19 9.96 10 11.68 1.68 16.80
4 20.60 10.36 9.86 10 10.34 0.34 3.45
5 23.55 9.79 9.83 26 14.15 -11.85 45.56
样点2 0 7.14
1 12.10 12.71 13.13 20 20.54 0.54 2.68
2 16.58 10.56 11.94 19 19.49 0.49 2.56
3 24.30 10.38 - 60 - - -
样点3 0 6.45
1 9.62 9.14 9.10 15 15.77 0.77 5.14
2 12.67 7.25 7.11 17 17.16 0.16 0.94
3 16.17 7.63 7.34 36 19.38 -16.62 46.17
4 25.00 9.78 13.82 13 35.62 22.62 174
样点4 0 8.89
1 11.58 13.07 13.00 11 11.19 0.19 1.73
2 20.51 12.49 12.22 37 38.31 1.31 3.53
3 30.48 12.13 12.74 39 41.93 2.93 7.50
4 35.20 11.15 11.49 14 20.87 6.87 49.06
本文方法对4个样点中深度在50 cm以内的土壤分层识别效果较好,绝对误差基本都限定在4 cm之内,随着土壤分层深度的增加,误差有所增大,最下面一层的层厚相对误差达到50%,对于样点2的通体黏型土体结构甚至未识别出其位于98 cm左右一层。总体来说,通过包络检波和STFT谱分析的方法进行土壤浅层地表分层厚度识别处理较可靠。

4.4 不确定性因素分析

实际农田环境中,存在着一些不确定性因素,从而为土层厚度测量带来一定误差,这些不确定性因素包括:
(1)农田环境的不确定性。由于农田环境中地表的复杂性和非均一性,例如地表植被、地上覆膜、土壤空洞以及土层盐分或者水分的不均匀等,使得探地雷达回波图像中出现异常的信号突变。数据采集时间为2016年4月,虽然我们选择了裸露地表或者植被稀疏地表,但杂草等地表覆盖物依然不可避免,有的地方有上年的地表残膜存在。另外曲周地区以前存在大量的盐碱地,经过多年的治理,虽然有了很大改善,但个别样点可能依然存在盐度较高的土块。这些都将带来探地雷达回波信号的不稳定。预处理中通过信号滤波降低了一部分由这些突变带来的噪声,但由土壤中较大空洞等造成的回波变化依然会引起其下部土层的层厚测定误差。
(2)土壤层的渐变性。由于土壤质地的变化在农田环境中通常是渐变的,介电常数也因此呈现出渐变的过程,但为了测算方便,本文将每层介电常数设为一常数进行计算,且实际通过介电常数对速度的估计也并非十分准确,土层厚度的计算是根据回波时延和所估算出的速度来确定的,因此速度估算的误差使得土层厚度的测算也存在一定误差。然而,这种局限性和约束性对土层厚度测算的影响在实际应用中经常处于可接受的范围内。
(3)土层深度的影响。电磁波信号随着穿透深度的增加,被散射、被吸收的程度也逐渐增加,因此探测深度达到一定程度后,回波强度减弱而噪声强度增加。因此,在实际农田环境中由于水分和黏土等较强吸收的影响,土壤深度超过70~80 cm后测算结果可能存在较大误差。对于农业应用中比较关心的耕层深度问题,由于一般其深度不超过50 cm,其识别精度是比较可靠的。

5 结论与讨论

5.1 结论

针对探地雷达在土壤分层信息识别上的应用,结合理想模拟实验环境和实际农田环境的探地雷达实测图像,研究了其分层特性及其厚度估算,通过采用包络检波的方法进行分层识别并采用对STFT谱分析得到的土壤层对应介电常数,进而计算土壤分层情况及各层厚度,得出以下主要结论:
(1)基于包络检波和Hilbert谱的时域分析方法,可以得到不同土体结构的时域分层图像,清晰地显示和表征出地下土壤分层特征,虽然可能存在一些检测误差,但总体来说分层层数和现场剖面实测结果基本相符。
(2)基于STFT算法对不同土壤分层进行介电常数及微波传播速度估算,结合基于包络检波的分层识别结果,得到实际农田环境土壤分层情况和层厚信息。对于深度在70 cm以内的土层,其层深计算的相对误差大都在10%以下;对于深度大于70 cm的层位,其层深识别误差则较大。在以后的工作中我们将考虑采用更高功率或更低频率的探地雷达进行数据获取,以提高信号的穿透性。此外,还将考虑尝试其它从弱信号中检测土层信息的方法。

5.2 讨论

本文在位于华北平原的河北曲周进行了实际农田环境实验,验证了本方法在华北平原地区的适用性,在华北平原其它地区使用时,由于雷达数据受下垫植被特征影响较大,因此采样时应注意避开作物覆盖的时间或样点,根据研究区土壤黏度和采样时土壤水分情况合理设置探地雷达参数。对于我国其它区域,由于土壤质地差异较大,雷达信号的穿透性和受干扰情况不尽相同,具体扩展到其它地区的实验将在下一步工作中进行。
[1]
解迎革, 李霞, 张风宝 , 等. 基于电阻率断层扫描技术探测林地土层厚度[J]. 农业工程学报, 2015,31(4):212-216.

[ Xie Y G, Li X, Zhang F B , et al. Detection of soil thickness in forest land based on electrical resistivity tomographic scanning technology[J]. Transactions of the Chinese Society of Agricultural Engineering, 2015,31(4):212-216. ]

[2]
郭秀军, 章光新, 谭笑平 . 物理探查方法在土壤改良中的应用研究[J]. 地理科学, 1999,19(5):470-474.

[ Guo X J, Zhang G X, Tan X P . Study on geophysical survey methods used in soil Improvement[J]. Scientia Geographica Sinica, 1999,19(5):470-474. ]

[3]
Power C, Gerhard J, Tsourlos P , et al. Improved time-lapse electrical resistivity tomography monitoring of dense non-aqueous phase liquids with surface-horizontal boreholes arrays[J]. Journal of applied geophysics, 2015,112(1):1-13.

[4]
何瑞珍, 胡振琪, 王金 , 等. 利用探地雷达检测土壤质量的研究进展[J]. 地球物理学进展, 2009,24(4):1483-1492.

[ He R Z, Hu Z Q, Wang J , et al. The progress of using ground penetrating to detect the soil quality[J]. Progress in Geophys, 2009,24(4):1483-1492. ]

[5]
胡振琪, 陈政, 陈星彤 . 应用探地雷达检测复垦土壤的分层结构[J]. 中国矿业, 2005,14(3):73-75.

[ Hu Z Q, Chen B Z, Chen X T . Study on layer structure of rehabilitated soil using ground penetrating radar[J]. China Mining Magazine, 2005,14(3):73-75. ]

[6]
彭亮, 徐清, 朱忠礼 , 等. 应用低频微波波段GPR测量土壤结构[J]. 北京师范大学学报, 2007,43(3):324-329.

[ Peng L, Xu Q, Zhu Z L , et al. Mapping the soil layers and texture with ground penetrating radar[J]. Journal of Beijing Normal University (Natural Science), 2007,43(3):324-329. ]

[7]
于秀秀, 马兴旺, 迪力夏提 , 等. 探地雷达在土层厚度调查中的试验研究[J]. 土壤学报, 2011,48(4):874-878.

[ Yu X X, Ma X W, Di L X T , et al. Using ground penetrating radar in determination of soil depth[J]. ActaPedologicaSinica, 2011,48(4):874-878. ]

[8]
Kataev S G . An approach of detection the structure of subsurface layers using reflected GPR signals[J]. Ultrawideband and Ultrashort Impulse Signals, 2008,12(3):176-176.

[9]
Ardekani M R, Druyts P, Lambot S , et al. Recovering the structure of a layered soil, including layer thickness and dielectric permittivity, using the interfaces and objects backscatter detected in GPR B-scans[C]. International Conference on Ground Penetrating Radar. IEEE, 2014: 397-400.

[10]
Léger E, Dafflon B, Soom F , et al. Quantification of arctic soil and permafrost properties using ground-penetrating radar and electrical resistivity tomography datasets[J]. IEEE Journal of Selected Topics in Applied Earth Observations & Remote Sensing, 2017,10(10):4348-4359.

[11]
Liu X, Cui X, Guo L , et al. Non-invasive estimation of root zone soil moisture from coarse root reflections in ground-penetrating radar images[J]. Plant and Soil, 2019,436(1-2):623-639.

[12]
Liu X, Dong X, Xue Q , et al. Ground penetrating radar (GPR) detects fine roots of agricultural crops in the field[J]. Plant & Soil, 2018: 1-15.

[13]
刘军, 刘敏, 程建川 . 国外探地雷达技术在测定沥青层厚度方面的应用[J]. 中外公路, 2006(2):102-105.

[ Liu J, Liu M, Chen J C , et al. The application of the ground penetrating radar in the measuring asphalt thickness[J]. Journal of China and Foreign Highway, 2006(2):102-105. ]

[14]
宗鑫, 王心源, 刘传胜 , 等. 探地雷达在地下考古遗存探测中的实验与应用[J]. 地球信息科学学报, 2016,18(2):272-281.

[ Zong X, Wang X Y, Liu C S , et al. Experiments and applications of ground penetrating radar in the investigation of subsurface archaeological interest[J]. Journal of Geo-information Science, 2016,18(2):272-281. ]

[15]
杨虎, 郭华东, 李新武 , 等. 极化雷达目标信息分解技术及其在古湖岸线探测中的应用[J]. 地球信息科学学报, 2003,5(2):109-114.

[ Yang H, Guo H D, Li X W , et al. Target decomposition using polarimetric imaging radar data and its application in identification of bank lines[J]. Journal of Geo-information Science, 2003,5(2):109-114. ]

[16]
张海如, 王国富, 张法全 . 联合瞬时参数分析用于探地雷达目标增强[J]. 微波学报, 2017,33(1):16-20.

[ Zhang H R, Wang G F, Zhang F Q . Conjoint instantaneous parameters analysis for GPR targets enhancement[J]. Journal of Mircrowaves, 2017,33(1):16-20. ]

[17]
张先武, 高云泽, 方广有 . Hilbert谱分析在探地雷达薄层识别中的应用[J]. 地球物理学报, 2013,56(8):2790-2798.

[ Zhang X W, Gao Z Y, Fang G Y . Applicaton of Hilbert spectrum analysis in ground penetrating radar thin layer recognition[J]. Chinese Journal of Geophysics. 2013,56(8):2790-2798. ]

[18]
朱军涛, 廖红建, 谢勇勇 . 采用短时傅里叶变换的铁路车载探地雷达数据解译方法[J]. 西安交通大学学报, 2012,46(7):108-114.

[ Zhu J T, Liao H J, Xie Y Y . Data interpretation of ground penetrating radar (GPR) via short-time Fourier transform for railway track detection[J]. Journal of Xi’an Jiaotong University, 2012,46(7):108-114. ]

[19]
戴前伟, 吴铠均, 张彬 . 短时傅里叶变换在GPR数据解释中的应用[J]. 物探与化探, 2016,40(6):1227-1231.

[ Dai Q W, Wu K J, Zhang B . A study of application of short-time Fourier transform to GPR data interpretation[J]. Geophysical and Geochemical Exploration, 2016,40(6):1227-1231. ]

[20]
于景兰, 王春和 . 探地雷达探测地下目标时的波速估计[J]. 地球物理学进展, 2003,18(3):477-480.

[ Yu J L, Wang C H . Estimation of EM wave velocity in detecting underground target by GPR[J]. Progress in Geophysics, 2003,18(3):477-480. ]

[21]
王新静, 赵艳玲, 胡振琪 . 不同水分条件下探地雷达电磁波波速估算方法与对比分析[J]. 煤炭学报, 2013,38(S1):174-179.

[ Wang X J, Zhao Y L, Hu Z Q , et al. Estimation method and comparative analysis of ground penetrating radar electromagnetic wave velocity based on different water contents[J]. Journal of China Coal Society, 2013,38(Suppl 1):174-179. ]

文章导航

/