基于电磁层析的地面目标内部结构图像重建

秦瑶　王其富　屈乐乐

(1.河南工业大学信息科学与工程学院,郑州 450001;2.河南工业大学 粮食信息处理与

控制教育部重点实验室,郑州 450001;3.河南省科学院应用物理研究所有限公司,郑州 450001;

4.沈阳航空航天大学电子信息工程学院,沈阳 110136)



基于电磁层析的地面目标内部结构图像重建

秦瑶1,2王其富3屈乐乐4

(1.河南工业大学信息科学与工程学院,郑州 450001;2.河南工业大学 粮食信息处理与

控制教育部重点实验室,郑州 450001;3.河南省科学院应用物理研究所有限公司,郑州 450001;

4.沈阳航空航天大学电子信息工程学院,沈阳 110136)

摘要针对电磁层析成像技术在工程物探中，采用钻孔方式进行地下介质内部结构探测有损的不足，为适应电磁层析成像技术无损探测的需求，研究了对地面目标内部结构探测的电磁层析成像技术，基于电磁层析的探测方式、联合成像方法和图像重建技术，为弥补跨孔探测方法应用于地面目标探测时，其边界条件与传统的地下目标探测不同，会对初至波的判断造成影响的不足，采用时域有限差分法仿真分析了波的传播状态，确定了初至波的提取方法，给出了初至零点的判定方法，并采用以上方法进行实验证实了对地上目标进行电磁层析成像的可行性和有效性.

关键词电磁波;层析成像;初至波;图像重建

DOI10.13443/j.cjors.2015021101

Electromagnetic tomography technology to image internal structure of overground target

QIN Yao1,2WANG Qifu3QU Lele4

(1.CollegeofInformationScienceandEngineering,HenanUniversityofTechnology,Zhengzhou450001,China;2.KeyLaboratoryofGrainInformationProcessingandControl,HenanUniversityofTechnology,MinistryofEducation,Zhengzhou450001,China;3.InstituteofAppliedPhysicsHenanAcademyofSciences,Zhengzhou450001,China;4.ShenyangAerospaceUniversityCollegeofElectronicandInformationEngineering,Shenyang110136,China)

Abstract Electromagnetic wave tomography technology is widely used in geophysical prospecting. Drilling method is the common way to detect the inside structures of underground targets. To enlarge electromagnetic wave tomography application range and raise its universality, this paper aimed to do the study of electromagnetic wave tomography on overground targets. At first, discuss the detecting method and unite tomography imaging algorithm. Then, according to the different boundary conditions of overground targets detecting, this paper uses FDTD method to simulate and analyze the distributing stage of electro-magnetic wave, makes clear the extract method of primary wave,and also gives out zero point exact method of primary arrival wave. At last, experiment results verify the effectiveness of electromagnetic wave tomography technology aimed on overground target.

Keywords electromagnetic wave; tomography; primary wave; image reconstruction

引言

电磁层析成像技术是20世纪末逐渐发展和成熟起来的一种基于电磁感应原理的层析成像技术．它通过在待测区域周围布置收发天线获取测量数据, 结合数学物理关系模型,采用图像重建算法反演物体内部电磁参数的空间分布图像,从而实现物体内部结构成像．美国康涅狄格大学、树城州立大学、俄罗斯莫斯科国立大学、国家地震局地球物理研究所、辽宁省地震局、吉林大学、南京大学等单位都曾开展过利用电磁层析成像技术进行矿藏探测、地下水的勘探、溶岩结构分析等工程物探方面的应用研究,并取得了相当丰富的成果[1-5]．上述应用研究都是针对地下目标内部结构的探测,收发天线通过钻孔的方式探入待测体内部获取测量数据,从而得到两钻孔之间区域的电磁特征参数分布图像．目前,各个行业对无损探测技术的需求越来越迫切,电磁波层析成像技术作为一种衰减小、穿透性强、探测精度高、成像效果好的无损探测手段,应该在更多的领域得到应用、发挥作用．针对电磁层析对地面目标内部结构探测成像的可行性研究,首先分析了跨孔探测方法应用于地面目标测量时由于边界条件不同所引发的新问题;接着构建了一个仿真模型,采用时域有限差分(Finite Difference Time Domain,FDTD)方法仿真分析了模型中波的传播状态,从而确定了初至波的提取方法,并给出了初至波零点的判定方法;最后通过实验证实了针对地上目标进行电磁层析成像的可行性和有效性．此项研究成果可以应用在固件的内部结构缺陷探测、桥墩和隧道顶部的坚固度检测,粮仓内部水分分布、异物探测等内部结构具有电磁参数差异的众多领域．

1理论基础

1.1探测方法

电磁层析成像的探测方法如图1所示．在待测体两侧的对应位置分别布置发射天线Tx和接收天线Rx．测量时,发射天线先固定在位置t1,接收天线从位置r1逐次移动至位置rn,每移动一个位置记录一道数据,共获得n道测量数据;再将发射天线固定至位置t2,接收天线仍从位置r1逐次移动至位置rn,又获取n道接收数据;按此规律接收天线逐次移动至tn,最终可获得n×n道接收数据[6]．

图1　层析成像测量方法示意图

1.2射线追踪

射线追踪原理是几何光学近似,是一种有效的波场近似计算方法,也是研究电磁层析成像正演计算的主要方法[7]．传统的射线方法包括初值问题的试射法和边值问题的弯曲法,但这两种方法都存在射线覆盖密度不足、难于处理介质中较强的速度变化、难于求取多值全局最小走时以及计算效率较低等缺陷[8]．目前常用的新型射线追踪算法不再局限于波的射线路径描述,而是直接从费马定理(Fermat)或惠更斯原理(Huygens)出发,采用等价波前来描述波场特征．其中最具代表性的是黃联捷等人提出的波前法[9]．

波前法射线追踪的计算是从源点出发,选取L×L(L小于单方向上的网格总数N)的计算方块,计算出波从源点到计算网格点的透射走时、射线路径和射线长度．然后把除源点之外的所有网格点相继当作次级源,再选取L×L的计算方块,计算出波从次级源点到计算网格点的透射走时、射线路径和射线长度．将每次计算出来的走时加上从主振源点的走时,作为波从主振源点到该网格点的走时,并记录下相应的射线路径位置及射线长度．当所有网格点均已当作次级源,即完成射线追踪的计算,最后从接收点出发选取最小走时,反向确定射线路径．

1.3联合成像算法

层析成像包括两种成像方法:波速层析和衰减层析．波速层析是利用接收天线接收波形的初至波时间进行反演成像;衰减层析是利用接收波形的幅值或频率信息进行反演成像[10]．波速层析可以对待测区域内部介质的波速特性分布成像,衰减层析可以对待测体内部介质的品质因子(或衰减系数)的分布特征成像．将两种算法相结合的联合成像算法能够较完整地反映待测区域内的波速、衰减等物性分布结构[11]．联合层析成像的实现流图如图2所示．首先从雷达接收天线的接收信号中提取初至波走时、幅度、频谱等信息,然后先利用初至波走时数据进行波速层析成像,再利用波速层析中得出的射线路径,结合初至波幅度、频谱信息进行衰减成像[12],最终通过图像重建算法得到待测区域的波速、品质因子分布图像．常见的图像重建算法有反投影法(Back Projection,BP)、代数重建法(Algebraic Reconstruction Technique,ART)、联合迭代重建法(Simultaneous Iterative Reconstruction Technique,SIRT)、奇异值分解正则化(Tikhonov Singular Value Decomposition,TSVD)以及最小二乘QR分解法(Least Squares QR Decomposition Algorithm, LSQR)等[13-15]．波速与衰减联合层析的具体实现方法可参考文献[6],在此不再赘述．

图2　波速与衰减联合层析实现流图

2初至波的提取

初至波的提取是层析成像的一个重要环节,初至波提取的精度在很大程度上影响着层析成像图像重建的精度．对于钻孔探测方式,收发天线通过钻孔探入地下,置于待测体内部,收发天线周围为无限大的同种介质,接收波形初至波走时的提取能够较严格地满足费马原理,因此,实际测量到的最先到达接收位置的波就是所需的初至波．而对于地面目标内部结构的电磁层析成像,由于待测体的高度、宽度范围有限,收发天线在待测体表面的测量过程中,待测体表面的绕射波会对真正的透射初至波的准确提取产生影响[16]．因此,地面待测体电磁层析初至波的提取,必须通过对待测体的仿真分析,观察各个测量位置上的初至波干扰情况,合理分析、剔除干扰,才能获取真正的初至波,从而提高图像重建的精度．

2.1初至波仿真分析

初至波的合理判断和分析是在对待测体进行电磁仿真的基础上实现的．仿真模型针对埋有目标的沙箱进行剖面的图像重建,沙箱和目标体的尺寸如图3(a)所示,电磁层析的测量剖面如图3中虚线所示．图3(b)为测量剖面的俯视图．需要注意的是,为了考察电磁波在测量剖面区域以外部分的传播对初至波提取的影响,仿真时需将图3(b)所示的测量剖面进行扩展,添加四周的空气区域．

(a) 内部结构图

(b) 测量剖面俯视图图3　仿真模型

FDTD法是模拟电磁波传播的常用方法[17]．FDTD法进行电磁仿真的具体实现方法可参考文献[18]．设雷达发射中心频率为3.5 GHz的高斯脉冲,沙子介电常数ε沙=4,电导率γ沙=10-3S/m,异物介电常数ε异物=8,电导率γ异物=10-1S/m．仿真区域为边长d=0.42 m的正方形,网格步长为0.003 m,即将此核心仿真区域剖分为140×140个网格．由于实际仿真中,需在核心区四周添加空气区域,故仿真图中实际显示的网格数为157×157．图4(a)～(c)分别为发射天线位于图4所示测量剖面的t1、t(n+1)/2和tn时,三个不同时间观察到的仿真波形．

(a) 发射天线位于t1位置

(b)发射天线位于t(n+1)/2位置

(c) 发射天线位于tn位置图4　FDTD仿真波形

由图4可以看出,沙箱外壁会对电磁波产生反射,同时部分电磁波会在沙箱表面发生绕射;穿过箱体的透射波,在介质中匀速向前传播,当波前遇到内部的异常体时,因异常体介电常数大于沙子,电磁波在异常体中的波速小于粮食中的波速,波前速度变慢．综上,在图3(b)所示的模型中,对初至波干扰最大的是沿箱体的绕射波．

空气中的绕射波可通过设置合适的时间阈值来消除．在测量过程中,测量系统的位置、移动步长等参数一旦确定,收发天线在各个位置对应的绕射波的到达时间也就可以通过计算得到准确的确定,之后可针对不同位置的接收波形,设置不同的时间阈值,屏蔽掉空气中初至波的绕射波之后将再次出现的次至波作为透射波初至波即可．在实际测试中,次至波初至波的准确判断也并不容易,尤其是当绕射波对波形进行干扰之后,加上接收信号微弱、仪器本身噪声过大等因素,更难于读准．可借鉴机械波中的基线控制线等方法进行判定[19]．

从上述过程中可以看出,地面无损层析初至波的提取,最重要的通过对待测模型的合理分析,消除初至波的干扰波,从而提取真正的初至波．

2.2初至波零点的判定

确定初至波后,需要判定初至波零点的位置．由于层析成像的测量结果中包含大量射线,因此对于初至波位置的选取必须适合自动提取．常用的初至波零点位置的选择如图5所示[20]．其中,A点为真正意义上的初至波,即波到达的最初起始位置,机械波初至波零点位置的判定大都采用此位置．而对于频率较高的电磁波,此位置的拾取存在困难,然而B、C、D位置的数据自动拾取在程序中的判定原则都十分简便,因而是电磁波探测信号采集中常用的判点．

图5　初至波零点判定位置

又由于走时层析需要提取的是真正初至波到达时间,即A点,而波在传播过程中的脉冲展宽无法确切计算．设原始脉冲宽度为τ0,传播过程中的脉冲展宽为τ′,则接收到的实际信号脉冲宽度为τ=τ0+τ′．令A、B、C、D点的走时分别为tA、tB、tC、tD,则真正初至波:

(1)

比较以上三个等式可知,tB所引入的脉冲展宽误差最小,因此选用B位置作为初至波的零点位置．

3实验

实验平台采用矢量网络分析仪和vivaldi天线搭建,实验中的待测体模型仿照图2搭建,沙箱中的异物为砖块．模仿图2,选取平行于地面的水平剖面作为测量剖面．为了减弱地面反射,将沙箱置于桌上,同时用四个支架架起两根贴有刻度的塑杆,将发射天线和接收天线分别固定于两个塑杆上．实验所用vivaldi天线的发射带宽为1～6 GHz,中心频率fc=3.5 GHz．定义垂直测线方向为宽度方向,宽度d=0.42 m;沿测线方向为长度方向,天线紧贴沙箱沿标注的测线位置逐点移动,步长0.03 m,测量点数15个,故测量区域总长度l=d=0.42 m．

实验中,固定好天线后,先将收发天线直接相对,即移开待测沙箱,此时测得的直达波走时为2.207ns,由此计算得发射延迟时间τ为: τ=t测-t计算=2.207ns-0.42m/c空气=0.807ns．之后,将收发天线对准沙箱,且均挪至同一长度位置,记录走时值,再将收发天线按刻度同方向移动一个步长,再记录走时值,共测15组,求得走时的均值为T=3.427ns,故可估算出沙子的介电常数为:

首先采用波前法对图3(b)所示的模型进行射线追踪,提取理论走时数据;再用2.1节方法对沙箱各个方向上的绕射波走时进行分析计算,排除干扰后进行测量数据走时的提取．理论走时和实验走时的对比结果如图6所示,其相对误差为5.27%．

图6　计算走时与实测走时对比图

图7(a)、7(b)为采用LSQR算法对走时和频谱质心偏移数据进行图像重建的结果．其中,调整因子取值为0.01,两种算法的迭代次数均为20次,成像结果中还做了中值滤波及聚类分析处理．图7中异物的位置信息基本准确,但由于测量误差系统噪声等因素的影响,非目标区域并不均匀,通过走时层析和衰减层析结果的对比,能够进一步确认异常区域的位置;同时从图7能够判断异常区域是低速区,同时品质因子Q小于沙,从而可以判断异常体目标为介电常数实部、虚部均大于沙的非金属介质．

(a) 走时层析

(b) 频谱质心偏移层析图7　图像重建结果

4结论

通过仿真和实验证实了电磁层析技术对地面目标内部结构探测成像的可行性,同时还得到以下几点经验:1) 层析成像能够分辨的异常区域的长宽尺寸,最小不能小于一个波长;2) 天线移动步长不能太小,至少应大于λ/2,否则,无法有效区别相邻位置间初至波走时的不同,成像结果与实际相差较大;3) 网格尺寸的划分应与天线移动步长相当(同一个量级),网格尺寸过大,图像重建结果粗糙,容易丢失目标;网格尺寸过小,图像重建质量并不会有明显提高,反而会大大增加测量时间．

参考文献

[1] LIU L B, GUO T S. Seismic non-destructive testing on a reinforced concrete bridge column using tomographic imaging techniques[J]. Journal of geophysics and engineering, 2005,2(1): 23-31.

[2] JOHNSON T C, ROUTH P S, KNOLL M D. Fresnel volume georadar attenuation-difference tomography [J]. Geophysical journal international, 2005, 162(1): 9-24.

[3] SUDAKOVA M S, KALASHNIKOV A Y, TERENTIEVA E B. GPR tomography feasibility study in the process of air-filled voids search in engineering structures[C]//Near Surface Geoscience 2014-20th European Meeting of Environmental and Engineering Geophysics. Athens, September 14-18, 2014.

[4]张旻舳, 师学明. 电磁波层析成像技术进展[J].工程地球物理学报,2009,6(4):418-425.

ZHANG M Z, SHI X M. Review of electromagnetic tomography technique[J]. Chinese journal of engineering geophysics,2009,6(4):418-425. (in Chinese)

[5]李才明, 张善法. 电磁波层析成像阻尼因子引入与应用[J].地球物理学进展,2005,20 (1):221-224.

LI C M, ZHANG S F.The application and estimate method of damping gene in the electromagnetic wave computerized tomography [J]. Progress in geophysics, 2005,20(1):221-224. (in Chinese)

[6] 秦瑶, 陈洁, 方广有, 等.雷达层析成像技术对粮仓储量信息的探测[J].电波科学学报, 2010, 25(1): 66-71.

QIN Y, CHEN J, FANG G Y, et al. Radar tomography technology to detect reserves of barn[J]. Chinese journal of radio science, 2010, 25(1): 66-71. (in Chinese)

[7] 孙方, 王红光, 康士峰, 等. 大气波导环境下的射线追踪算法[J].电波科学学报, 2008, 23(1): 179-183.

SUN F, WANG H G, KANG S F. A ray tracing algorithm for duct environment[J]. Chinese journal of radio science, 2008, 23(1): 179-183. (in Chinese)

[8] 和锐, 杨建思, 张翼. 地震层析成像方法综述[J]. CT理论与应用研究, 2007, 16(1): 35-48.

HE R, YANG J S, ZHANG Y. A review on the technology of seismic tomograph[J]. CT theory and applications, 2007, 16(1): 35-48.(in Chinese)

[9] 黄联捷, 李幼铭, 吴如山. 用于图像重建的波前法射线追踪[J].地球物理学报, 1992, 35: 223-233.

HUANG L J, LI Y M, WU R S. Wave-front ray tracing method for image rebuild [J]. Chinese journal of geophysics, 1992, 35: 223-233. (in Chinese)

[10]王辉, 常旭, 刘尹克. 时间域相邻道地震波衰减成像研究[J]. 地球物理学报, 2001, 44(3): 396-403.

WANG H, CHANG X, LIU Y K. Seismic neighboring traces attenuation tomography in time domain[J]. Chinese journal of geophysics, 2001, 44(3): 396-403. (in Chinese)

[11]ZHOU C, LIU L, LANE J W, et al, Nonlinear inversion of borehole-radar tomography data to reconstruct velocity and attenuation distribution in earth materials [J]. Journal of applied geophysics, 2001, 47: 271-284.

[12]LIU L B, LANE J W, QUAN Y L. Radar attenuation tomography using the centroid frequency downshift method[J].Journal of applied geophysics, 1998, 40(1/3):105-116.

[13]赵连锋. 井间地震波速与衰减联合层析成像方法研究[D]. 成都: 成都理工大学,2003.

ZHAO L F. Study on crosswell seismic tomography combining velocity and attenuation[D]. Chengdu: Chengdu: University of Technology, 2003. (in Chinese)

[14]欧明, 甄卫民, 於晓,等.一种基于截断奇异值分解正则化的电离层层析成像算法[J]. 电波科学学报,2014,29(2):345-352.

OU M, ZHEN W M, YU X, et al. A computerized ionospheric tomography algorithm based on TSVD regularization[J]. Chinese journal of radio science,2014,29(2):345-352. (in Chinese)

[15]刘劲松, 刘福田, 刘俊,等.地震层析成像LSQR算法的并行化[J].地球物理学报,2006,49(2):540-545.

LIU J S, LIU F T, LIU J, et al. Parallel LSQR algorithms used in seismic tomography[J].Chinese journal of geophysics, 2006, 49(2):540-545. (in Chinese)

[16]FU L, LIU S X, LIU L B. Internal structure characterization of living tree trunk cross-section using GPR: numerical examples and field data analysis[C]//15th International Conference on Ground Penetrating Radar-GPR 2014, December 1, 2014: 155-160.

[17]方广有, 张忠治, 汪文秉. 无载频脉冲探地雷达性能分析及其FD-TD法数值模拟计算[J].电波科学学报,1997,12(1):44-50.

FANG G Y, ZHANG Z Z, WANG W B. Study about the properties of a carrier-free pulse ground penetrating radar and their numerical simulation by 2.5D FD-TD method[J].Chinese journal of radio science, 1997,12(1):44-50. (in Chinese)

[18]秦瑶,陈洁,方广有.粮仓储量信息无损探测技术研究[J]. 中国粮油学报, 2010, 25(4):51-55.

QIN Y, CHEN J, FANG G Y. Non-destructive detection of barn reserves information [J]. Journal of the chinese cereals and oils association,2010, 25(4):51-55.(in Chinese)

[19]吴新璇.混凝土无损检测技术手册[M]. 北京: 人民交通出版社, 2003: 172-177.

WU X X. Concrete NDT technical manual[M]. Beijing: People Transportation Press, 2003, 172-177. (in Chinese)

[20]YELF R, YELF D. Where is true time zero[J]. Electromagnetic phenomena, 2006, 7(18): 158-163.

秦瑶(1981-),女,河南人,博士,河南工业大学信息学院副教授,主要从事电磁波探测及成像技术研究．

王其富(1981-),男,广东人,博士,河南省科学院应用物理研究所副研究员,主要从事信号处理、智能用电及物联网相关技术的研究工作．

屈乐乐(1983-),男,河南人,博士,沈阳航空航天大学电子信息工程学院副教授,主要从事雷达信号处理及成像技术研究．

IEEE国际地球科学与遥感大会IGARSS

2016年将在北京召开

IEEE国际地球科学与遥感大会(IGARSS，International Geoscience And Remote Sensing Symposium)是地球科学与遥感(GRS)领域最具影响力的国际学术会议，现与会人数达2700人，我国科学工作者出席IGARSS的人数已达数百人。该会议每年轮流在北美、欧洲和亚太地区举办，至今已有30年历史。但是，从未在我国举办过。

2012年7月21日，由中国科学院国家空间科学中心主任吴季率领的中国申办“IGARSS 2016”团队，在德国慕尼黑举行的IEEE GRSS理事会上做了申办IGARSS2016的申请报告，有力地阐述了中国在空间遥感科学领域所取得的进展以及作为IGARSS2016主办城市-北京的申办条件。最后，经过IEEE GRSS理事会成员激烈的讨论，中国获多数票表决通过。

IGARSS2016大会主席由中国科学院国家空间科学中心吴季和复旦大学金亚秋共同担任，会议将得到国家各主管部门、研究所、高等院校的大力支持。

期望我国从事遥感科学领域研究与应用的专家学者和研究生，对IGARSS 2016筹备工作给予高度关注与大力支持。IGARSS2016在北京举办将促进中国在该重要的科学技术领域有更深入广泛的国际合作和交流，增强中国遥感科学的国际知名度，提升我国科学研究的水平。

咨询地址： IGARSS2016@nssc.ac.cn

作者简介

中图分类号TN959

文献标志码A

文章编号1005-0388(2016)01-0143-07

收稿日期：2015-02-11

秦瑶, 王其富, 屈乐乐. 基于电磁层析的地面目标内部结构图像重建[J]. 电波科学学报,2016,31(1):143-149. DOI: 10.13443/j.cjors.2015021101

QIN Y, WANG Q F, QU L L. Electromagnetic tomography technology to image internal structure of overground target[J]. Chinese journal of radio science,2016,31(1):143-149. (in Chinese). DOI: 10.13443/j.cjors.2015021101

资助项目: 国家自然科学基金(No.61201389;No.61302172;No.61201390;No.11201120)； 河南省高等学校重点科研项目(No.16A510003); 粮食信息与控制部级重点实验室开放基金课题(No.KFJJ-2015-102); 粮食公益性行业科研专项(No.201413003)； 河南省科技厅重点科技攻关项目(No.152102310083)； 河南工业大学省属高校基本科研业务费专项资金(No.2014YWQQ11)

联系人： 秦瑶 E-mail:eqinyao@163.com