基于共中心点道集约束的探地雷达波阻抗反演_戴前伟.pdf
第 48 卷 第 3 期 煤田地质与勘探 Vol. 48 No.3 2020 年 6 月 COAL GEOLOGY 2. Key Laboratory of Metallogenic Prediction of Nonferrous Metal and Geological Environment Monitoring, Ministry of Education, Central South University, Changsha 410083, China Abstract GPR impedance inversion is an effective to obtain accurately the intrinsic parameters of the subsurface medium. This relies on low frequency ination provided by logging data, but it is rarely ac- companied by drilling in practical applications. To solve the problem, the common midpoint velocity analysis is employed to provide more low frequency component ination for impedance inversion of GPR to finely recon- struct the dielectric parameter of subsurface in the framework of velocity-constrained inversion technique. Firstly, a layer model is specifically set up as an example to verify the feasibility of the developed velocity- con- strained inversion to regulate the initial model. Then, the impedance inversion test is pered on two random media models, result of mean relative error from the true model is 8.73, which show that the overall structure is more consistent with the real model, and more importantly, the microstructure is also finely depicted. The proposed is more efficient and economical in impedance inversion of ground penetrating radar with random medium model, with more detailed ination contained in imaging results, the is very feasible and applicable in the estimation of other physical parameters in soil investigation. Keywords ground penetrating radarGPR; impedance inversion; common midpointCMP; velocity analysis; pa- rameter estimation 准确探测土壤的介电常数对水文地质调查、环境领域、工程领域、考古研究具有重要意义,例如, ChaoXing 212 煤田地质与勘探 第 48 卷 地下水含量、污染物运移扩散、管线渗漏、古墓挖 掘等[1-4]。地球物理勘探是现代探测地下相关参数的 一种先进科学技术,常规的钻孔取心和地球物理测 井方法能够提供详细地质信息,但探测范围仅限井 眼附近的有限区域[5]。探地雷达Ground Penetrating Radar,简称 GPR是通过用高频电磁波来确定地下 埋藏物内部或结构的无损探测技术,电磁波在介质 中传播时,其传播路径、波场强度和波形随所通过 介质的电性及几何形态而变化[6]。在探地雷达探测 时,地下介质的电性参数影响着电磁波的传播,从 而控制了探地雷达的回波响应。因此,获取这些介 质的电性参数信息对探地雷达的研究具有重要价 值。共中心点Common Midpoint,简称 CMP速度 分析是从地面探地雷达多偏移距数据中估计目标大 尺度介电常数速度分布的常用方法,但其成像分 辨率有限,难以刻画目标结构细节信息[7-9]。近年, 利用走时、振幅和相位等雷达波信息的全波形反演 Full Wave Inversion, FWI能够获取高分辨率的 地下物性参数分布[10-12],但 FWI 计算花费的时间太 多[13]。新兴的神经网络、遗传算法、粒子群算法等 完全非线性算法,也因其需要在较庞大的解空间中 寻找最优解,由此导致大量的计算量,而不适用于 大范围探测数据的解释[14-16]。 在地震勘探中,狭义的地震反演是指从有限频 带宽度的地震数据中恢复出宽带的波阻抗参数模 型。R. O. Lindseth[17]开发的阻抗反演方法,其表示 地震资料缺少低频数据,必须通过其他方法补充低 频数据才可行; D. A. Cooke 等[18]详细阐述了地震资 料的广义线性反演方法,拉开了地震阻抗反演的序 幕;随后,周竹生等[19]利用地质、地震和测井资料 进行联合反演,克服了单一线性反演的缺陷;R. J. Ferguson 等[20]基于声波测井技术提出有限带宽阻抗 反演,该方法结合了测井资料纵向的强约束和地震 横向分辨率高的优点,实现了全宽带资料的反演; C. Schmelzbach 等[21]通过地震阻抗反演流程,实现 了探地雷达高分辨率地下含水量的估计,并在数值 模拟及实测数据中都验证了其方法的可行性;李 静 [22]利用阻抗反演方法实现了三维探地雷达复杂 随机层状介质含水量等目标本征属性参数的提取; Li Jing 等[23]利用探地雷达波阻抗反演估算了月球风 化层的相对介电常数,得到了地下详细的结构,为 了解月球地下结构提供了可靠的方法;刘钰[24]采用 基于模型的宽带约束探地雷达阻抗反演方法进行古 墓探测,经过苏州木渎古城的探测工作被验证是一 种有效的方法,但其表示基于模型的反演需要建立 较好的初始模型。初始模型的低频信息是否准确关 乎重建探测目标本征属性参数的成败,地震勘探常 常通过测井资料和叠后速度分析建立准确的初始速 度模型。 由于探地雷达在工程项目中钻孔资料极少, 共中心点速度分析成为了探地雷达效仿地震勘探获 取二维速度场,建立初始模型的常用方法。S. Busch 等[25]利用 GPR 共中心点道集进行速度分析, 将速度 分析结果定义为全波形反演的初始模型,得到了可 靠的反演结果。张彬[26]指出共中心点道集能够较方 便地利用叠加速度谱求取不同深度的平均波速,能 够为偏移处理提供初始速度模型,但在地下横行非 均匀介质中效果不佳。J. H. Bradford 等[27]提出利用 偏移后的共中心点数据做反射层析成像估算雷达波 速方法,提高了解译近地表水分布的能力。 探地雷达波阻抗反演是了解浅地表层状介质物 性参数的有效方法,但是探地雷达工程项目中稀缺 的钻孔资料无法为波阻抗反演建立准确的初始模 型,补充所需的低频信息。为此,笔者开展基于共 中心点道集约束的地面探地雷达波阻抗反演研究, 依据采集的地面雷达数据, 通过 CMP 速度分析技术 获取地层的速度分布, 为波阻抗反演建立初始模型, 补充低频信息,获取地下目标介质的中尺度本征属 性参数。 1 方法原理 1.1 速度分析原理 地下介质速度分布图是 GPR 数据解释的关键, 其中,获取地下介质速度分布通过时深转换可以定 位异常体和反射界面的位置, 也可以对 GPR 剖面做 偏移处理提高探测精度[28]。 CMP 速度分析的前提是 假设地层水平,介质均匀且各向同性,则当偏移距 不大时,反射信号满足双曲线时距方程 2 22 0 2 x tt v 1 式中t 为反射信号的双程走时;x 为偏移距; 0 t为 零偏移距处的双程走时;v为均方根速度。 则正常时差t可以表示为 2 000 2 x ttttt v 2 速度分析的原理是根据沿着反射同相轴方向使相 关系数、叠加能量或相似性系数最大。在野外工作中, x 是已知的, 0 t可通过自激自收方式获取,依据地质条 件给定一系列试验速度 minminmax ,,, i vvvvv 代 替式1的均方根速度,便可得到对应的反射双曲 ChaoXing 第 3 期 戴前伟等 基于共中心点道集约束的探地雷达波阻抗反演 213 线,倘若地层速度在试验速度区间内,那么在这反 射双曲线中必然存在一条双曲线与反射轴重合,使 得叠加的能量在速度谱上有一个最大的峰值。则计 算该双曲线所对应的试验速度为反射波的叠加速度 v。选择合适的振幅能量准则叠加能量是获取最佳 速度谱的关键[29],常用的振幅能量准则有叠加能量 法、互相关振幅能量法和归一化振幅能量法,本次 选取互相关振幅能量法[30]制作速度谱。 互相关振幅 能量法的公式为 1 0 1 1 , NN iijj ij i C tv ftft - 3 式中 0, C tv为互相关函数; ii ft为双程走时 t 第 i 道的振幅值;N 为 CMP 多偏移距道数。 当偏移距远小于地下目标体埋深时,叠加速度 近似等于均方根速度,根据速度谱拾取叠加速度后 可以使用迪克斯dix公式[31]求取层速度, 表达式为 22 0,,0,1,1 0,0,1 nnnn n nn tvtv v tt 4 式中 n v为第 n 层层速度; ,n v为第一层至第n层 的均方根速度; 0,n t为第一层至第n层的时间。 得到各层层速度后,结合各层的双程走时便可 计算得到各层的层厚度, 第n层层厚度计算公式为 int1nnn dvtt 5 式中 n d为第n层层厚度; int v为第n层层速度; n t 为第一层至第n层的时间; 1n t 为第一层至第n1 层的时间。 1.2 探地雷达波阻抗反演原理 探地雷达的响应一般由介电常数、电导率和磁 导率控制,一般情况下假设地下介质是无磁性的, 电导率的影响通常只考虑对电磁波的衰减和损耗, 所以在低损耗的介质中介电常数决定了电磁波波速 0 GPR r c v 。电磁波在地下传播时,当相邻地层存 在阻抗差异则会发生反射,电磁波阻抗可以用相对 介电常数估计[32] 0 r Z Z 6 式中 0 377Z,为自由空间阻抗值; r 为相对 介电常数。 因此,阻抗反演可以通过获得的相对介电常数 转化成其他参数,如采用 Topp,CRIM 公式估计含 水量等[33]。 基于上述假设,自激自收的地面探地雷达反射 系数可以表示为[34] 1 1 rr 1 1 rr nn nn nn n nn ZZ R ZZ 7 1 1 1 exp 2 n ni i ZZR 8 式中 n R为第n层和第n1 层界面的反射系数; n Z 为探地雷达数据的反射系数通过式8递推估计的 第n层阻抗值; rn 为第n层相对介电常数。 全频带阻抗值可以根据以下公式获得 b Z tZ tZt 9 式中 b Zt为初始模型通过低通滤波器得到的约 束阻抗值,即所需的低频分量; Z t为通过雷达 资料计算反射系数递推得到的阻抗值;为权重 因子。 一维有限带宽阻抗方法反演流程步骤[35]可以 表示为 ① 雷达反射波能量在地下传播时会随着时间、 几何扩散、吸收和传输损耗等因素而衰减,所以在 计算反射系数前需要对振幅进行真振幅恢复。如若 子波振幅未知,则可以利用先验信息通过式8估算 振幅补偿因子; ② 计算反射系数前需要通过预测反褶积和偏 移等预处理技术提升 GPR 剖面的分辨率; ③ 通过稀疏脉冲反褶积从步骤 b 雷达数据中 计算反射系数; ④ 利用步骤③的反射系数通过式8获得递推 的阻抗值,利用傅里叶变化将递推的阻抗值转至频 率域,并通过高通滤波器滤去不可靠的低频信息, 保留中高频信息,记为式9中的 Z t; ⑤ 利用测井或 CMP 速度函数通过式6得到 先验阻抗信息 Z t,利用傅里叶变化将 Z t转至 频率域并通过低通滤波器滤去中高频信息,得到 符合该测点一定范围内的低频信息,记为式9中 的 b Zt; ⑥ 可以利用先验阻抗信息通过式9求取, 随后通过式9将整个剖面的 Z t和 b Zt在频率 域中相加得到全频带的阻抗值,最后通过式6求取 相应的相对介电常数。 2 算例分析 为了验证基于共中心点道集约束的探地雷达波 阻抗反演的效果,设置了 3 组试验。试验一为共中 ChaoXing 214 煤田地质与勘探 第 48 卷 心点速度分析作为波阻抗反演的初始模型有效性试 验;试验二利用一维模型正演模拟高信噪比信号, 试验基于共中心点道集约束的 GPR 波阻抗反演对 高信噪比资料的反演效果,并与测井约束的波阻抗 反演对比;试验三以二维随机介质为模型,试验基 于共中心点道集约束的 GPR 阻抗反演在随机介质 中的反演效果。 2.1 速度分析 为了验证 CMP 速度分析程序的有效性,建立 了如图 1a 所示的深 1.5 m,水平距离 2 m 的三层水 平均匀模型。其中,第一层介质的埋深为 0.5 m,其 相对介电常数为 10;第二层介质的埋深为 0.9 m, 其相对介电常数为 15;第三层介质的相对介电常数 为 20;假设三层介质的电导率均为 0.001 S/m。 图 1 层状模型正演速度分析 Fig.1 Velocity analysis experiment of layered model forward 采用时域有限差分算法[36-37]对模型正演计算, 空间步长为 0.1 m,网格数为200 150,激励源为中 心频率 500 MHz 的雷克子波,采样率为 0.02 ns,时 窗长度为 40 ns。零偏移距位于地表测线 1.0 m 处, 每次各向两侧对称地移动 0.01 m,最大偏移距为 1.0 m。对上述模型正演计算结果如图 1b 所示,从 图中看到强能量的地滚波掩盖了其他反射信号。为 了更好地突出有效信息,对雷达剖面的地滚波进行 切除,消除地滚波的影响,并选择合适的增益函数 突出深部的有效信号。数据处理结果如图 1c 所示, 图 1b 中 10 ns 和 40 ns 之间被掩盖的两条双曲线在 数据处理后清晰可见。然后选取一系列试验速度对 数据处理后的雷达剖面进行速度扫描,根据互相关 振幅能量法制作速度谱得到如图 1d 所示两个能量 团,假设叠加速度近似为均方根速度,提取能量团 峰值点对应双程走时和叠加速度。最后通过 dix 公 式和式5即可得到层速度及相应的层厚度。从表 1 所示 CMP 速度分析相对误差分析可知, 地层层速度 的最大相对误差为 6.41 ,层厚度最大相对误差为 4.04 。 通过模拟算例分析可知, 本文所设计的 CMP 速度分析算法反演结果误差小,能够很好地提取水 平地层层速度及相应层厚度。 表 1 速度分析结果的相对误差分析 Table 1 Relative error analysis of velocity analysis results 层号 模型层速度/mns –1 计算层速度/mns–1层速度相对误差/模型层厚度/m 计算层厚度/m 层厚度相对误差/ 1 0.095 0.090 5.21 0.500 0.520 4.04 2 0.078 0.073 6.41 0.400 0.413 3.30 ChaoXing 第 3 期 戴前伟等 基于共中心点道集约束的探地雷达波阻抗反演 215 2.2 一维随机介质模型 为了验证 CMP 资料约束下的探地雷达阻抗反 演方法的有效性,在大尺度相对介电常数背景下加 入随机噪声建立如图 2a 所示的一维随机介质模型, 模型的深度为 2.6 m,分为五层,其相对介电常数分 别为 2.9、7、9、11、13,网格数为 260,空间步长 为 0.01 m,假设所有地层的电导率为 0.001 S/m。 采用时域有限差分算法对该模型进行正演模 拟,采用中心频率为 500 MHz 的雷克子波作为激励 源,采样率为 0.02 ns,时窗长度为 52 ns,采用共偏 移距的方式记录雷达数据。 图 2b 为正演模拟预处理 后的雷达波形数据。图 2c 为稀疏脉冲反褶积迭代 20次计算真振幅数据得到的反射系数。 图2d为CMP 速度分析初始模型、常规测井初始模型和理论模型 相对介电常数的对比; 图 2e 黑线为理论模型相对介 电常数, 蓝线为基于 CMP 为初始模型的阻抗反演结 果,红线为基于测井低频信息为初始模型的阻抗反 演结果。分别将两种方法的相对介电常数反演结果 与理论模型对比,其中,测井约束的阻抗反演估计 相对介电常数结果的相对误差平均值为 6.81 , CMP 约束的阻抗反演估计相对介电常数结果的相 对误差平均值为 8.73 。从分析结果看出,虽然测 井约束的阻抗反演结果比速度分析约束下的阻抗反 演结果与理论模型更为接近,但速度分析约束下的 阻抗反演的细节信息与模型有着相近的趋势,其结 果能够为地下介质细节的刻画提供较好的分辨率; 而且,运用该技术反演,采集成本更经济,处理快 捷。通过算例分析,本文提出的基于共中心点道集 约束的阻抗反演方法能够高效、准确地获取层状随 机介质中尺度参数。 a 模型相对介电常数;b 预处理后仿真数据;c 反射系数;d 两种方法的初始模型与理论模型;e 两种方法的反演估计结果与理论模型 图 2 两种探地雷达波阻抗反演方法流程 Fig.2 Process of two wave impedance inversion s of GPR 2.3 二维随机介质模型 建立倾斜层状随机介质模型,如图 3a 所示,模 型大小为 2.30 m1.55 m,设定所有地层电导率均 为 0.001 S/m。 采用时间域有限差分法对模型进行正 演模拟,离散步长为 0.01 m,网格数为230 155, 激励源采用中心频率为 500 MHz 雷克子波,采样率 为 0.02 ns,时窗为 40 ns,采用共偏移距的方式记录 雷达数据。 图 3a 中模型的相对介电常数值呈随机分 布,整体上分成三层。图 3b 为三层随机介质模型有 限差分正演模拟及预处理后的雷达剖面,从图 3 中 可以直观地反映探测区域的地下结构分布,地层分 为三层且层位清晰,在测区内沿层位开展速度分析 工作,本次速度分析零偏移位置于测线 0.8 m 处, 每次各向两侧对称地移动 0.02 m,最大偏移距 0.8 m。截取有效深度的信号拾取地下介质速度, 得到三层介质的层速度分别为 0.116、0.079 7、 0.063 5 m/ns,相应的反射界面埋深分别为 48.7、 70.2 cm,层厚度的最大相对误差为 8.2。分析结果 说明CMP 拾取速度的误差在允许范围内,拾取的 地层层速度和层厚度是准确的,满足作为初始模型 的精度要求。 随后利用 CMP 获取的速度沿层位插值 得到二维速度分布,结合该雷达资料依照上述介绍 的步骤进行阻抗反演。根据反演流程获得如图 4 所 示的反演结果,可以从图 4 观察到反演估计结果整 ChaoXing 216 煤田地质与勘探 第 48 卷 体趋势继承了速度分析初始模型的准确性, 细节信息 通过雷达资料得到了补充, 剖面的分辨率得到了一定 的提升。 综上所述, 该方法在横向连续性较好的层状 介质中有着不错的反演效果。 图 3 三层随机介质模型及有限差分正演结果 Fig.3 Three layer random soil medium model and FDTD forward modeling result 图 4 三层随机介质模型阻抗反演结果 Fig.4 Impedance inversion results of a three layer random medium model 3 结 论 a. 层状模型测试结果显示,与实际测试模型的 参数相比,共中心点道集速度分析的上下层速度的 相对误差分别为 5.21 和 6.41 ,上下层厚度的相 对误差分别为 4.04 和 3.30 ,表明共中心点道集 速度分析能较准确地获取地下大尺度的速度结构信 息,证明了 CMP 速度分析作为 GPR 波阻抗反演初 始模型的策略可行性。 b. 通过随机介质模型的 GPR波阻抗反演测试, 获取了随机介质模型的相对介电常数信息,结果与 测井约束的阻抗效果一致, 为随机介质中 GPR 波阻 抗反演提供了一种经济高效的大尺度约束方法。 c. 复杂随机介质模型测试结果显示, 相对介电 常数的相对误差值为 8.73 , 随机介质细微结构得 到了较好的重构, 证明了基于共中心点道集速度分 析的 GPR 波阻抗反演方法的精确性,在复杂的土 壤随机介质物理参数估计和反演中极具可行性和 适用性。 d..在实际工作中,地层结构复杂多变,如何利 用 CMP 速度分析建立准确的二维速度初始模型, 以期提高反演效果,有待研究。 致谢感谢吉林大学地球探测科学与技术学院李静 副教授在阻抗反演理论及计算上给予的帮助和讨 论;同时感谢两位匿名审稿专家提出的宝贵意见。 在此一并致谢 请听作者语音介绍创新技术成果 等信息,欢迎与作者进行交流 参考文献References OSID 码 [1] GREAVES J R,LESMES P D,LEE M J,et al. Velocity variations and water content estimated from multi offset,ground penetrating radar[J]. Geophysics,1996,613683–695. [2] 王仙丽. 地下油类污染区探地雷达GPR探测能力研究[D]. 青岛中国海洋大学,2011. WANG Xianli. Research on detectability of ground penetrating radar for oil contaminated site[D]. QingdaoOcean University of China,2011. [3] 金鑫. 管线渗漏异常探地雷达数据的电场分量成像分析[J]. ChaoXing 第 3 期 戴前伟等 基于共中心点道集约束的探地雷达波阻抗反演 217 煤田地质与勘探,2018,462159–163. JIN Xin. Study on electric field component imaging of leakage in pipeline using GPR[J]. Coal Geology Exploration,2018, 462159–163. [4] 吴秋霜,王齐仁,皮海康. 水泥混凝土路面脱空的探地雷达图 像特征分析[J]. 煤田地质与勘探,2018,464181–185. WU Qiushuang,WANG Qiren,PI Haikang. Analysis on the image features of ground penetrating radar for cavity of concrete pavement[J]. Coal Geology Exploration,2018,464 181–185. [5] KOBR M,MAREŠ S,PAILLET F. Geophysical well logging Borehole geophysics for hydrogeological studiesPrinciples and applications in hydrogeophysics[M]. BerlinSpringer,2005. [6] 李大心. 探地雷达方法与应用[M]. 北京地质出版社, 1994. LI Daxin. and application of Ground Penetrating Ra- dar[M]. BeijingGeological Publishing House,1994. [7] 刘四新,蔡佳琪,傅磊,等. 利用探地雷达精确探测铁路路基 含水率[J]. 地球物理学进展,2017,322878–884. LIU Sixin,CAI Jiaqi,FU Lei,et al. Accurate detection of moisture content of subgrade by GPR[J]. Progress in Geophys- ics,2017,322878–884. [8] 董泽君,鹿琪,冯晅,等. 探地雷达测量土壤含水量的应用研 究[J]. 地球物理学进展,2017,3252207–2213. DONG Zejun,LU Qi,FENG Xuan,et al. Estimation of soil water content using ground penetrating radar[J]. Progress in Geophysics,2017,3252207–2213. [9] 王洪华,戴前伟. 探地雷达有限元正演及介电参数反演[M]. 长沙中南大学出版社,2016. WANG Honghua, DAI Qianwei. Ground Penetrating Radar finite element numerical simulation and dielectric parameter inver- sion[M]. Changsha Central South University Publishing House, 2016. [10] 林朋, 彭苏萍, 卢勇旭, 等. 基于共轭梯度法的全波形反演[J]. 煤田地质与勘探,2017,451131–136. LIN Peng,PENG Suping,LU Yongxu,et al. Full wave inversion based on the conjugate gradient [J]. Coal Geol- ogy Exploration,2017,451131–136. [11] 孟旭,刘四新,吴俊军,等. 时间域跨孔雷达全波形反演及实 际应用[J]. 世界地质,2016,351256–263. MENG Xu,LIU Sixin,WU Junjun,et al. Full wave in- version of time-domain cross-hole radar and its field applica- tion[J]. Global Geology,2016,351256–263. [12] 冯德山,王珣. 基于 GPU 并行的时间域全波形优化共轭梯度 法快速 GPR 双参数反演[J]. 地球物理学报,2018,6111 4647–4659. FENG Deshan,WANG Xun. Fast ground penetrating radar dou- ble-parameter inversion based on GPU-parallel by time-domain full wave optimization conjugate gradient [J]. Chi- nese Journal of Geophysics,2018,61114647–4659. [13] MELES G A, GREENHALGH S A, GREEN A G, et al. GPR full wave sensitivity and resolution analysis using an FDTD ad- joint [J]. IEEE Transactions on Geoscience Remote Sensing,2012,5051881–1896. [14] 刘敦文,徐国元,黄仁东,等. 一种基于神经网络的探地雷达 信号解释研究[J]. 地球物理学进展,2004,191179–182. LIU Dunwen,XU Guoyuan,HUANG Rendong,et al. Study on a signal interpretation of GPR based on BP neural network[J]. Progress in Geophysics,2004,191179–182. [15] 袁克阔. 粒子群算法改进及内变量本构模型参数反演[J]. 煤 田地质与勘探,2017,452112–117. YUAN Kekuo. Improved particle swarm optimization and pa- rameter inversion in internal variable constitutive model[J]. Coal Geology Exploration,2017,452112–117. [16] 王升,陈洪松,付智勇,等. 基于探地雷达的典型喀斯特坡地 土层厚度估测[J]. 土壤学报,2015,5251024–1030. WANG Sheng, CHEN Hongsong, FU Zhiyong, et al. Estimation of thickness of soil layer on typical karst hillslopes using a ground penetrating radar[J]. Acta Pedologica Sinica,2015, 5251024–1030. [17] LINDSETH R O. Synthetic sonic logsA process for strati- graphic interpretation[J]. Geophysics,1979,4413–26. [18] COOKE D A,SCHNEIDER W A. Generalized linear inversion of reflection seismic data[J]. Geophysics, 1983, 486 665–676. [19] 周竹生,周熙襄. 宽带约束反演方法[J]. 石油地球物理勘探, 1993,285523–536. ZHOU Zhusheng,ZHOU Xixiang. Band-constrained inver- sion[J]. Oil Geophysical Prospecting,1993,285523–536. [20] FERGUSON R J,MARGRAVE G F. A simple algorithm for band-limited impedance inversion[R]. CREWES Research Re- port,1996,8211–10. [21] SCHMELZBACH C , TRONICKE J , DIETRICH P. High-resolution water content estimation from ground-penetrating radar reflection data by impedance inversion[J]. Water Resources Research,2012,4888505. [22] 李静. 随机等效介质探地雷达探测技术和参数反演[D]. 长 春吉林大学,2014. LI Jing. Ground penetrating radar detection and parameter in stochastic effective medium[D]. ChangchunJilin University, 2014. [23] LI Jing, ZENG Zhaofa, LIU Cai, et al. A study on lunar regolith quantitative random model and lunar penetrating radar parameter inversion[J]. IEEE Geoscience and Remote Sensing Letters, 2017,14111953–1957. [24] 刘钰. 探地雷达数据波阻抗反演方法及其应用研究[D]. 杭 州浙江大学,2018. LIU Yu. The study of ground penetrating radar impedance inver- sion and its application[D]. HangzhouZhejiang Uni- versity,2018. [25] BUSCH S,KRUK J V D,BIKOWSKI J,et al. Quantitative permittivity and conductivity estimation using full wave in- ChaoXing 218 煤田地质与勘探 第 48 卷 version of on ground GPR data[J]. Geophysics,2012,776 79–91. [26] 张彬. 探地雷达中的逆时偏移及速度估计[D]. 长沙中南大 学,2010. ZHANG Bin. Reverse time migration and velocity estimation of ground penetrating radar[D]. ChangshaCentral South Univer- sity,2010. [27] BRADFORD J H. Measuring water content heterogeneity using multifold GPR wit