基于电磁计算的非平整表层回波抑制方法

全息穿透成像雷达(Holographic Subsurface Imaging Radar，HSIR)利用电磁波的穿透特性和微波全息成像工作原理，通过二维合成孔径扫描获取回波数据，利用相位信息进行相干成像，以二维图像呈现介质内埋藏目标或者材料缺陷的投影图像，具有高空间分辨率、高穿透能力、快速、无损等优点，可广泛应用于无损检测、建筑、反恐、安检等领域[1-3]。

全息穿透成像雷达向介质发射电磁波时，除目标外，不同来源的反射波也被接收，如天线间的直接耦合波、空气-介质交界面的表层回波、介质内腔隙颗粒等杂散物体造成的散射波[4]等。对于非平整表面介质，尤其对于高程差可与发射波长相比拟的介质而言，表层回波随介质表面高程变化而变化，没有一致的规律，难以校准和去除，且介质表面相较于介质内埋藏目标距离天线更近，电磁波传输过程中损耗更小，因此一般来说非平整表层回波远大于埋藏目标回波，这使得目标回波淹没于表层回波等杂波之中[5-6]。

为抑制非平整表面效应对穿透雷达成像产生的影响，现有研究采取了多种物理方式。典型的有在估算介质介电常数的基础上，获得与此介电常数相似的柔性填充物，以填平非平整面[7-8]，然而此种方法难以获取合适的填充物；另一种方法是设计共形天线阵列或天线的扫描面与介质面共形，这对于形状规则的介质可以做到，比如圆柱、球形等，这方面的应用可以追溯到对乳腺癌的检测—半球形天线阵列对乳房探测成像[9]，以及圆柱扫描系统安检设备对人体衣物下隐藏武器的探测等[10-12]，这类方法只能应用于特定场景。基于介质的表层高程及其属性以估计表层回波并加以消除是一类更通用的方法：文献[13]提出对于特定形状的介质，利用其形成的回波特征来识别并消除表层回波，实测实验验证该方法对特定介质有效；文献[14]提出在获得非平整表面高程数据后，基于几何光学法对表层回波拟合，通过背景对消抑制实测数据中的表层回波。文献[14]的方法能够抑制表层回波，然而因其拟合精度不高的问题，杂波抑制过后，仍有过量杂波影响成像质量。

本文提出了一种基于电磁计算的非平整表层回波抑制方法，该方法依据雷达散射特性加性模型，在已知非平整表面高程信息条件下，利用电磁计算工具拟合表层回波，并与实测数据相减消除表层回波，得到近似目标回波，从而达到抑制表层回波的目的。

1 非平整表层回波抑制方法原理

全息穿透成像雷达扫描示意图如图1所示，对该扫描空间建立三维直角坐标系，其中雷达扫描面平行于XY平面，收发天线以固定扫描间隔依次沿X向Y向二维扫描，扫描面与介质之间为空气层。电磁波离开发射天线后，在介质表面处发生后向散射与折射。假设介质的介电参数处处均匀，则电磁波在介质内传输时，除遇目标外，不发生反射。此外，由于收发天线无法完全隔离，部分电磁波经空间直达接收天线，形成直接耦合波。因此，雷达在扫描位置(x,y)处接收到的回波能量分别是收发天线之间的直接耦合波So(x,y)、介质表面反射回波Ss(x,y)以及目标回波St(x,y)三者的矢量叠加[4]：

S(x,y)= Ss(x,y)+St(x,y)+So(x,y)

对于直耦波，因收发天线相对空间位置保持固定，其值不随扫描位置变化，可作为直流分量去除，不影响对目标的成像。此时回波可简化为

S′(x,y)= Ss(x,y)+St(x,y)

式中，介质表面反射波是关于空间扫描位置的变量，当介质表层不平整时，表层回波随扫描位置变化，使得表层回波与目标回波叠加后，难以直接分离。考虑到介质表层靠近天线，表层回波幅度往往强于目标回波，使得目标回波杂糅于杂波之中，甚至淹没于杂波之中，因此必须将表层回波抑制才能有效成像。

考虑到实际情况下无法测得不包含目标的介质表层雷达回波，可通过电磁数值计算拟合非平整介质表层回波，设拟合回波为Cs(x,y)，基于式(2)，用拟合的表层回波对消扫描数据中的表层回波，即可得近似目标回波

再对S′t(x,y)进行成像处理，即可实现对表层不平整介质内隐藏目标的成像。

在拟合回波Cs(x,y)时，还需附加幅度和相位校准因子，以补偿计算模型与实际雷达系统的误差。总结而言，可得算法流程图如图2所示。

2 表层杂波电磁计算方法

用电磁数值计算方法对表层回波进行拟合时，应根据发射频率、天线特性、介质大小等选择不同数值计算方法，例如时域法、频域法以及高频法。同时，电磁计算工具可对天线精确建模，相比于几何光学法[14]提高了仿真精度。不仅如此，可对回波的最低接收能量进行设置，从而提高拟合精度。

利用电磁计算工具拟合雷达回波的基本步骤为：

1) 结构建模，主要包括：

① 对已知高程数据的非平整表面介质三维建模，如图3(a)所示；

② 量测实际天线尺寸，对收发天线三维建模并设置端口，如图3(b)所示。

2) 设置激励信号模拟雷达发射信号。

3) 计算并获取回波数据。

获取回波数据前，雷达系统的参数：天线类型、天线三维数据、激励源、天线距介质表面的高度、扫描点的位置等需针对具体的雷达系统进行设置。同时，介质的表面高程及介电属性也需与实际尽可能一致。

3 表层回波自适应对消处理

电磁数值计算与实际探测场景设置不可避免地存在一定误差，如介质介电参数、介质表面高程、背景电磁环境，以及雷达系统的建模误差。因此，拟合回波与实测回波在幅度与相位上存在差异。造成差异的原因包括但不限于以上因素难以精确补偿，但对于给定的探测场景，表现为固定差异的部分，可通过相位与幅度校准去除。

3.1 幅度校准

理想条件下，拟合回波与实测无目标背景回波对应扫描点的幅度之比应相同。为应对上述因素造成的干扰，且为了避免随机干扰带来的误差，对实测回波处理，将疑似目标区域去除得到背景回波，对其幅度累加求得与拟合回波相差倍数，用该倍数对拟合回波校准，可用以下公式表达：

式中S′p(x,y)为实测回波中除去疑似目标区域的部分，相应的Cs(x,y)中该部分区域也应去除得到Csp(x,y)。当目标回波完全淹没于杂波中时，可认为S′(x,y)≈S′p(x,y)，Csp(x,y)≈Cs(x,y)。

3.2 相位校准

幅度校准后，采用以下准则估计相位误差，计算准则[14]为

式中C′sp(x,y)为C′s(x,y)中去除疑似目标区域的背景回波，F为代价函数。当代价函数F在0≤t≤1范围内值最小时，此时的t=t0值即为所求，即对拟合回波的相位偏移2πt0值，得到C″s(x,y)= C′s(x,y)ej2πt，可对拟合回波实现校准。

4 实验验证

4.1 场景设置

实验的介质材料为图4(a)、(b)所示的两块非平整表面介质，其中方块尖顶介质为40 cm×40 cm×3.5 cm ABS(Acrylonitrile Butadiene Styrene)塑料(3.5 cm为最大厚度)，最大高程差1.5 cm，双余弦介质为40 cm×40 cm×3.3 cm ABS塑料(3.3 cm为最大厚度)，最大高程差1.3 cm，分别模拟非平滑突变曲面和平滑渐变曲面。如图4(d)所示实测场景，介质在XY扫描平台控制下沿固定距离移动，收发天线连接矢量网络分析仪(图4(c))并固定在介质上方，发射波为7 GHz正弦连续波，扫描方式为停-走-停，扫描范围18 cm×18 cm，扫描间隔6 mm。采用有限积分法(Finite Integration Technique，简称FIT)数值计算方法计算拟合回波，该方法将待求空间网格化，每个网格具有不同的电磁属性，将麦克斯韦积分方程离散化，可采用计算机进行求解。利用电磁计算工具拟合回波的计算模型如图4(e)、(f)所示。

实验中以垂直交叉布置的铜条为目标，如图5所示，铜条紧贴于介质下表面，长10 cm，宽1 cm。为检查天线极化对成像的影响，以及非平整表面对目标回波传播路径影响下的算法性能，图5(a)所示铜条与水平方向成45°角，图5(b)所示铜条置于水平方向。实验中，方块尖顶介质与双余弦介质均使用这两组目标。

4.2 幅相校准及数据处理结果

图6(a)、(b)为两介质电磁数值计算的拟合回波，图6(c)、(d)为无目标的实测背景回波。对比可知，两者所反映的介质表层对电磁波的散射趋势相似。

结合公式(4)、(5)、(6)，幅相校准步骤表示如下：

1) 对实测回波成像，以确定疑似目标区域，若无法确定，即目标回波淹没于杂波中时，跳至步骤3)；

2) 去除实测回波和拟合回波中疑似目标区域；

3) 求得实测回波、拟合回波各点的幅度并累加求和；

4) 求得拟合回波、实测回波幅度累加和之比k；

5) 原始拟合回波各点数据除以k，得到C′s(x,y)。

幅度校准后，进行相位校准：

1) 去除实测回波、拟合回波C′s(x,y)中疑似目标区域；

2) 对拟合回波各点数据附加一固定相位，该相位在[0,2π]范围内等间隔获取，即二维数据增加了一维Z向，得到一三维拟合回波；

3) 实测回波与每一附加固定相位的拟合回波相减，得到一新的三维数据；

4) 遍历该三维数据中的Z向，对每一Z向二维数据幅度累加求和，得到一维向量F；

5) 求向量F中的最小值；

6) 当F取最小时，获取此时拟合回波附加的相位φ；

7) 对拟合回波C′s(x,y)附加相位φ。

图7所示为方块尖顶介质与45°角铜箔的处理结果。其中，图7(a)为原始回波，图7(b)为原始回波成像，图7(c)为原始回波经所提算法抑制杂波的成像，图7(d)为原始回波减背景回波的成像，其中背景回波为实验室条件下得到的数据，在现实应用中该值无法获得，然而此时杂波抑制最为彻底，成像效果最好，本文用做图7(c)成像的对照组。

由图7可知，原始回波成像(图7(b))的表层回波相较于图7(c)、(d)较为严重，杂波形成伪像，且目标像不完整，无法辨别目标的形状大小。抑制杂波后的成像(图7(c))，杂波被有效抑制，其强度远小于目标回波，且目标像得以完整呈现。图7(c)与图7(d)对照，图7(c)的杂波稍强于图7(d)。对图7(c)目标长度测算为100 mm(见图中白色矩形框)，与目标实物的大小、形状一致。

图8除具有图7的实验结论外，可以看到目标像水平方向与垂直方向铜箔的长度不一致，这与实物不同，分析知与天线的极化方式有关，本实验使用角锥喇叭天线，为线极化，会致使某一扫描方向的接收回波相对弱小，使得目标像形状与实际形状存在差异。

采用双余弦介质进行实验，如图9、图10所示。由图9可知，原始回波成像图中，目标本身回波较强，因而可直接成像，但受介质表层杂波的影响，直接成像结果中目标严重失真。抑制杂波后(图9(c))，目标像的形状大小均得以呈现。图10结果与图9结果相似，不再赘述。

若无幅相校准，实测回波与拟合回波直接对消，四组实验的成像如图11所示。分别与图7～10(c)幅相校准后成像图对照，可见图11中杂波较强，目标像被掩埋。证明了幅相校准有效地抑制了拟合回波与实测回波幅相的差异性。

4.3 处理结果评估

为了评估图7～10中各成像图(原始回波成像，抑制杂波后成像，减背景回波成像)的成像质量，采用信杂比评估成像质量提高的程度。为提高数据间比较的有效性，计算信杂比时，同一原始回波的三幅成像图中，选择同一区域作为目标区。计算结果如表1所示。

观察表1中数据可知，信杂比反映的图像质量改善情况与视觉观察一致。应用算法抑制杂波后的成像相较直接成像的信杂比提高明显，尤其方块尖顶介质表层回波较强的情况，所提算法也能很好将其抑制。两介质抑制杂波后的成像图与减背景回波成像图的信杂比相比仍有一定差距，但抑制杂波后的成像已能够对目标检测、判别，对此类存在弱杂波的成像，可通过非线性增益控制算法[15]抑制，以获得高质量成像。

5 结束语

非平整表面介质产生的表层回波是制约全息穿透成像雷达发展的重要因素。本文为抑制介质表面高程差可与波长相比拟的表层回波，提出了一种基于电磁计算的表层回波拟合方法，该方法利用电磁计算工具对已知介质表面高程的表层回波进行拟合，并对拟合数据的幅度、相位自适应校准，之后与实测数据对消，从而抑制其影响。实验结果表明，拟合的表层回波与实测中的表层回波匹配良好，实现了对表层回波的有效拟合，进而抑制了表层回波对成像的影响，改善了成像质量。

[1] IVASHOV S, ZHURAVLEV A, CHIZH M, et al. High Resolution MW Holographic System for NDT of Dielectric Materials and Details[C]∥Proceedings of 16th International Conference on Ground Penetrating Radar, Hong Kong, China:IEEE ,2016:1-5.

[2] 毛强,晋良念,刘庆华.穿墙雷达多维参数人体姿态识别方法[J].雷达科学与技术,2021,19(1):40-47.

[3] CHIZH M, ZHURAVLEV A, RAZEVIG V, et al. Defects Investigation in Thermal Insulation Coatings with Microwave Imaging Based on a 22 GHz Holographic Radar[J]. NDT & E International, 2020, 109:102191.1-102191.8.

[4] 陈诚，刘涛，刘雨，等．强杂波干扰下全息穿透雷达目标分辨方法[J]．现代雷达，2021，43(9)：31-36.

[5] HUANG C L, LIU T. The Impact of an Uneven Medium Surface in Holographic Penetrating Imaging and a Method to Eliminate the Interference[C]∥Proceedings of 2013 IEEE International Conference of IEEE Region 10, Xi’an, China:IEEE, 2014:1-4.

[6] CHEN C, HE Z H, SONG X J, et al.A Subspace Projection Approach for Clutter Mitigation in Holographic Subsurface Imaging[J]. IEEE Geoscience and Remote Sensing Letters, 2022, 19:1-5.

[7] IVASHOV S I, RAZEVIG V V, VASILIEV I A, et al. Holographic Subsurface Radar of RASCAN Type: Development and Applications[J]. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing. 2011, 4(4):763-778.

[8] IVASHOV S I, BUGAEV A S, ZHURAVLEV A V, et al. Holographic Subsurface Radar Technique for Nondestructive Testing of Dielectric Structures[J]. Technical Physics, 2018, 63(2):260-267.

[9] HAGNESS S C, TAFLOVE A, BRIDGES J E. Three-Dimensional FDTD Analysis of a Pulsed Microwave Confocal System for Breast Cancer Detection: Design of an Antenna-Array Element[J]. IEEE Trans on Antennas and Propagation, 1999, 47(5):783-791.

[10] MCMAKIN D L, SHEEN D M, HALL T E, et al. Cylindrical Holographic Radar Camera[J]. Proceedings of SPIE - the International Society for Optical Engineering, 1998, 3575:79-88.

[11] MCMAKIN D L, SHEEN D M, GRIFFIN J W, et al. Extremely High-Frequency Holographic Radar Imaging of Personnel and Mail[J]. Proceedings of SPIE- the International Society for Optical Engineering,2006, 6201:580-591.

[12] SHEEN D M, MCMAKIN D L, HALL T E. Cylindrical Millimeter-Wave Imaging Technique and Applications[C]∥Proceedings of Passive Millimeter-Wave Imaging Technology IX, Kissimmee, FL:SPIE ,2006:1-10.

[13] LI M T, HUANG C L, SU Y. A Method of Removing Interference Fringes on Spherical Subsurface Imaging with Continuous Wave Penetrating Radar[C]∥ Proceedings of 16th International Conference on Ground Penetrating Radar, Hong Kong, China: IEEE:1-4.

[14] RAZEVIG V V, ZHURAVLEV A V, BUGAEV A S, et al. Imaging Under Irregular Surface Using Microwave Holography[J]. 2017 Progress in Electromagnetics Research Symposium-Fall (PIERS-FALL), 2017:172-177.

[15] HUANG Y W, LIU T, HE Z H, et al. Improvement of Surface Penetrating Radar Imaging by Suppressing Clutter Using Nonlinear Gain Control[J].Proceedings of SPIE, 2020, 11519:289-293.