CN109343022A - 估测层间土壤含水量的方法 - Google Patents

估测层间土壤含水量的方法 Download PDF

Info

Publication number
CN109343022A
CN109343022A CN201811485381.1A CN201811485381A CN109343022A CN 109343022 A CN109343022 A CN 109343022A CN 201811485381 A CN201811485381 A CN 201811485381A CN 109343022 A CN109343022 A CN 109343022A
Authority
CN
China
Prior art keywords
wave
velocity
interlayer
soil
radar
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201811485381.1A
Other languages
English (en)
Other versions
CN109343022B (zh
Inventor
崔喜红
刘新波
曹鑫
陈学泓
陈晋
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Normal University
Original Assignee
Beijing Normal University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Normal University filed Critical Beijing Normal University
Priority to CN201811485381.1A priority Critical patent/CN109343022B/zh
Publication of CN109343022A publication Critical patent/CN109343022A/zh
Application granted granted Critical
Publication of CN109343022B publication Critical patent/CN109343022B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/12—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation operating with electromagnetic waves
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G01S7/414—Discriminating targets with respect to background clutter
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/885—Radar or analogous systems specially adapted for specific applications for ground probing
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G01S7/411—Identification of targets based on measurements of radar reflectivity
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G01S7/418—Theoretical aspects
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/38—Processing data, e.g. for analysis, for interpretation, for correction
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30—Assessment of water resources

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Electromagnetism (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Geology (AREA)
  • Geophysics (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geophysics And Detection Of Objects (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明提供了一种估测层间土壤含水量的方法。所述方法包括如下步骤：利用固定收发天线间隔的探地雷达测量方式获取粗根反射数据；利用基于RHT的双曲线自动识别算法从粗根双曲线信号中提取层间波速；将所得到的层间波速转化为层间土壤含水量。该方法首次将植物根系作为反射体估测土壤含水量，实现了在未知粗根深度信息的情况下，准确获取层间波速。这为快速、无损获取土壤含水量垂直变化情况提供了一种行之有效的方法，对研究植物根系与土壤之间相互作用关系具有重要意义。

Description

估测层间土壤含水量的方法
技术领域
本发明涉及探地雷达图像的处理方法，特别是利用基于RHT的双曲线自动识别算法，在雷达图像上自动识别双曲线并获取平均波速的过程，属于图像处理领域。
背景技术
土壤含水量(soil water content，SWC)是土壤的重要属性之一，是理解自然生态系统和生物多样性形成的关键变量。土壤含水量也是水循环中的一个重要环节，它通过控制着地面与大气之间能量、水的交换，影响土壤水渗透，渗流以及再分配等水文过程。另外，在农业的精细化管理，高速公路维护以及大气水文模式研究中，土壤含水量都是一个关键因子。
由于土壤含水量具有较高的时间变化和空间变异的特点，对其在中尺度上快速、无损、准确测量一直是一个巨大挑战。探地雷达作为一种无损的地球物理探测工具，在过去的二十年已经被证明有能力在中尺度上获取土壤含水量的分布。目前，经常被使用的方法主要有两种：(1)基于雷达地面波分析的地面波法；(2)基于雷达反射波分析的反射波法。
雷达地面波是传播在雷达发射天线和接收天线之间土壤表层的电磁波，能够反映土壤表层土壤含水量信息。在实地测量应用中，地面波法是目前最被广泛使用的方法之一。然而，地面波法能够获取的土壤含水量深度有限，尤其对于地面波影响深度的定义目前还没有统一，参见Galagedara L.W.,Parkin G.W.,Redman J.D.et al.Field studies ofthe GPR ground wave method for estimating soil water content duringirrigation and drainage[J].J HYDROL.2005,301(1-4):182-197.(《灌溉和排水过程中利用GPR地面波方法估测土壤含水量的实地研究》，Galagedara L.W.等，Journal ofHydrology杂志，2005年)，这不利于地面波法所获取的土壤含水量数据被进一步推广使用。另外，地面波易受杂波干扰，在雷达图像中不容易识别，这也限制了该方法的测量精度和使用范围。相比而言，反射波法操作简单、稳定性好，能够根据反射体位置获取不同深度范围(这里指反射体与地表面之间)的平均土壤含水量，许多研究已经证明了反射波法在实地准确测量和监测土壤含水量变化的能力。
反射波法能够被利用反演土壤含水量主要是基于雷达反射波波速在土壤中对土壤水分非常敏感这一特点。因为土壤一般被认为是由水、空气和固相组成的三相混合物，在非饱和土壤中，水的介电常数为80，土壤固相介电常数为3～10，空气介电常数为1。所以水的介电常数远大于土壤其它属相介电常数，土壤含水量也就决定了土壤的介电常数大小，参见Huisman J.A.,Hubbard S.S.,Redman J.D.et al.Measuring soil water contentwith ground penetrating radar[J].VADOSE ZONE J.2003,2(4):476-491.(《探地雷达测量土壤含水量》，Huisman J.A.等，VADOSE ZONE杂志，2003年)。在土壤介质中，雷达波速v与土壤介电常数ε近似存在关系：
这里c代表电磁波在真空中的传播速度(3×108m·s-1)，参见Davis J.L.,AnnanA.P..Ground penetrating radar for high-resolution mapping of soil and rockstratigraphy.Geophysical Prospecting[J].Geophysical Prospecting,1989,37(5):531-551(《探地雷达用于土壤和岩石地层的高分辨率测绘》，Davis J.L.，地球物理勘探杂志，1989年)。因此，如果已知雷达反射波波速，就可以根据上述公式计算出土壤介电常数，然后通过建立土壤介电常数与含水量之间的关系，最终得到土壤含水量。目前，有许多经验性的土壤介电常数与含水量关系可以被直接用来计算土壤含水量，正如使用最广的Topp公式，参见Topp G.C.,Davis J.L.,Annan A.P.Electromagnetic determination of soilwater content:measurements in coaxial transmission lines[J].WATER RESOURRES.1980,16(3):574-582(《土壤含水量的电磁测定：同轴传输线的测量》，Topp G.C.等，水资源研究杂志，1980年)：
θ＝-5.3×10-2+2.92×10-2ε-5.5×10-4ε2+4.3×10-6ε3
为了计算探地雷达反射波波速，首先需要获取反射波信号数据。目前，反射波信号可以通过固定天线间隔(FO，fixed-offset)和共中心点(CMP，common mid-point)两种测量方式获取。其中，CMP测量是以侧线上某一点为中心，按照一定间隔依次增大发射天线和接收天线之间的距离来获取反射波数据，这种测量方式要求接收天线和发射天线可分离的雷达系统，且测量耗时、数据时间和空间分辨率比较低；FO测量是以固定的天线间距，沿侧线拖动雷达获取数据，该测量方式方便快捷，具有较好的时间分辨率和空间分辨率，容易通过汽车等移动平台进行大范围获取雷达数据，非常适合中尺度土壤含水量测量，并且，Huisman等(《探地雷达测量土壤含水量》，Huisman J.A.等，VADOSE ZONE杂志，2003年)认为使用FO测量的反射波法(后面简称为FO反射波法)，将有潜力在土壤异质性强的实地条件下获取土壤含水量的空间变化。
但是，使用FO反射波法的首要条件是选择合适的反射体，并且反射体的深度信息通常是快速计算平均波速的必要条件，这限制了该方法的推广使用。就目前研究而言，土壤中连续的反射界面(如岩性转换界面、潜水面)，或者被埋藏的人造反射体(如管线、铝板)都可以作为反射体来估计平均土壤含水量。反射体的深度可以通过土钻、已知的土壤剖面信息或埋藏反射体时提前测量获得，但这些方法费时费力，会破坏测量点且仅能提供有限的深度和空间范围的信息，这不利于在中尺度上快速、重复地测量土壤含水量。另外，目前反射波法所获取的土壤含水量代表的是地面与反射体之间的平均含水量，对于反射体以上不同土壤层的含水量变化，还不能进行有效估计。因此，为了能够无损、快速和重复获取平均土壤含水量和土壤层间土壤含水量信息，有必要选择合适的反射体，以及克服需要提前已知反射体深度信息才可以计算平均波速这一局限。
近二十年，探地雷达已经在实地成功用于植物粗根的无损测量，主要包括植物根系形态绘图，根系三维结构恢复，以及根茎和生物量等参数的定量估计。这些研究充分证明探地雷达可以有效探测不同深度的粗根。
发明内容
为此本发明提供了一种利用探地雷达粗根反射估测层间土壤含水量的方法，该方法可以解决前面所提到的问题。
为解决以上问题，本发明提供的估测层间土壤含水量方法，首先提出使用植物粗根作为反射体估测土壤含水量，其理由主要有四点：(1)不同深度的粗根可以被探地雷达所探测，且信号在雷达图像中可以识别；(2)植物不同深度的侧根一般以辐射状的方式分布在植株的周围，这一生物特点有助于估计土壤层间土壤含水量；(3)粗根具有特定双曲线形状的雷达信号，这为从信号形状中自动提取波速提供有利条件；(4)由粗根反射信号估测出的土壤含水量，代表粗根以上的土壤含水量情况，这对于研究植物根系与土壤相互作用关系具有重要意义。
本发明方法主要包括以下步骤：
A、利用FO(固定天线间隔，fixed-offset)探地雷达测量方式获取粗根反射数据：该步骤使用探地雷达FO测量方式获取粗根反射数据，雷达设备可为收发天线分置装置，也可为收发天线一体的便携式雷达。
B、利用基于RHT的双曲线自动识别算法从粗根双曲线信号中提取平均波速(有关基于RHT的双曲线自动识别算法可以参见Xu,L.,Oja,E..Randomized hough transform(RHT):basic mechanisms,algorithms,and computational complexities[J].ComputerVision&Image Understanding.1993,57(2),131-154(《随机Hough变换(RHT)：基本机制，算法和计算复杂性》，Xu,L.等，计算机视觉和图像理解杂志，1993年)和A.Simi,S.Bracciali,G.Manacorda.Hough transform based automatic pipe detection for array GPR:Algorithm development and on-site tests[C].2008.(《基于阵列GPR自动探测管道的Hough变换：算法发展和实地验证》，Simi A.等，IEEE雷达会议，2008年))：该步骤利用基于RHT的双曲线自动识别算法，在雷达数据图像中自动识别粗根双曲线信号，并通过所建立的双曲线方程，求解出粗根反射波在土壤中的平均波速等参数，所述双曲线方程形式为：
其中x表示雷达测量位置与粗根之间的水平距离，d表示粗根的深度，tw表示雷达测量位置所接收到的信号的双程走时，vsoil表示雷达电磁波在粗根以上土壤中的平均波速，代表地表面和粗根之间的平均速度。
C、基于步骤B获取的平均波速计算层间波速：该步骤针对深度不同的相邻两粗根，根据步骤B得到它们各自的平均波速和对应深度的信号双程走时，计算得到两根之间的层间波速，层间波速计算公式如下：
其中和t12表示深度深的粗根的平均波速和信号双程走时，和t1表示深度浅的粗根的平均波速和信号双程走时，表示两粗根之间的层间波速。
D、将所得到的层间波速转化为土壤含水量：该步骤首先根据步骤B和步骤C获得的层间波速，利用波速与介电常数之间的关系，将它们转化为对应的土壤介电常数，波速v与介电常数ε之间的关系如下：
其中，c代表电磁波在真空中的传播波速，值为3×108m·s-1；v为层间波速；ε代表土壤介电常数。然后，进一步利用介电常数与土壤含水量的经验关系(Topp公式)，计算得到层间土壤含水量，介电常数与土壤含水量的关系如下：
θ＝-5.3×10-2+2.92×10-2ε-5.5×10-4ε2+4.3×10-6ε3
其中，θ为土壤层间含水量。
本发明具有的有益效果是：
本发明方法可以利用土壤中天然存在的粗根作为反射体，在未知反射体深度的情况下，使用基于RHT的双曲线自动识别算法从粗根反射信号中自动提取粗根以上平均土壤含水量，以及不同土壤层含水量的变化情况，可以达到在较大范围内快速、无损、重复获取土壤含水信息的目的。
本发明的其他特征和优点将在随后的说明书中阐述，并且，部分的从说明书中变得显而易见，或者通过实施本发明而了解。本发明的目的和其他优点可通过在所写的说明书、权利要求书、以及附图中所特别指出的结构来实现和获得。
附图说明
以下附图仅旨在于对本发明做示意性说明和解释，并不限定本发明的范围。其中，
图1为根据本发明的一个具体实施例的利用探地雷达粗根反射数据估测平均土壤含水量和层间土壤含水量方法的流程示意图；
图2(a)为根据本发明提供的利用探地雷达粗根反射数据估测平均土壤含水量和层间土壤含水量方法中，在步骤A中使用探地雷达FO测量方式对粗根进行测量示意图；图2(b)显示的是图2(a)测量获得的信号示意图；
图3为根据本发明提供的利用探地雷达粗根反射数据估测平均土壤含水量和层间土壤含水量方法，在步骤B中一个从粗根双曲线信号提取平均波速的具体过程实例：(a)探地雷达FO测量方式采集的原始图像，(b)预处理后的雷达图像，(c)边缘提取后二值化显示的边缘图像，(d)使用随机Hough变换在粗根附近识别的双曲线，(e)确定最终的目标双曲线和平均波速；
图4为根据本发明提供的利用探地雷达粗根反射数据估测平均土壤含水量和层间土壤含水量方法，在步骤C中计算层间波速原理示意图；
图5(a)为根据本发明所设计的控制实验的实地场景，图5(b)为实地场景的示意图；
图6(a)为根据本发明提供的利用探地雷达粗根反射数据估测平均土壤含水量和层间土壤含水量方法，在验证实验中所获取的粗根雷达图像，以及被基于RHT双曲线自动识别算法识别的双曲线信号和对应各粗根以上的平均波速，图6(b)为平均土壤含水量与土钻数据结果对比；
图7(a)和图7(b)分别显示的为根据本发明提供的利用探地雷达粗根反射估测平均土壤含水量和层间土壤含水量方法，在实地控制实验中所获取的平均土壤含水量和层间土壤含水量与土钻数据随深度变化趋势对比。
具体实施方式
为了对本发明的技术特征、目的和效果有更加清楚的理解，现说明本发明的具体实施方式。但本领域的技术人员应该知道，以下实施例并不是对本发明技术方案作的唯一限定，凡是在本发明技术方案精神实质下所做的任何等同变换或改动，均应视为属于本发明的保护范围。
图1为根据本发明的一个具体实施例的利用探地雷达粗根反射数据估测平均土壤含水量和层间土壤含水量方法的流程示意图；参照图1所示，下面详细说明根据本发明提供的利用探地雷达粗根反射数据估测平均土壤含水量和层间土壤含水量方法的原理，所述方法包括如下四个步骤：
步骤A，利用FO(固定天线间隔，fixed-offset)雷达测量方式获取粗根反射数据；
步骤B，利用基于RHT的双曲线自动识别算法从粗根双曲线信号中提取平均波速(有关基于RHT的双曲线自动识别算法可以参见Xu,L.,Oja,E..Randomized houghtransform(RHT):basic mechanisms,algorithms,and computational complexities[J].Computer Vision&Image Understanding.1993,57(2),131-154(《随机Hough变换(RHT)：基本机制，算法和计算复杂性》，Xu,L.等，计算机视觉和图像理解杂志，1993年)和A.Simi,S.Bracciali,G.Manacorda.Hough transform based automatic pipe detection forarray GPR:Algorithm development and on-site tests[C].2008.(《基于阵列GPR自动探测管道的Hough变换：算法发展和实地验证》，Simi A.等，IEEE雷达会议，2008年))；
步骤C，基于步骤B获取的平均波速计算层间波速；
步骤D，将所得到的平均波速和层间波速转化为土壤含水量。
对于上述四个步骤，步骤A属于雷达数据收集和获取的过程，需在待测量地点完成；步骤B，C和D属于数据处理过程，可根据实时情况和需求进行安排完成。在步骤A中，本发明选择采用FO测量方式，保证了本发明方法可以方便、快速收集雷达粗根反射数据，并为以后方法在大范围高效地测量土壤含水量提供了可能。在步骤B中，利用基于RHT的双曲线自动识别算法从粗根信号中提取平均波速，这保证了本发明方法可以在未知反射体深度情况下，自动获取平均波速，克服了前人相关方法的局限，既达到了无损测量的目的，也提高了数据解译速度的目的。在步骤C中，本发明可以获取层间波速，这成为本发明的一个重要优势。在步骤D中，本发明简单地采用经验的Topp公式将土壤介电常数转化为土壤含水量，旨在初步测试本发明方法的可行性，对于特殊的测量地点，可采用其它公式或实地定标公式。
下面分别对步骤A，B，C和D予以详细介绍：
步骤A，利用FO雷达测量方式获取粗根反射数据
本发明使用探地雷达FO测量方式获取粗根反射数据，雷达设备可为收发天线分置装置，也可为收发天线一体的便携式雷达系统，测量过程如图2(a)和2(b)所示，其中，所述粗根优选直径大于等于5mm的植物根系，以利于提供更易分辨的反射数据。紧贴地面的雷达收发天线以固定间隔沿着测线移动，由于发射天线对地下发射的电磁波脉冲具有一定的波束宽度，所以在天线未到达根系正上方时，接收天线已经可以接收到粗根的反射信号。这样，随着雷达天线沿侧线移动，天线所获取的粗根反射信号在雷达图像中呈现双曲线形状，并且，根据雷达移动与粗根位置的几何关系(参见图2(b))，双曲线信号方程可表示为如下形式：
其中，vsoil表示雷达电磁波在粗根以上土壤中的平均波速；x表示雷达测量位置与粗根之间的水平距离；tw表示雷达测量位置所接收到的信号的双程走时；d表示粗根的深度；a表示接收天线和发射天线之间的间隔。一般来说，对于使用FO反射波法反演土壤含水量的数据收集过程，为了保证具有较高的分辨率和便捷的测量方式，收发天线间隔a都比较小，可以在上述方程中不被考虑，因此，方程可简化表达为：
由于tw可从图像中直接获取，则该方程中只存在三个未知参数vsoil，x和d。
本发明所提供的方法在步骤B中将采用这种简化的双曲线方程表达形式。
步骤B，利用基于RHT的双曲线自动识别算法从粗根双曲线信号中提取平均波速
本发明在利用基于RHT的双曲线自动识别算法从粗根信号中提取平均波速之前，首先需要对步骤A所获取的雷达原始数据进行预处理操作，主要包括初至时间订正、背景去除、带通滤波和信号增益，效果如图3(b)所示。初至时间订正步骤是将雷达信号记录起始时间调整为信号在地面发生反射的时刻；背景去除、带通滤波和信号增益三个操作具有去除信号杂波和突出双曲线信号在雷达图像中显示的作用。上述预处理操作是探地雷达数据处理的基本操作，可以采用本领域公知的处理方法进行处理。
然后，利用图像边缘检测梯度算子(sober或canny)(canny和sober算子是图像处理中非常常见的处理操作，目前matlab等许多软件都有现成的函数直接可以利用)对雷达灰度图像进行处理，得到二值化(0和1)显示的边缘图像(图3(c))，并将每一个连通的边缘作为一个感兴趣区(即有可能出现双曲线的地方)。
再然后，对每一个感兴趣区内的边缘进行随机Hough变换(RHT，Randomized houghtransform)，如果检测出存在双曲线形状，则根据检测时所建立的目标双曲线方程，返回该方程的三个未知参数vsoil，x和d。随机Hough变换算法与一般曲线拟合的目的相同，都是求解曲线参数，它是基于从变量空间到参数空间的转换过程，通过不断迭代而实现的，在每一次迭代过程，需要从感兴趣区中被提取的边缘上随机选取一组点，根据这组点求解出双曲线方程中的未知参数vsoil，x和d，包括平均波速等，并将参数记录在参数空间中。当完成所有迭代次数后(迭代次数可以自行设定，次数越高，得到结果精度越高，但所花费的时间越多，这里设为10000)，选择参数空间中频次出现最多的一组参数所对应的双曲线为目标双曲线，这样在感兴趣区中既识别了双曲线信号，也获得了平均波速等参数，有关基于RHT的双曲线自动识别算法可以参见Xu,L.,Oja,E..Randomized hough transform(RHT):basicmechanisms,algorithms,and computational complexities[J].Computer Vision&ImageUnderstanding.1993,57(2),131-154(《随机Hough变换(RHT)：基本机制，算法和计算复杂性》，Xu,L.等，计算机视觉和图像理解杂志，1993年)A.Simi,S.Bracciali,G.Manacorda.Hough transform based automatic pipe detection for array GPR:Algorithm development and on-site tests[C].2008.(《基于阵列GPR自动探测管道的Hough变换：算法发展和实地验证》，Simi A.等，IEEE雷达会议，2008年)。
最后，因为对于每一个粗根的双曲线信号，以上步骤会在信号的不同位置自动识别出多条目标双曲线(如图3中d部分所示)，而不同位置双曲线对应的波速大小是不同的，仅仅识别出双曲线是不满足本发明的需要的，因为本发明的目的是为了获得平均波速，因此，本发明采用了进一步的优选步骤：选择位于最亮带与最暗带之间位置(即信号正负最大振幅之间的零点位置)的双曲线作为最终所识别的目标双曲线(如图3中e部分所示)，也就是选择信号最强和最弱位置之间的双曲线作为目标双曲线，其对应方程的波速即为所要求的平均波速，并可进一步用于计算平均土壤含水量和层间波速。
步骤C，基于步骤B获取的平均波速计算层间波速
对于深度不同的相邻两粗根(如图4所示)，根据步骤B的双曲线自动识别算法，可以获取地面与深度较深的粗根之间的平均波速为以及对应粗根深度的信号双程走时t12。同样对于深度较浅的粗根，也可以获取平均波速和信号双程走时分别为和t1。因此，假设两粗根之间的层间波速为对应两粗根垂直间隔的信号双程走时为t2，那么根据两粗根位置的几何关系，可得到数学公式：
和
t2+t1＝t12
然后通过对以上两公式进行化简，得到层间波速计算公式：
这样，只要将步骤B所获取的结果t1和t12代入上式，就可以计算得到两粗根之间的平均波速
步骤D，将所得到的平均波速和层间波速转化为土壤含水量
根据步骤B和步骤C获得的平均波速和层间波速，利用波速与介电常数之间的关系和介电常数与土壤含水量的关系，就可以直接计算得到平均土壤含水量和层间土壤含水量。其中，波速与介电常数之间的关系如下：
介电常数与土壤含水量的关系这里采用Topp公式，正如下所示：
θ＝-5.3×10-2+2.92×10-2ε-5.5×10-4ε2+4.3×10-6ε3
为了更好地说明本发明的技术效果，通过在实地设计一个控制实验，利用本发明提供的方法，从探地雷达粗根反射数据中获取平均土壤含水量和层间土壤含水量，并和相应位置的土钻数据进行对比。控制实验实地场景和设计示意图如图5(a)和5(b)所示。实验时间为2016年7月，地点选择在中国内蒙古阿巴嘎旗(43°55′55″N,114°41′32″E)。实地地点土壤相对均一，含水量较低(通常小于15％)，属于沙质土壤类型(大约80％沙土和约20％粘土)，非常适宜探地雷达测量。在控制实验中，首先挖一个宽1.5米，长8米和深度为1.2米的沙槽(如图5(a))。在沙槽一侧的槽壁上，以水平间隔为1米，在深度分别为10厘米，20厘米，30厘米，40厘米，50厘米，60厘米，70厘米和80厘米位置处，垂直槽壁插入相互平行的八条灌木粗根，并在地表面标记根的位置。这些粗根长度都为50厘米，直径为1.5厘米。然后，将沙槽填埋、平整地面，根据根系位置确定雷达测线位置，使用探地雷达沿侧线垂直根系长轴方向进行测量，获取粗根反射数据。最后，在雷达测量后，沿侧线在每一个粗根位置附近，用两个1米长的土钻分别以不同间隔获取土样，其中一个土钻以20厘米间隔取土样，另一个土钻只取对应粗根深度的土样，这样就获取了不同粗根以上的平均土壤含水量，以及每个粗根附近20厘米间隔的土壤含水量扩线。需要指出的是，雷达数据所获取的土壤含水量是体积含水量，而土钻获取的含水量属于质量含水量，，且体积含水量等于质量含水量与容重的乘积，因此为了将土壤质量含水量转化为体积含水量，在挖沙槽时，同时也采集了不同深度的土壤容重数据。
图6(a)为利用步骤B的双曲线自动识别算法从控制实验雷达粗根数据所获取的平均波速。可以发现，对应八个粗根以上的平均波速，随着深度范围增大，整体存在减小的趋势。将这些平均波速按照步骤D所提供的公式，转化为平均土壤含水量，并和土钻数据比较，结果如图6(b)所示。可以明显发现，本发明所计算出的平土壤含水量与土钻数据十分接近，均方根误差(RMSE)仅为0.0088m3·m-3。
图7(a)为在控制实验中，本发明所获取的平均土壤含水量和土钻数据随深度范围增大的趋势变化。可以明显看出，两组结果变化趋势基本一致，相关系数可达0.939。并且，由土钻数据分析可知，实验地点在深度80厘米以上，平均土壤含水量变化范围较小，变化梯度基本小于0.015m3·m-3，但本发明方法得到了较好结果，基本抓住了这些微小的变化。为了匹配以20厘米间隔采样的土钻数据，选取深度为20，40，60和80厘米处粗根以上的平均波速和对应深度信号双程走时(可由步骤B获取)，代入步骤C的转换公式，则计算出它们之间的层间波速，然后经步骤D得到层间土壤含水量，结果和土钻数据对比如图7(b)所示。可以发现，在控制实验中本发明所计算出的层间土壤含水量略大于土钻数据，RMSE为0.012m3·m-3，仍具有较高的准确度。并且土壤层间含水量扩线与土钻数据变化趋势一致，相关系数可达0.975，很好地反映了实验地点的土壤含水量垂直变化情况。
由上述图实例可看出，本发明所提供的方法可以从探地雷达粗根反射数据中有效估测平均土壤含水量和层间土壤含水量。本发明以自然存在的粗根作为反射体，使用基于RHT的双曲线自动识别算法，能够在未知粗根深度的情况下，准确获取土壤平均含水量和层间土壤含水量，这为快速、无损获取土壤含水量水平分布和垂直变化情况提供了一种行之有效的方法，对研究植物根系与土壤之间相互作用关系具有重要意义。
本领域技术人员应当理解，虽然本发明是按照多个实施例的方式进行描述的，同时在所设计的控制实验中和土钻数据进行了对比验证，但是并非每个实施例仅包含一个独立的技术方案。说明书中如此叙述仅仅是为了清楚起见，本领域技术人员应当将说明书作为一个整体加以理解，并将各实施例中所涉及的技术方案看作是可以相互组合成不同实施例的方式来理解本发明的保护范围。
以上所述仅为本发明示意性的具体实施方式，并非用以限定本发明的范围。任何本领域的技术人员，在不脱离本发明的构思和原则的前提下所作的等同变化、修改与结合，均应属于本发明保护的范围。

Claims (3)

1.一种估测层间土壤含水量的方法，其特征在于，所述方法包括如下步骤：
步骤A，利用固定收发天线间隔的探地雷达测量方式获取粗根反射数据；
该步骤利用固定收发天线间隔的探地雷达测量方式对植物根系进行测量，从而获取粗根反射数据，雷达设备为收发天线分置装置或收发天线一体的便携式雷达；
步骤B，利用基于RHT的双曲线自动识别算法从粗根双曲线信号中提取平均波速；
该步骤根据步骤A获得的粗根反射雷达数据，利用基于随机Hough变换的双曲线自动识别算法，在雷达图像上自动识别双曲线，建立双曲线方程，进而求解出雷达电磁波在土壤中的平均波速的参数，所述双曲线方程形式为：
其中x表示雷达测量位置与粗根之间的水平距离，d表示粗根的深度，tw表示雷达测量位置所接收到的信号的双程走时，vsoil表示雷达电磁波在粗根以上土壤中的平均波速，即地表面和粗根之间的平均速度；
步骤C，基于步骤B获取的平均波速计算层间波速；
该步骤针对深度不同的相邻两粗根，根据步骤B得到两粗根以上的平均波速和对应深度的信号双程走时，计算得到两根之间的层间波速，层间波速计算公式如下：
其中和t12表示深度深的粗根的平均波速和对应深度的信号双程走时，和t1表示深度浅的粗根的平均波速和对应深度的信号双程走时，表示两粗根之间的层间波速；
步骤D，将所得到的层间波速转化为层间土壤含水量；
该步骤根据步骤C获得的层间波速，利用波速与介电常数之间的关系和介电常数与土壤含水量的关系，得到层间土壤含水量，其中，波速v与介电常数ε之间的关系如下：
其中，c代表电磁波在真空中的传播速度，值为3×108m·s-1；v为层间波速；ε代表土壤介电常数；
然后，利用介电常数ε与土壤含水量θ的关系，计算得到层间土壤含水量：
θ＝-5.3×10-2+2.92×10-2ε-5.5×10-4ε2+4.3×10-6ε3；
其中所述步骤B的具体实现方法为：
首先，利用图像边缘检测梯度算子对雷达灰度图像进行处理，得到二值化显示的边缘图像，并将每一个连通的边缘作为一个感兴趣区；
然后，对每一个感兴趣区内的边缘进行随机Hough变换，检测出目标双曲线，并得到该双曲线方程参数，包括平均波速；
最后，对每一个粗根的双曲线信号，上一步骤会在信号的不同位置识别出多条目标双曲线，需要选择位于信号最亮带与最暗带之间位置的双曲线作为最终所识别的目标双曲线，其对应方程的波速参数即为所要求的平均波速。
2.根据权利要求1所述的方法，其特征在于，所述方法还包括在所述步骤B之前对雷达数据进行信号初至时间订正、背景去除、带通滤波和增益的预处理过程。
3.根据权利要求1所述的方法，其特征在于，所述方法还包括在所述步骤C中，使用深度不同的相邻两粗根以上的平均波速和对应深度的信号双程走时，推导层间波速计算公式的过程，具体推导过程为：
根据双曲线自动识别算法，能够获取地面与深度较深的粗根之间的平均波速为以及对应粗根深度的信号双程走时t12；同样对于深度较浅的粗根，可获取平均波速和信号双程走时分别为和t1；因此，假设两粗根之间的层间波速为对应两粗根垂直间隔的信号双程走时为t2，那么根据两粗根位置的几何关系，可得到数学公式：
和
t2+t1＝t12
然后通过对以上两公式进行化简，得到层间波速计算公式：
CN201811485381.1A 2016-12-30 2016-12-30 估测层间土壤含水量的方法 Active CN109343022B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811485381.1A CN109343022B (zh) 2016-12-30 2016-12-30 估测层间土壤含水量的方法

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201611250710.5A CN106772641B (zh) 2016-12-30 2016-12-30 估测平均土壤含水量的方法
CN201811485381.1A CN109343022B (zh) 2016-12-30 2016-12-30 估测层间土壤含水量的方法

Related Parent Applications (1)

Application Number Title Priority Date Filing Date
CN201611250710.5A Division CN106772641B (zh) 2016-12-30 2016-12-30 估测平均土壤含水量的方法

Publications (2)

Publication Number Publication Date
CN109343022A true CN109343022A (zh) 2019-02-15
CN109343022B CN109343022B (zh) 2022-11-15

Family

ID=58927951

Family Applications (2)

Application Number Title Priority Date Filing Date
CN201611250710.5A Expired - Fee Related CN106772641B (zh) 2016-12-30 2016-12-30 估测平均土壤含水量的方法
CN201811485381.1A Active CN109343022B (zh) 2016-12-30 2016-12-30 估测层间土壤含水量的方法

Family Applications Before (1)

Application Number Title Priority Date Filing Date
CN201611250710.5A Expired - Fee Related CN106772641B (zh) 2016-12-30 2016-12-30 估测平均土壤含水量的方法

Country Status (1)

Country Link
CN (2) CN106772641B (zh)

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107831177A (zh) * 2017-09-08 2018-03-23 山东农业大学 一种测定滨海盐渍土壤中水份含量与盐渍化程度的方法
CN110133643A (zh) * 2019-05-22 2019-08-16 北京林业大学 植物根系探测方法及装置
CN110376584A (zh) * 2019-06-28 2019-10-25 浙江大学 基于探地雷达图像特征信号识别的供水管道漏损检测方法
CN111880176A (zh) * 2020-07-24 2020-11-03 安徽理工大学 一种煤矸石充填复垦地土壤含水量的测定方法
CN111929678A (zh) * 2020-07-24 2020-11-13 安徽理工大学 一种浅地表土壤含水量的测定方法
CN111982969A (zh) * 2020-07-02 2020-11-24 中国科学院地理科学与资源研究所 分层土壤的含水量的测量方法及测量设备
CN112816663A (zh) * 2021-02-02 2021-05-18 陆相荣 防洪工程的黄河堤坝土质含水量监测方法及装置
CN113985403A (zh) * 2021-10-26 2022-01-28 中国科学院西北生态环境资源研究院 一种流动沙丘克隆植物根茎生长监测设备
CN116577418A (zh) * 2023-07-14 2023-08-11 山东省科学院激光研究所 一种基于分布式光纤传感的声速测量反演方法

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107727016A (zh) * 2017-10-11 2018-02-23 山东农业大学 一种滴灌灌溉方式下土壤湿润体的探测方法
CN107782640B (zh) * 2017-10-20 2020-07-31 东北大学 一种进行岩石试件含水均匀性检测和扩散系数计算的方法
CN107977720A (zh) * 2017-12-20 2018-05-01 云南电网有限责任公司文山供电局 一种接地极状态分析和异常预警系统
CN108387893A (zh) * 2018-03-13 2018-08-10 曲阜师范大学 一种探地雷达电磁波传播速度的提取方法

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101813775A (zh) * 2010-04-16 2010-08-25 北京师范大学 一种利用探地雷达测量植物地下根的直径以及生物量的方法
CN101887116A (zh) * 2010-04-14 2010-11-17 电子科技大学 一种用于微弱目标的改进型随机Hough变换检测方法
CN101915771A (zh) * 2010-09-03 2010-12-15 中国科学院南京土壤研究所 基于低频探地雷达地波法的测定砂质土壤含水量的方法
CN104062654A (zh) * 2014-07-17 2014-09-24 电子科技大学 一种基于超宽带雷达的土壤含水量测量方法
US8849523B1 (en) * 2013-05-20 2014-09-30 Elwha Llc Systems and methods for detecting soil characteristics
CN104569972A (zh) * 2015-01-20 2015-04-29 中国科学院成都生物研究所 一种植物根系三维构型无损检测方法

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5384715A (en) * 1993-08-27 1995-01-24 The Texas A&M Univeristy System System identification and analysis of subsurface radar signals
CN101832952A (zh) * 2010-05-19 2010-09-15 北京航空航天大学 一种采用探地雷达检测混凝土结构含水率分布的方法
CN103913777B (zh) * 2014-04-28 2016-06-22 中国科学院亚热带农业生态研究所 喀斯特非连续性土壤与岩石分布的地质雷达探测方法
CN104020495B (zh) * 2014-06-24 2015-05-06 中国矿业大学(北京) 一种基于探地雷达的地下管线参数自识别方法
CN105044785B (zh) * 2014-09-19 2017-05-31 中国矿业大学(北京) 一种基于模糊聚类与Hough变换的探地雷达地下管线探测方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101887116A (zh) * 2010-04-14 2010-11-17 电子科技大学 一种用于微弱目标的改进型随机Hough变换检测方法
CN101813775A (zh) * 2010-04-16 2010-08-25 北京师范大学 一种利用探地雷达测量植物地下根的直径以及生物量的方法
CN101915771A (zh) * 2010-09-03 2010-12-15 中国科学院南京土壤研究所 基于低频探地雷达地波法的测定砂质土壤含水量的方法
US8849523B1 (en) * 2013-05-20 2014-09-30 Elwha Llc Systems and methods for detecting soil characteristics
CN104062654A (zh) * 2014-07-17 2014-09-24 电子科技大学 一种基于超宽带雷达的土壤含水量测量方法
CN104569972A (zh) * 2015-01-20 2015-04-29 中国科学院成都生物研究所 一种植物根系三维构型无损检测方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
康学远等: "土壤介电常数-含水量关系模型研究", 《中国农村水利水电》 *
杨鑫: ""基于Hough变换的林木根系探地雷达图像处理与波速估计"", 《中国优秀硕士学位论文全文数据库 农业科技辑》 *
王前锋等: "基于探地雷达快速测定土壤含水量试验研究", 《自然资源学报》 *
郭立等: ""探地雷达应用于植物粗根探测的研究进展"", 《中国科技论文》 *

Cited By (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107831177A (zh) * 2017-09-08 2018-03-23 山东农业大学 一种测定滨海盐渍土壤中水份含量与盐渍化程度的方法
CN110133643A (zh) * 2019-05-22 2019-08-16 北京林业大学 植物根系探测方法及装置
CN110133643B (zh) * 2019-05-22 2021-08-20 北京林业大学 植物根系探测方法及装置
CN110376584A (zh) * 2019-06-28 2019-10-25 浙江大学 基于探地雷达图像特征信号识别的供水管道漏损检测方法
CN110376584B (zh) * 2019-06-28 2021-09-14 浙江大学 基于探地雷达图像特征信号识别的供水管道漏损检测方法
CN111982969B (zh) * 2020-07-02 2022-02-18 中国科学院地理科学与资源研究所 分层土壤的含水量的测量方法及测量设备
CN111982969A (zh) * 2020-07-02 2020-11-24 中国科学院地理科学与资源研究所 分层土壤的含水量的测量方法及测量设备
CN111880176A (zh) * 2020-07-24 2020-11-03 安徽理工大学 一种煤矸石充填复垦地土壤含水量的测定方法
CN111929678A (zh) * 2020-07-24 2020-11-13 安徽理工大学 一种浅地表土壤含水量的测定方法
CN112816663A (zh) * 2021-02-02 2021-05-18 陆相荣 防洪工程的黄河堤坝土质含水量监测方法及装置
CN112816663B (zh) * 2021-02-02 2023-01-10 陆相荣 防洪工程的黄河堤坝土质含水量监测方法及装置
CN113985403A (zh) * 2021-10-26 2022-01-28 中国科学院西北生态环境资源研究院 一种流动沙丘克隆植物根茎生长监测设备
CN113985403B (zh) * 2021-10-26 2022-08-12 中国科学院西北生态环境资源研究院 一种流动沙丘克隆植物根茎生长监测设备
CN116577418A (zh) * 2023-07-14 2023-08-11 山东省科学院激光研究所 一种基于分布式光纤传感的声速测量反演方法
CN116577418B (zh) * 2023-07-14 2023-09-22 山东省科学院激光研究所 一种基于分布式光纤传感的声速测量反演方法

Also Published As

Publication number Publication date
CN106772641A (zh) 2017-05-31
CN106772641B (zh) 2019-01-25
CN109343022B (zh) 2022-11-15

Similar Documents

Publication Publication Date Title
CN106772641B (zh) 估测平均土壤含水量的方法
Liu et al. Detection and localization of rebar in concrete by deep learning using ground penetrating radar
Guo et al. Subsurface lateral preferential flow network revealed by time‐lapse ground‐penetrating radar in a hillslope
Liu et al. Measurement of soil water content using ground-penetrating radar: A review of current methods
Cui et al. Modeling tree root diameter and biomass by ground-penetrating radar
Zenone et al. Preliminary use of ground-penetrating radar and electrical resistivity tomography to study tree roots in pine forests and poplar plantations
Lunt et al. Soil moisture content estimation using ground-penetrating radar reflection data
CN104020495B (zh) 一种基于探地雷达的地下管线参数自识别方法
CN104471443B (zh) 使用近场与远场超低频和极低频干涉测量合成孔径雷达以用于地下成像的方法
Cui et al. Estimating tree-root biomass in different depths using ground-penetrating radar: Evidence from a controlled experiment
CN109031290A (zh) 一种土地工程的质量检测方法
US20150051883A1 (en) Process to identify and classify oil seep areas at the seabed through inverse modeling
Angelis et al. Combined application of GPR and ERT for the assessment of a wall structure at the Heptapyrgion fortress (Thessaloniki, Greece)
Lorenzo Cimadevila et al. Forestry applications of ground-penetrating radar
Merkle et al. Fusion of ground penetrating radar and laser scanning for infrastructure mapping
CN111551927A (zh) 一种基于三维探地雷达的地下管道管径测量方法
Molon et al. Quantitative, nondestructive estimates of coarse root biomass in a temperate pine forest using 3‐D ground‐penetrating radar (GPR)
JP5629840B1 (ja) 空洞厚探査方法
Cui et al. Pairing dual‐frequency GPR in summer and winter enhances the detection and mapping of coarse roots in the semi‐arid shrubland in China
Cui et al. Chinese radioglaciological studies on the Antarctic ice sheet: Progress and prospects
Wang et al. Qualitative research: The impact of root orientation on coarse roots detection using ground-penetrating radar (GPR)
CN113468804B (zh) 一种基于矩阵束和深度神经网络的地下管道识别方法
Illawathure et al. Evaluating soil moisture estimation from ground‐penetrating radar hyperbola fitting with respect to a systematic time‐domain reflectometry data collection in a boreal podzolic agricultural field
Zhou et al. A novel reconstruction method based on changes in phase for subsurface large sloped dielectric target using GPR
Li et al. Identification of geo-bodies in borehole radar image based on Curvelet transform

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant