CN110176046B - 收缩系数改进Tikhonov正则化参数的电阻层析成像方法 - Google Patents

Abstract

本发明涉及一种收缩系数改进Tikhonov正则化参数的电阻层析成像方法，步骤如下：对被测场进行测量，获取重建图像所需的相对边界测量向量b和反映其与内部电特性关系的灵敏度矩阵A；求解Tikhonov正则化标量参数λ；根据电场衰减分布特性，预先设置正则化参数的收缩系数β向量；收缩系数改进正则化参数后的Tikhonov正则化算法的求解；对解x进行非负约束和归一化处理；电阻层析成像。

Description

收缩系数改进Tikhonov正则化参数的电阻层析成像方法

技术领域

本发明属于电学层析成像技术领域，涉及一种用收缩系数改进Tikhonov正则化参数的电阻层析成像重建算法。

背景技术

电学层析成像技术(Electrical Tomography，ET)是自上世纪80年代后期逐渐发展起来的过程层析成像技术，基于不同的电特性(电导率/介电系数/复导纳/磁导率)敏感机理，通过敏感场内电特性分布推断场内煤质分布情况。电学层析成像技术主要包括电阻层析成像(Electrical Resistance Tomography，ERT)、电容层析成像(ElectricalCapacitance Tomography，ECT)、电阻抗层析成像(Electrical Impedance Tomography，EIT)和电场层析成像(Electrical Magnetic Tomography，EMT)。电学层析成像在多相流和生物医学领域有广泛的应用前景，可以实现长期、持续监测。

电阻层析成像逆问题是非线性的，通过线性化处理，转化成求解线性方程的问题，但是该线性逆问题具有严重的不适定性。对此，通常选择正则化方法进行求解，其中Tikhonov正则化是应用比较广泛的一种，它基于最小二乘算法，通过先验信息对解的2-范数施加约束，以求得逼近真实解的稳定解，例如M·Vauhkoen等人在1998年发表于《IEEE医学成像》(Medical Imaging，IEEE Transactions on Medical Imaging)第17卷，第2期，第285-293页，题为《基于电阻抗层析成像的Tikhonov正则化和先验信息》(Tikhonovregularization and prior information in electrical impedance tomography)的文章。正则化参数平衡残差范数和解范数的大小，如果参数过大，则正则化解过光滑，残差过大；如果参数过小，正则化影响很小，正则化解将不稳定。所以选择合适的正则化参数具有重要意义，L-curve方法是选择正则化参数常用的方法，例如Per Christian Hansen等人在1993年发表于美国SIAM《科学计算杂志》(SIAM Journal on Scientific Computing)第14卷，第6期，第1487-1503页，题为《L-curve在离散不适定问题的正则化中的应用》(The useof the L-curve in the regularization of discrete ill-posed problems)文章及D·Calvetti等人在2000年发表于《计算与应用数学杂志》(Journal of Computational andApplied Mathematics)第123卷，第1期，第423-446页，题为《Tikhonov正则化和L-curve解决大规模离散不适定问题》(Tikhonov regularization and the L-curve for largediscrete ill-posed problems)的文章。

上述研究中选择的正则化参数是标量参数，没有考虑物体分布场域的空间特性，所以根据物体空间分布特性，有人提出了自适应正则化参数选择的方法，例如D.Watzenig等人在2004年发表于IEEE Transactions on Magnetics第40卷，第2期，第1116-1119页，题为《重构问题的自适应正则化参数调整》(Adaptive regularization parameteradjustment for reconstruction problems)的文章及X.Z.Song等人在2016年发表于《流量测量与仪表设备》(Flow Measurement and Instrumentation)第50卷，第1-12页，题为《一种用于电阻层析成像的自适应Tikhonov正则化参数选择方法》(An adaptive Tikhonovregularization parameter choice method for electrical resistance tomography)的文章。这些文章中都通过迭代算法进行自适应参数选择。

在电阻层析成像的重建场内，从激励电极发出的激励电流产生的电场是衰减分布的，所以距离激励电极越远的边界数据对电导率变化的敏感性最低；在圆形阵列场中，成像平面的中心部分具有最低的灵敏度，从而使得选定标量正则化参数后，中心内含物的成像效果比其他位置差。例如David S.Holder在2004年出版于CRC press，题为《电阻抗层析成像：方法、历史与应用》(Electrical Impedance Tomography：Methods,History andApplications)，第188页内容。

发明内容

本发明的目的是利用电场衰减分布的特性，基于L-curve选择参数，提出一种收缩系数对Tikhonov正则化参数作改进，并结合电导率分布的非负特性，提高电阻层析成像逆问题中小圆形内含物的成像质量，尤其是中心物体成像。

本发明的技术方案如下：

一种收缩系数改进Tikhonov正则化参数的电阻层析成像方法，适用于电阻层析成像算法重建图像。该方法将电阻层析成像逆问题近似看成线性逆问题Ax＝b，其中A是灵敏度矩阵，b是相对边界测量值向量，x是与成像灰度值对应的待重建的内部电导率分布向量，被称为解向量。基于Tikhonov正则化算法对图像进行重建，步骤如下：

(1)对被测场进行测量，获取重建图像所需的相对边界测量向量b和反映其与内部电特性关系的灵敏度矩阵A；

(2)利用L-curve法求解Tikhonov正则化标量参数λ；

(3)Tikhonov正则化算法的目标函数：

其中，||·||是欧几里得范数；

(4)Tikhonov正则化算法的求解公式：

x＝(A^TA+λI)^-1A^Tb

(5)根据电场衰减分布特性，预先设置正则化参数的收缩系数β向量，步骤如下：

将被测的圆形的场域，划分成16层圆域，其圆心也是所画圆的圆心，最小圆半径为1/16，以1/16等量地递增半径画圆，最外环内像素点对应的正则化参数的收缩系数设为1，依次往内一环，对应的正则化参数的收缩系数按照指数衰减，中心圆内对应点的参数收缩系数最小，设为0.1，具体的选取方式为：1和0.1之间等距选取14个数，和1,0.1一起从大到小排列，用ξ代表各个数，通过u＝10^-ξ求解，得到对应的16个收缩系数；这16个收缩系数，依次分配到所划分的16层圆域，然后得到与各像素点对应的正则化参数的收缩系数向量β；

(6)收缩系数改进正则化参数后的Tikhonov正则化算法的求解公式记为：

x＝(A^TA+diag(λβ))^-1A^Tb

其中，向量β与解x的维数相同，diag(·)是将向量λβ构成对角阵的函数，向量各个元素构成对角阵的对角元素，对角阵的非对角元素均为零。

(7)对解x进行非负约束和归一化处理，公式如下：

非负约束：

归一化处理：

其中x_k表示解向量的第k个元素，max(·)和min(·)分别表最大值和最小值；

(8)根据所求解的灰度值进行成像。

本发明根据电场衰减特性对重建算法中的正则化参数进行改进，提出的预先设置的参数收缩系数的使用，使小圆形物体重建图像的质量得到改善，具有更高的重建精度，同时本发明降低了算法对噪声的敏感性。

附图说明

图1为本发明的一种收缩系数改进Tikhonov正则化参数的电阻层析成像算法的流程图；

图2为物体摆放位置及变化方向示意图；

图3(a)为边界电压数据随物体位置的变化情况图，图3(b)为归一化边界数据后的半对数坐标图；

图4为本发明提出的收缩系数分布场域的分层示意图；

图5为本发明中选取的部分仿真对比图像：5％噪声下得到的各模型重建图像。

具体实施方式

下面结合附图对本发明中改进Tikhonov正则化参数的电阻层析成像重建算法的具体实施方式进行详细说明。

本发明是一种用收缩系数改进Tikhonov正则化参数的电阻层析成像重建算法，基于L-curve选参，用根据电场衰减特性预先设置的收缩系数对参数进行优化，得到空间分布的正则化参数，从而使小圆形内含物的重建图像更接近真实分布。

如图1所示，为本发明的一种收缩系数改进Tikhonov正则化参数的电阻层析成像算法的流程图。如图2所示，为小圆形内含物在场域中的各个位置及变化方向，通过改变物体位置，得到边界测量数据平均值随所处位置的变化情况，如图3(a)所示；图3(b)所示为将边界电压平均值数据归一化后的半对数坐标图。从而反映出电场衰减特性对重建数据的影响及衰减的规律。为了量化衰减规律，将之近似看作指数衰减。

如图4所示，为正则化参数的收缩系数分布场域的分层示意图，最外环内像素点对应的正则化参数的收缩系数设为1，依次往内一环，对应的正则化参数的收缩系数按照指数衰减，中心圆内对应的参数收缩系数最小，设为0.1，具体的选取方式为：1和0.1之间等距选取14个数，和1,0.1一起从大到小排列，如果用ξ代表各个数，则通过u＝10^-ξ求解，得到对应的16个收缩系数；这16个收缩系数，依次分配到图4所示划分的16层区域，然后得到与各像素点对应的正则化参数的收缩系数向量β。

为了更好地体现本发明中收缩系数的使用对重建图像质量的影响，分别选择了五个模型，包含不同数量的小圆形内含物，利用Tikhonov正则化算法及改进后的算法进行模型重建，得到如图5所示的仿真重建结果的对比图，第一行是真实分布图，第二行是Tikhonov正则化算法重建图像，第三行是改进后的算法的重建图像，实施例子包括以下具体步骤：

一种收缩系数改进Tikhonov正则化参数的重建算法，适用于电阻层析成像，该方法将电学成像问题近似为一个线性不适定问题Ax＝b，其中A是灵敏度矩阵，b是相对边界测量值向量，x是与成像灰度值对应的待重建的内部电导率分布向量，被称为解向量。

Tikhonov正则化的目标函数为：

其中λ＞0是正则化参数，由L-curve方法求得相对适宜的参数值，||·||是欧几里得范数；

Tikhonov正则化算法的求解公式为：

x＝(A^TA+λI)^-1A^Tb

收缩系数改进正则化参数后的求解公式记为：

x＝(A^TA+diag(λβ))^-1A^Tb

其中向量β是根据衰减特性预先设置的收缩系数，与解x的维数相同，diag(·)是将向量λβ构成对角阵的函数，向量各个元素构成对角阵的对角元素，对角阵的非对角元素均为零。

重建模型的具体步骤如下：

(1)针对五个含有不同数量小圆形内含物的模型，分别获取各自的边界测量值和灵敏度矩阵：

被测场采用16电极模型，边界电极的分布在图2中体现出来；采用相邻激励相邻测量的方式，循环激励循环测量，采集边界电压数据，每个模型获得208个测量值，分别构成向量；逆问题求解中的b为模型测量数据向量值与为内含物时测得的数据向量值之差；因为实际中噪声不可避免，所以加入5％的随机噪声。

灵敏度矩阵是根据被测场没有内含物时的边界测量数据，结合灵敏度理论计算得到的，公式记为：

其中A_ij是灵敏度系数矩阵A的(i,j)位置元素，

分别为第i个电极和第j个电极对在激励电流为I_i、I_j时的电势分布；

(2)使用Tikhonov正则化算法，利用L-curve法求出五个模型对应的重建的正则化参数；

(3)使用Tikhonov正则化算法求各个模型的解，再进行非负约束和归一化处理，然后进行成像：

成像结果如图5第二行所示。

(4)对正则化参数使用收缩系数，使用改进后的求解公式：

x＝(A^TA+diag(λβ))^-1A^Tb

(5)求得各个模型的内部相对电导率值，对各个解向量进行非负约束，然后归一化处理；

(6)根据处理后的解得到的灰度值，进行成像:

求解得到的各个模型的成像结果如图5中第三行所示。

从用改进前和改进后的算法对五个模型进行成像，发现改进后的成像质量更好，改进的算法对于噪声的灵敏度变低了。

以上所述实例选择了部分成像结果，本发明不限于该实例和附图所述内容。凡基于本发明所公开的基本思想和方法下完成的等效和修改，均在本发明保护的范围。

Claims (2)

1.一种收缩系数改进Tikhonov正则化参数的电阻层析成像方法，适用于电阻层析成像算法重建图像；该方法将电阻层析成像逆问题近似看成线性逆问题Ax＝b，其中A是灵敏度矩阵，b是相对边界测量值向量，x是与成像灰度值对应的待重建的内部电导率分布向量，被称为解向量；基于Tikhonov正则化算法对图像进行重建；步骤如下：

(1)对被测场进行测量，获取重建图像所需的相对边界测量向量b和反映其与内部电特性关系的灵敏度矩阵A；

(2)求解Tikhonov正则化标量参数λ；

(3)Tikhonov正则化算法的目标函数：

其中，||·||是欧几里得范数；

(4)Tikhonov正则化算法的求解公式：

x＝(A^TA+λI)^-1A^Tb

(5)根据电场衰减分布特性，预先设置正则化参数的收缩系数β向量，步骤如下：

将被测的圆形的场域，划分成16层圆域，其圆心也是所画圆的圆心，最小圆半径为1/16，以1/16等量地递增半径画圆，最外环内像素点对应的正则化参数的收缩系数设为1，依次往内一环，对应的正则化参数的收缩系数按照指数衰减，中心圆内对应点的参数收缩系数最小，设为0.1，具体的选取方式为：1和0.1之间等距选取14个数，和1,0.1一起从大到小排列，用ξ代表各个数，通过u＝10^-ξ求解，得到对应的16个收缩系数；这16个收缩系数，依次分配到所划分的16层圆域，然后得到与各像素点对应的正则化参数的收缩系数向量β；

(6)收缩系数改进正则化参数后的Tikhonov正则化算法的求解公式记为：

x＝(A^TA+diag(λβ))^-1A^Tb

其中，向量β与解x的维数相同，diag(·)是将向量λβ构成对角阵的函数，向量各个元素构成对角阵的对角元素，对角阵的非对角元素均为零；

(7)对解x进行非负约束和归一化处理，公式如下：

非负约束：

归一化处理：

其中x_k表示解向量的第k个元素，max(·)和min(·)分别表最大值和最小值；

(8)电阻层析成像。

2.根据权利要求1所述的方法，其特征在于，步骤(2)利用L-curve法求解Tikhonov正则化标量参数λ。

Images