2017, Vol. 60
2. 近地面探测技术重点实验室, 江苏无锡 214035
2. Science and Technology on Near-Surface Detection Laboratory, Jiangsu Wuxi 214035, China
近地表地球物理(Near-surface Geophysics)应用范围涵盖自然灾害监测、工程质量探测、环境污染监测、地下管网与未爆弹药(Unexplored ordnance UXO)探查、考古等众多领域(Butler，2005；Butler et al.，2006，2012；程红杰等，2006；Miller，2008；Beran et al.，2013)，其中针对未爆弹药、金属管道等有限导体的探测与公共利益、公共安全甚至军事需求密切相关(Butler et al.，2006，2012；Beran et al.，2013).
不同地球物理探测方法对近地表有限导体的探测与识别能力各不相同：磁法探测未爆弹药具有成本低，探测深度大的优点，但是探测结果易受环境噪声干扰，且磁异常响应对目标体形状区分能力弱，因此探测虚警率高达90%(Billings，2004；Beran et al.，2015)；探地雷达对目标体形状、结构区分能力强，但是探测效果受地表杂波、低阻覆盖等干扰严重(Chen et al.，2001；Carin et al.，2002)；瞬变电磁法被广泛用于地质勘探以及工程、灾害探测等领域(薛国强等，2008；丁艳飞等，2012；陈卫营和薛国强，2013；张爽等，2014；齐彦福等，2015)，由于其抗干扰能力强、对目标体形状有较强的分辨力，近年来被广泛用于近地表有限导体的探测并已经取得丰硕成果.
Grimm(2001) 通过仿真对比磁性与非磁性目标体的瞬变电磁响应，重点研究炮弹中涡旋电流的分布规律，这种侧重于机理的定性研究有助于更好的认识有限导体的瞬变电磁响应规律，却难以对各种有限导体的响应进行定量描述，也难以实现各种有限导体的区分、识别.Khadr等(1998) 研究有限导体的磁极化张量与其长径比的关系，指出未爆弹药与金属弹片长径比差异大，可用于目标体识别，而论文并未对目标体瞬变电磁响应随时间衰减特性进行充分利用，这些衰减特性可以更加准确的实现有限导体的识别.Fernández等(2008) 将用于频率域有限导体数据处理的标准化激励方法(SEA)用于时间域电磁响应分析，以实现对结构复杂、材质各异的目标体瞬变电磁响应的建模与分析.Barrow和Nelson(2001) 、Pasion等(2007) 、Billings等(2010) 以及Song等(2012) 建立有限导体响应的三维感应偶极子模型，正演轴对称目标体瞬变电磁响应，根据实测响应估计目标体的位置、轴向等信息，实现对目标体的形状、结构的识别.
本文在三维感应偶极子模型基础上定义了目标体的特征响应，并将有限导体特征响应进行参数化建模，以进一步研究有限导体瞬变电磁响应与形状、结构的关系.为此，本文首先通过理论计算与实验测量得到球体特征响应，分析其响应规律从而建立有限导体特征响应的参数化模型；其次选择几种典型轴对称有限导体，以研究目标体的外径、长度、壁厚、体积等对响应的影响规律；最后分析非对称有限导体的特征响应.在此基础上揭示轴对称目标体与普通目标体瞬变电磁响应的特性，为进一步建立精确的目标体响应模型以及目标体识别提供理论指导与实验依据.
2 有限导体响应的三维感应偶极子模型有限导体的三维感应偶极子模型是将有限导体等效为一组中心重合、方向正交的三维感应偶极子，以三维感应偶极子响应表示有限导体的瞬变电磁响应.
如图 1所示，在一次场Bp的激励下，目标体的总感应偶极子m可以分解为3个独立的感应偶极子m1、m2、m3，总感应偶极子m的三分量表示为
|
(1) |
|
图 1 有限导体的三维感应偶极子模型 Fig. 1 Three-dimensional induction dipole model of a finite conductor |
根据一次场BP计算总感应偶极子m为
|
(2) |
其中E为感应偶极子m1、m2、m3的方向矩阵，公式为
|
(3) |
其中L为感应偶极子m1、m2、m3随时间衰减的磁矩，公式为
|
(4) |
若空间中点P到原点O的距离为r，方向为er，感应偶极子m在P点产生二次场Bs为
|
(5) |
将二次场Bs对时间t求导，得到单位面积接收线圈的感应电压为
|
(6) |
由式(2) 和(4) 可得：
|
(7) |
由式(6) 和(7) 可见：目标体响应不仅取决于一次场、接收线圈与目标体的相对位置、距离等，还取决于目标体感应偶极子的时域衰减特征，即目标体的特征响应.根据(6) 式，本文定义偶极子m1的特征响应l1(t)：偶极子m1在单位磁场激励下，在BP、e1、er同向时，单位面积接收线圈在单位距离上的感应电压为
|
(8) |
由(8) 式可得，三维感应偶极子构成的有限导体的特征响应l(t)为
|
(9) |
(9) 式表明，目标体特征响应l(t)由三维偶极子特征响应l1(t)、l2(t)、l3(t)构成，由于各个偶极子特征响应由目标体形状、结构决定，因此目标体的特征响应也直接反映其形状、结构信息，相反，根据目标体形状、结构信息可以获得目标体特征响应.本文研究有限导体特征响应与其形状、结构等之间的关系，需要获得有限导体特征响应的定量描述，所以首先借助球体的瞬变电磁响应完成有限导体特征响应的定量表述，同时研究球体的特征响应与其形状的关系.
3 球体的特征响应一般有限导体的响应不能给出明确的数学表达，只有球体响应可以给出解析结果.为此本文首先分析球体理论与实测响应，并以此为依据建立一般目标体特征响应的参数化模型，同时研究球体特征响应与其半径的关系.
3.1 球体的特征响应如图 2所示：自由空间中球体半径为a，磁导率为μ，相对磁导率为μr，电导率为σ，一次场BP沿z方向.
|
图 2 磁性球体响应 Fig. 2 Response of a magnetic sphere |
一次场BP阶跃关断，从0时刻开始，在球体正上方r点处，沿z方向的二次场BS可以表示为
|
(10) |
LB(t)为二次场的衰减函数，公式为
|
(11) |
其中τ=σμa2，s为响应模式，从1开始直到无穷，qs为超越方程的解，公式为
|
(12) |
(10) 式经求导可得：
|
(13) |
对比(7) 、(9) 与(13) 式可以看出，球形目标体的特征响应为
|
(14) |
由(14) 式可见，球体特征响应幅度与其半径成正比，电导率成反比，衰减速度则与半径的平方、电导率、磁导率均成反比.理论上从球体的特征响应中直接可以获取球体半径、磁导率、电导率等信息.由于球体完全对称，特征响应具有各向同性，其三维表达式为
|
(15) |
根据(14) 式，计算在μr=40，σ=3.73×106s的条件下，半径32 mm铁球的特征响应，并与实测结果进行对比.
由图 3可见，铁球体的理论响应在早期按照幂函数进行衰减，衰减速率从早期t-1/2增加到中期t-3/2，在晚期则以指数函数进行衰减；实际测量结果受到系统影响，不能准确的反映目标体早期t-1/2衰减特性，但在中期、晚期与理论响应高度吻合.因此，本文只针对目标体中、晚期响应分段建立其特征响应的参数化模型.
|
图 3 铁球体理论响应与实征响应 Fig. 3 Theoretical and measured responses of an iron sphere |
根据上面的分析，针对有限导体中、晚期响应分段建立特征响应模型为
|
(16) |
其中kl、ke分别表示中期幂函数衰减、晚期指数衰减的幅度；α为幂函数衰减的指数，λ为指数衰减的时间常数.
根据公式(16) ，利用最小二乘拟合算法，对图 3中实测信号进行参数拟合，算法为
|
(17) |
其中，l1(t(i)，kl，α，ke，λ)为参数化的响应模型在t(i)时刻的值，dat(i)为t(i)时刻的测量值，Φ(kl，α，ke，λ)为最小二乘目标函数，当Φ取最小值时的参数即为拟合参数.
依据(17) 式对实测的球体响应进行最小二乘拟合，结果如表 1所示.从表 1可知，由于拟合数据包含部分早期t-1/2衰减与中期t-3/2衰减，因此参数α的拟合结果略小于1.5；参数λ为球体晚期指数衰减时间常数，拟合结果为11 ms.
|
表 1 球体参数拟合结果 Table 1 Fitting results of sphere parameters |
根据拟合结果构建球体参数化特征响应，并与实测数据比较(图 4).由图 4可见，在40 μs到1 ms、20 ms到70 ms两段区间内，实测信号与拟合曲线吻合度较高，误差在5%以内.在40 μs之前、70 ms之后的误差超过10%.早期响应受到发射电流关断延时的影响，实测值高于理论值导致误差偏大；晚期响应受到噪声干扰严重导致误差较大.总体上，利用分段参数化响应模型可以准确描述目标体的特征响应.
|
图 4 参数拟合结果 (a)球体的实测响应与拟合曲线；(b)拟合的相对误差. Fig. 4 Results of parameter fitting (a)Measured response and fitted curves of sphere；(b)Relative errors of fitting results. |
本文将以此模型为基础研究各种典型有限导体的特征响应，分析各种目标体形状、结构对其特征响应的影响规律.
4 轴对称目标体的二维特征响应复杂目标体特征响应一般不能通过理论计算获得，本文通过实验方法，利用上文中提出特征响应的参数化模型，研究复杂目标体的特征响应与其形状的关系.为了研究方便本文按照结构的对称性将有限导体分为两类：一类为具有轴对称结构的有限导体，如金属管道、未爆炸物等；另一类则不具有对称结构.轴对称有限导体作为一类典型的探测目标，其特征响应可分为两类：沿主轴方向(l1(t))和垂直于主轴方向(l2(t))，公式(9) 中的l3(t)与l2(t)相等.
本文选择12种典型目标体，这12种目标体的参数与编号如表 2所示，其中9种为轴对称结构目标体，3种为复杂目标体，如图 5所示.
|
图 5 12种目标体照片 Fig. 5 Pictures of 12 targets |
|
|
表 2 目标体编号 Table 2 Number of targets |
实测上述9种目标体特征响应，考虑到早期信号畸变，晚期受到噪声影响，本文只取50 μs到20 ms之间的信号进行抽道，提取目标体特征响应参数进行分析.管状目标体参数包括长度、外径、壁厚、体积等，本文分别研究这几种参数对响应的影响.
4.1 长度对特征响应影响实测表 2中编号为n1—n4铸铁管特征响应，获得如图 6所示的目标体特征响应.由图 6可见，同一目标体沿轴向特征响应l1(t)的幅度比垂直轴向特征响应l2(t)的幅度高出一个量级，且衰减更慢，这一规律可用于区分沿轴向特征响应与垂直轴向特征响应；随着目标体长度减小，特征响应l1(t)、 l2(t)的响应幅度均不断减小，可见较短的目标体在实际探测中更难被发现；图中300 mm与200 mm铁管特征响应l1(t)在晚期重合，可见为管状目标体沿轴向的特征响应在目标体长度较大时，晚期衰减特性相同，可以根据这一特性区分不同类型的管道.
|
图 6 不同长度铸铁管的特征响应 (a)沿主轴方向特征响应；(b)垂直主轴方向特征响应. Fig. 6 Characteristic responses of cast iron pipes with different lengths (a)Characteristic responses along major axis direction；(b)Characteristic responses vertical to major axis direction. |
根据公式(17) ，利用最小二乘拟合算法提取这些特征响应参数，得到如表 3所示结果.表 3中的拟合结果显示，由于铁管管壁较薄，晚期衰减的时间常数一般小于300 μs，且l1(t)是l2(t)的1.5~2倍，这些特征可用于初步判断目标体的类型；随着管状目标体长度减小，特征响应l1(t)早期幂函数的指数α逐步变大，衰减加快；l2(t)早期衰减体现的目标体半径，由于目标体半径不变，l2(t)早期幂函数的指数α维持不变，因此衰减速率也不变；目标体l1(t)、 l2(t)的早期响应幅度(kl)比值随着长径比减小而减小，因此根据比值可以大致确定目标体的长径比，而长径比对于区分目标体至关重要.
|
|
表 3 n1—n4号目标体模型参数 Table 3 Model parameters of targets No.n1—n4 |
实测表 2中编号为n1、n6铸铁管特征响应，获得如图 7所示的目标体特征响应.由图 7可见：在目标体长度、壁厚相同的条件下，目标体特征响应早期幅度相近，外径越大，目标体晚期衰减时间常数越大，衰减速率越低，目标体垂直轴向特征响应l2(t)的早期衰减速率越低.所以管状目标体外径影响目标体晚期衰减速率以及l2(t)的早期衰减速率，对目标体早期响应幅度影响有限.
|
图 7 不同直径目标体特征响应 Fig. 7 Characteristic responses of targets with different diameters |
根据公式(17) ，利用最小二乘拟合算法提取这些特征响应参数，得到如表 4所示结果.表 4中拟合结果对目标体特征响应的早期响应幅度与晚期衰减速率均有较好的反应；在晚期，90 mm外径目标体指数衰减的时间常数为75 mm外径的10倍，因此晚期信号衰减更慢，信号幅度更大，即管状目标体直径主要影响其晚期衰减的时间常数，决定晚期信号幅度.
|
|
表 4 目标体模型参数 Table 4 Model parameters of targets |
实测表 2中编号n2、n5铸铁管特征响应，获得如图 8所示的目标体特征响应.由图 8可见，目标体长度、外径相同的情况下，目标体早期响应的幅度与衰减速率完全相同，目标体壁厚越大，感应电流形成环路的等效电阻越小，因此导体晚期响应的时间常数越大，响应衰减越缓慢.
|
图 8 不同壁厚目标体特征响应 Fig. 8 Characteristic responses of targets with different thicknesses |
根据公式(17) ，利用最小二乘拟合算法提取这些特征响应参数，得到如表 5所示结果.表 5的拟合结果显示，随着目标体壁厚由3 mm增加到5 mm，目标体沿轴向的特征响应l1(t)与垂直轴向的特征响应l2(t)的晚期衰减的时间常数λ，分别增加了38%与36%，考虑测量误差，可以认为壁厚对l1(t)与l2(t)的晚期时间常数影响一致；l1(t)与l2(t)的早期响应幅度与衰减速率基本一致，即壁厚的增加仅影响目标体晚期指数衰减的时间常数.
|
|
表 5 目标体模型参数 Table 5 Model parameters of targets |
由此可见：目标体特征响应幅度以及衰减速率主要取决于目标体长度，垂直轴向特征响应的早期衰减速率主要取决于目标体的直径，晚期指数衰减时间常数主要取决于目标体直径与壁厚；目标体外形相同则早期响应相同，目标体壁厚越大晚期指数衰减时间常数越大.
4.4 体积对特征响应影响实测表 2中编号n7、n8、n9的特征响应，这3个目标体形状相似、大小不同，实测结果如图 9所示.通过图 9可知，在早期，目标体特征响应l1(t)的幅度比l2(t)的幅度高出约一个量级，且衰减更为缓慢，120 mm迫击炮弹的特征响应幅度最大，晚期衰减最缓慢；82 mm、60 mm迫击炮弹体积相近，响应幅度以及晚期衰减速率等均较为接近.
|
图 9 不同型号迫击炮弹的特征响应 (a)沿主轴方向特征响应;(b)垂直主轴方向特征响应. Fig. 9 Characteristic responses of different shells (a)Characteristic responses parallel to major axis direction;(b)Characteristic responses vertical to major axis direction. |
根据公式(17) ，利用最小二乘拟合算法提取这些特征响应参数，得到如表 6所示结果.表 6中，迫击炮弹的指数衰减时间常数λ在1~10 ms之间，远大于管状目标体，且目标体体积越大，响应幅度与时间常数越大；120 mm与82 mm迫击炮弹幂函数的指数α基本相同，即早期衰减速率相同，60 mm迫击炮弹幂指数α最大，衰减最快；晚期响应中，82 mm与60 mm迫击炮弹的时间常数接近，远小于120 mm迫击炮弹，因此晚期信号衰减快，信噪比低.总体看，迫击炮弹弹体厚度决定了其晚期衰减的时间常数远大于铸铁管，因此衰减速度慢，在整个测量区间中信号幅度大，与铸铁管有较明显的差异.
|
|
表 6 目标体模型参数 Table 6 Model parameters of targets |
大部分地下有限导体不具有轴对称结构，这些有限导体对轴对称目标体探测会产生强烈的干扰，导致系统探测虚警率居高不下.本文选取50 mm铸铁管的T型头(编号n10) ，沿轴向开槽的铸铁管(编号n11) 以及金属碎片(编号n12) 作为典型的非对称结构目标体，实测这些目标体响应结果如图 10所示.
|
图 10 复杂目标体的特征响应 (a)T型接头特征响应；(b)开缝铁管特征响应；(c)金属弹片特征响应. Fig. 10 Characteristic responses of complex targets (a)Characteristic responses of T-joint；(b)Characteristic responses of slotted iron pipe；(c)Characteristic responses of metal shrapnel. |
由图 10可见，复杂目标体的三维特征响应l1(t)、l2(t)、l3(t)互不相同，衰减过程甚至相互交叉；T型头三维特征响应在早期幅度相近，但由于目标体三维结构差异大，导致晚期特征响应差异显著；沿轴向开槽的铸铁管在垂直轴向方向仍具有较高的对称性，因此l2(t)、l3(t)的早期响应完全一致，晚期则由于缝隙导致衰减速率出现差异；金属碎片的三维特征响应之间无任何关联，特征响应甚至相互交叉.
根据公式(17) ，利用最小二乘拟合算法提取这些特征响应参数，得到如表 7所示结果.
|
|
表 7 复杂目标体特征响应模型参数 Table 7 Model parameters of complex targets |
表 7的拟合结果定量反映了复杂目标体三维特征响应之间的差异以及差异的来源：T型头的三维特征响应l1(t)、l2(t)、l3(t)幂函数衰减指数α分别为1.06、1.40、1.46，而晚期指数衰减的时间常数λ分别为0.316 ms、0.310 ms、0.287 ms，可见块状目标体特征响应的分离主要源于早期幂函数衰减过程，晚期指数函数的时间常数则比较接近；轴向对称性遭到破坏的目标体，其垂直轴向响应在早期响应幅度与衰减速率都非常接近，仅晚期指数衰减的时间常数(1.49 ms、1.61 ms)有轻微差异，这种目标体和完好的轴对称目标体响应目标体响应极为相似，实际探测中极难与完整的轴对称目标体进行区分；金属碎片的特征响应四个拟合参数没有特定的规律，三维特征响应早期幂函数衰减速率与晚期指数函数衰减速率不一致，导致特征响应出现交叉现象.
6 结论在近地表有限导体探测中，目标体的形状、结构等信息对于目标体的区分、识别至关重要.本文以三维感应偶极子模型为基础，根据球体的理论响应构建有限导体特征响应的参数化模型，分析有限导体的形状、结构对其瞬变电磁特征响应影响规律.
目标体形状的对称性决定了目标体特征响应的相似性：结构越对称，特征响应的相似程度越高.结构高度对称的球体所有特征响应完全相同；轴对称的有限导体具有沿轴向与垂直轴向的两组特征响应；普通有限导体的特征响应互不相同.由于探测目标多为轴对称有限导体，因此对称性是判断目标体有效性的重要依据.轴对称目标体沿轴向特征响应的响应幅度比垂直与轴向特征响应幅度高一个量级；早期幂函数衰减响应的幅度以及衰减速率主要取决于目标体长度，长度越长，响应幅度越大，衰减越慢；晚期指数衰减响应的时间常数主要取决于目标体的外径与壁厚，外径越大、壁厚越厚，时间常数越大，衰减越慢；目标体沿轴向特征响应幅度与垂直轴向响应幅度之比与目标体的长径比成正比.复杂目标体三维特征响应之间差异显著.
本文的研究表明有限导体的瞬变电磁响应能够准确其反映形状、结构信息，还能以参数化的响应模型精确描述，从而在目标体形状、结构与模型参数之间建立起对应关系.随着这些关系的不断丰富与完善，系统对有限导体的识别能力将得到极大的提高，虚警率将显著降低.
| Barrow B, Nelson H H. 2001. Model-based characterization of electromagnetic induction signatures obtained with the MTADS electromagnetic array. IEEE Transactions on Geoscience and Remote Sensing, 39(6): 1279-1285. DOI:10.1109/36.927450 | |
| Beran L, Zelt B, Billings S. 2015. Detecting and classifying UXO. Journal of ERW and Mine Action, 17(1): Article 15. | |
| Beran L, Zelt B, Pasion L, et al. 2013. Practical strategies for classification of unexploded ordnance. Geophysics, 78(1): E41-E46. DOI:10.1190/geo2012-0236.1 | |
| Billings S D. 2004. Discrimination and classification of buried unexploded ordnance using magnetometry. IEEE Transactions on Geoscience and Remote Sensing, 42(6): 1241-1251. DOI:10.1109/TGRS.2004.826803 | |
| Billings S D, Pasion L R, Beran L, et al. 2010. Unexploded ordnance discrimination using magnetic and electromagnetic sensors:Case study from a former military site. Geophysics, 75(3): B103-B114. DOI:10.1190/1.3377009 | |
| Butler D K. Near-Surface Geophysics.Tulsa: Society of Exploration Geophysicists, 2005. | |
| Butler D K, Bennett H H, Ballard J H. 2006. Overview of multimethod geophysical system development for enhanced near-surface target detection, discrimination, and characterization. The Leading Edge, 25(3): 352-356. DOI:10.1190/1.2184105 | |
| Butler D K, Simms J E, Furey J S, et al. 2012. Review of magnetic modeling for UXO and applications to small items and close distances. Journal of Environmental and Engineering Geophysics, 17(2): 53-73. DOI:10.2113/JEEG17.2.53 | |
| Carin L, Sichina J, Harvey J F. 2002. Microwave underground propagation and detection. IEEE Transactions on Microwave Theory and Techniques, 50(3): 945-952. DOI:10.1109/22.989977 | |
| Chen C C, Higgins M B, O'Neill K, et al. 2001. Ultrawide-bandwidth fully-polarimetric ground penetrating radar classification of subsurface unexploded ordnance. IEEE Transactions on Geoscience and Remote Sensing, 39(6): 1221-1230. DOI:10.1109/36.927444 | |
| Chen W Y, Xue G Q. 2013. Application on coal-mine voids detection with multi-device TEM technology. Progress in Geophys. (in Chinese), 28(5): 2709-2717. DOI:10.6038/pg20130554 | |
| Cheng H J, Hu X Y, Zhang R, et al. 2006. The development of the near surface geophysics. Chinese Journal of Engineering Geophysics (in Chinese), 2(1): 73-76. | |
| Ding Y F, Bai D H, Xu C. 2012. A rapid algorithm for calculating time domain transient electromagnetic responses of a large fixed rectangular loop on the half space. Chinese J. Geophys. (in Chinese), 55(6): 2087-2096. DOI:10.6038/j.issn.0001-5733.2012.06.030 | |
| Fernández J P, Barrowes B, O'Neill K, et al. 2008. A data-derived time-domain SEA for UXO identification using the MPV sensor.//Proc. SPIE 6953, Detection and Sensing of Mines, Explosive Objects, and Obscured Targets XIII. Orlando, FL:SPIE. | |
| Grimm R E. 2001. Time domain measurements and modeling of UXO.//UXO Forum. | |
| Khadr N, Barrow B J, Bell T H. 1998. Target shape classification using electromagnetic induction sensor data.//Proceedings of UXO Forum'98. | |
| Miller R. 2008. Introduction to this special section:Near-surface geophysics. The Leading Edge, 27(11): 1423-1423. DOI:10.1190/1.3036100 | |
| Pasion L R, Billings S D, Oldenburg D W, et al. 2007. Application of a library based method to time domain electromagnetic data for the identification of unexploded ordnance. Journal of Applied Geophysics, 61(3-4): 279-291. DOI:10.1016/j.jappgeo.2006.05.006 | |
| Qi Y F, Yin C C, Wang R, et al. 2015. Multi-transient EM full-time forward modeling and inversion of m-sequences. Chinese J. Geophys. (in Chinese), 58(7): 2566-2577. DOI:10.6038/cjg20150731 | |
| Song L P, Oldenburg D W, Pasion L R. 2012. Estimating source locations of unexploded ordnance using the multiple signal classification algorithm. Geophysics, 77(4): WB127-WB135. DOI:10.1190/geo2011-0379.1 | |
| Xue G Q, Li X, Di Q Y. 2008. Research progress in TEM forward modeling and inversion calculation. Progress in Geophysics (in Chinese), 23(4): 1165-1172. | |
| Zhang S, Liu Z X, Chen S D. 2014. Calibration study on the damping characteristics of TEM sensor. Chinese J. Geophys. (in Chinese), 57(2): 662-670. DOI:10.6038/cjg20140229 | |
| 陈卫营, 薛国强. 2013. 瞬变电磁法多装置探测技术在煤矿采空区调查中的应用. 地球物理学进展(5): 2709–2717. DOI:10.6038/pg20130554 | |
| 程红杰, 胡祥云, 张荣, 等. 2006. 近地表地球物理发展状况综述. 工程地球物理学报, 2(1): 73–76. | |
| 丁艳飞, 白登海, 许诚. 2012. 均匀半空间表面大定源瞬变电磁响应的快速算法. 地球物理学报, 55(6): 2087–2096. DOI:10.6038/j.issn.0001-5733.2012.06.030 | |
| 齐彦福, 殷长春, 王若, 等. 2015. 多通道瞬变电磁m序列全时正演模拟与反演. 地球物理学报, 58(7): 2566–2577. DOI:10.6038/cjg20150731 | |
| 薛国强, 李貅, 底青云. 2008. 瞬变电磁法正反演问题研究进展. 地球物理学进展, 23(4): 1165–1172. | |
| 张爽, 刘紫秀, 陈曙东. 2014. 瞬变电磁传感器阻尼特性的标定研究. 地球物理学报, 57(2): 662–670. DOI:10.6038/cjg20140229 | |