基于IHT算法的EMT金属探伤稀疏成像方法

孙春光， 何 敏， 曾星星， 冯肖维

(1.上海海事大学物流工程学院，上海 201306; 2.华北电力大学控制与计算机工程学院，北京 102206)

0 引 言

金属探伤技术的目的是探测金属材料表面或内部的缺陷，是工业生产中大型金属结构安全运行的重要保障措施，目前常见的探伤方法有射线探伤、磁粉探伤、超声波探伤、电容层析成像等[1-3]，这些技术各有优缺点。电磁层析成像[4](electromagnetic tomography,EMT)是一种基于电磁感应原理、可将被测对象电磁特性的分布进行图像重建的技术，具有非接触性、低成本、无辐射、可在线检测、图像化等优势，在冶金、轨道交通等方面应用前景广阔。该技术应用于金属结构探伤中，金属缺陷部分主要由空气构成，空气电导率和金属电导率相差较大，因此可以根据电导率分布判断缺陷情况。首先施加一定的激励电磁场，缺陷的存在会导致激励电磁场发生变化，采集变化后的电磁场信息并采用图像重建算法可以得到被测金属结构的缺陷分布图像，通过信息提取进一步可以得到缺陷的长度、宽度等信息。其中，图像重建是EMT探伤成像过程中的重要一环，其精度直接影响整个检测过程的准确度。

裂纹是金属结构缺陷的最常见形式，实际检测中发现，裂纹在金属结构构件上的分布是稀疏的[5]，说明常见的带裂纹金属结构可以被归类为稀疏分布特性对象。目前在层析成像技术研究领域，多数研究从l1、l2范数约束的角度进行求解。张玲玲等[6]研究了基于1范数的电阻层析成像图像重建算法。王丕涛等[7]将l1范数作为图像重建中的数据项和正则化项，对凸优化问题进行求解，同时也对l1范数和l2范数的4种组合进行了仿真分析。Ye等[8]设计了一种扩展灵敏度矩阵，将介电常数分布进行稀疏表示。但在EMT用于金属结构探伤的应用研究中，有关裂纹的稀疏性分析及其图像重建很少。文献[7]虽然对缺陷的稀疏性有一定介绍，但采用l1正则化约束，而研究表明l0正则化的稀疏性高于l1正则化[9]，对稀疏信号的成像效果更好。

针对金属结构裂纹的稀疏分布特性和l0范数更适合求解稀疏信号的特点，本文使用l0范数最小化约束。但直接求解l0范数是NP-hard问题，本文研究了稀疏信号重建理论中的贪婪类算法，该类算法在求解问题时，每一步都选择最好或最优的计算，使最终结果接近最优解。其中迭代硬阈值(iterative hard thresholding,IHT)算法是一种较流行的解决非线性逆问题的方法，该算法通过l0范数将图像重建的不适定逆问题正则化，通过IHT进行求解得到被测量的分布。因此，本文将IHT算法用于EMT金属探伤的图像重建中。

1 EMT金属探伤的图像重建

1.1 EMT金属探伤的理论基础

EMT金属探伤技术的数学模型[10]如下式所示：

式中：σ—被测物体的电导率分布矩阵，S/m；

v———检测线圈的感应电压列向量，V；

F金属缺陷情况与检测线圈上感应电压之间的映射关系。

式(1)可以用佩亚诺型泰勒公式表示，如下式所示：

σ0被测金属没有裂纹时的电导率分布矩阵，S/m；

o( )——表示泰勒展开的余项。

当电导率变化比较小时，泰勒公式的佩亚诺型余项部分可以省略，式(2)可以简化为：

式中：B=v-F(σ0)——被测金属电导率改变导致检测线圈电压变化的矩阵；

G——表征裂纹分布的被测金属体电导率分布矩阵。

1.2 EMT金属探伤中的常见图像重建算法

根据先验信息灵敏度矩阵S和测量到的电压列向量B求解方程（4），得到表征裂纹位置和大小信息的被测金属体电导率分布矩阵，并由图像形式将裂纹信息表现出来，这就是EMT探伤中的逆问题求解，即图像重建。因为测量数据有限，而被测区域内的像素点数远远多于检测的数据量，即方程未知数的个数远远多于方程数，方程（4）的求解具有病态性。求解这种病态方程，常见的方法包括线性反投影（LBP）算法、Landweber迭代算法、Tikhonov正则化算法。

2 基于稀疏分析的迭代硬阈值图像重建算法

2.1 EMT探伤图像重建的稀疏性分析

传统图像重建算法往往根据问题的病态性，单纯的从数学角度分析求解过程，如Landweber迭代算法由最速下降法演变而来，以最小二乘为依据，在数据残差的负梯度方向修正计算结果，进而求得该病态方程的优化解，一些改进的Landweber迭代算法也仅仅是改善解的限制条件[11]；Tikhonov正则化则是将不适定问题转化为近似适定问题，将非凸解集转化为凸优化问题[12]。这些传统算法都是建立在被测对象数理模型的既定基础上，即：模型一旦确定，则不再考虑被测对象本身，图像重建往往转化为数学领域的方程求解问题，因此传统图像重建算法纷纷从方程求解范围、求解方向、求解速度等方面去做文章，以提高重建质量。但正是传统图像重建过程中忽视了被测对象本身的特性，因此反而图像重建质量遭遇了瓶颈。如Tikhonov正则化算法，它基于l2正则化约束条件，所得到的解不是稀疏的，会造成重建的图像在缺陷与非缺陷交界处过度平滑，图像重建质量受限。

本文对待测对象（带裂纹的金属结构）进行深入分析，发现金属结构特别是大型金属结构（比如港机），裂纹出现的位置通常为受力集中部位，且该部位往往出现单条裂纹，如图1所示。因此作者认为，港机金属结构上的裂纹分布并不是密集的，具备稀疏特性。

图1 常见港机裂纹

2.2 基于l0正则化的IHT图像重建算法

信号的稀疏性在信号处理和逼近理论中有悠久的历史，在数据压缩[13]、去噪[14]等方向有重要作用。从数学上来讲，如果离散信号u中最多有K个非零元素，则该离散信号u是K阶稀疏的，即 ‖u‖0≤K，K阶稀疏信号的集合可以用下式表示：

稀疏信号的重建即为基于测量信号y，求解最具稀疏性的目标信号，即得到稀疏分布的缺陷图像，也就是求解下式：

式中：I⊆ {1,···,N}——索引集；

φi——矩阵 Φ 的第i列。

因使用组合的方法直接求最稀疏的解是NP-hard问题，特别是N比较大时更不现实。但是，在经典的稀疏逼近方法中，贪婪算法逐步选择矩阵Φ的列来逼近y，进而确定索引集I，其中迭代硬阈值(iterative hard thresholding,IHT)算法[15-16]常用于解决非线性逆问题。IHT计算结构简单，设初始重建目标为xˆ0，通过下式得到预测值：

其中T(x,K)为保留x中K个较大值而将其他值置零的结果，循环计算中会收敛到最优的固定点ˆ。

目前，求解非线性不适定问题的主要方法有正则化和迭代两种方法，研究发现，正则化是求解该问题的主要方法。正则化过程中常用的有l2范数、l1范数、l0范数，作者分析了l2、l1、l0范数约束下的最小化问题求解，发现方程解的稀疏度按照l2、l1、l0范数约束的顺序增加[17]；虽然l2和l1范数约束下的求解方法比较多，但是解的精度和稀疏性有待提高，而l0范数约束下的方程解的稀疏性较好、计算结果比较准确。

因此，在IHT图像重建算法中采用l0正则化，病态方程（4）的求解问题可以转化为解的优化问题，如下式所示：

优化式(8)的目标函数为：

直接优化式(9)比较难，这里引入替代目标函数，如下式所示：

如果 ‖S‖2＜ 1，则替代目标函数是目标函数的一种优化，而且也是替代函数的极小化，即为Majorization-Minimization优化算法。式(10)可以改写为：

式中：G*=G′+ST(B-SG′)，为计算得到的电导率分布矩阵即，即为式(8)要解决的优化问题，求解得到：

其中H为非线性算子。

也可以表示为：

式(14)中的G′=Gn,表示上一次迭代计算的结果，反复迭代式(15)即可求出像素矩阵G。

3 仿真实验与分析

为了对比使用不同算法时金属缺陷的图像重建效果，文中使用AnsoftMaxwell软件进行了仿真实验，如图2所示为6线圈传感器结构，所有线圈的轴线平行，检测线圈可以得到较大的感应电压，同时线圈轴线垂直于被测金属面，穿过被测金属的磁通量最大，提高了测量精度。根据雷登变换，该传感器结构使用较少的线圈数，可以得到整个检测区域的物场信息。传感器的线圈匝数均为500，正弦激励信号电压的峰峰值为10 Vpp，仿真实验中选择1 010低碳钢，电导率为2 000 000 S/m，相对磁导率为1，实际应用中一般检测深度1.5 mm深度的缺陷，根据趋肤效应公式计算得到激励频率为50 kHz，传感器线圈外圆以内作为有效检测范围。

图2 EMT传感器仿真模型

仿真实验中使用模型扰动法测量得到灵敏度矩阵，然后分别使用Landweber迭代算法、Tikhonov正则化算法、IHT算法进行图像重建，如图3所示（图中右侧0～1的图例表示像素矩阵G的元素），半径R分别为4 mm、3 mm、2 mm的圆形缺陷图像重建效果。从图3中可以看出Landweber迭代算法重建的图像，缺陷边缘部分比较模糊，干扰伪影较多，Tikhonov正则化算法的图像重建效果比Landweber迭代算法稍微改善，而本文介绍的IHT算法图像重建效果比Tikhonov正则化算法的图像清晰度高，缺陷形状的准确度增加。

图3 不同算法对不同尺寸缺陷重建效果

为了定量比较不同图像重建算法下，不同尺寸的金属缺陷图像重建效果，这里使用图像重建相对误差RE作为评价指标，其定义式为：

将图3中设定的电导率分布和不同算法计算得到的电导率分布代入式(16)，即可计算得到使用不同图像重建算法时，不同缺陷的图像重建相对误差，如图4所示。

图4 不同算法对不同尺寸缺陷重建的相对误差

从图4可以看出，Landweber迭代算法重建图像的相对误差最大，Tikhonov正则化算法的图像重建相对误差比Landweber迭代算法小一点，IHT算法的图像重建相对误差最小。相同算法情况下，随着缺陷尺寸的减小，图像相对误差都有一定的增加。

4 实验分析

为了验证仿真实验结果，本文采用了自行开发的EMT金属探伤系统进行了实验。该探伤系统由硬件电路、传感器、上位机三部分组成，其中硬件电路包括控制器、激励源、信号调理电路、通道选择电路，如图5所示。

图5 EMT系统结构图

系统工作原理为：控制器使激励源产生正弦激励信号，经过信号调理电路和通道选择电路后施加到激励线圈上，检测线圈的感应电压经过信号调理电路后送入控制器，处理后传送给上位机进行图像重建。传感器在整个系统中具有重要作用，根据电磁感应原理，线圈轴线平行时，通过线圈的磁通量最大，因此本文使用了如图6所示的6线圈传感器结构，激励线圈和检测线圈复用，圆形虚线内为有效检测区域。

图6 六线圈传感结构

实验中圆形缺陷半径分别是5 mm、3 mm，长条缺陷宽度为2 mm，长度为6 mm，传感器通入的激励信号是50 kHz、10 V的正弦信号，传感器线圈匝数为500，线圈内半径4 mm，外半径6 mm。首先使用扰动法获得灵敏度矩阵，然后使用不同算法进行图像重建，得到的图像如图7所示。

图7 不同算法对不同缺陷的重建效果

根据图7中实验的不同算法对比不同尺寸缺陷重建效果，可以看出Tikhonov正则化算法比Landweber迭代算法图像干扰少一点，但总体效果都没有IHT算法的图像重建效果好，符合仿真结果。使用相对误差作为标准，得到的结果如图8所示，可以看出，实际实验结果与仿真结果相符。

图8 不同算法对不同缺陷重建图像的相对误差

为比较不同算法的收敛速度，对同一缺陷用相同算法计算20次，然后求平均计算时间，得到3种算法的计算时间如图9所示。虽然重建图像时间都大于3 s，但对于实时性要求不高的情况，如港机金属在线监测中，因裂纹变化是准静态情况，故满足实时性要求。对于不知道缺陷位置的情况，可以将传感器安装在可移动平台上，传感器在被测金属表面移动，即可快速检测缺陷。

图9 不同算法重建图像时间

5 结束语

本文分析了金属探伤中缺陷分布的稀疏性，结合l2、l1、l0范数最小化约束解的特点，采用了l0范数约束下的IHT算法进行EMT金属探伤的图像重建。使用AnsysMaxwell软件进行仿真实验，并将IHT算法与Landweber迭代算法、Tikhonov正则化算法进行对比，结果表明l0范数约束下的IHT算法图像重建清晰度高、准确性更好；在根据仿真模型搭建的实验平台上进行探伤实验，实验结果与仿真结果相符，为EMT探伤的实际应用具有现实意义。