一种基于模型细化的改进牛顿-拉夫逊ERT图像重建法

摘要

本发明公开了一种基于模型细化的改进牛顿-拉夫逊ERT图像重建法，针对电阻层析成像技术中灵敏度分布的不均匀性，通过采用在每个三角形有限元的形心位置增加一个节点的方法，细化有效提高正问题计算精度的有限元模型b，并在算法重建过程中，遵循“计算正问题时采用细化前有限元模型b，修正电阻率分布时采用细化后的有限元模型c及其对应的灵敏度矩阵”的原则。本发明充分利用了优化后有限元模型b与细化后的有限元模型c的优点，在有效提高了正问题计算精度的基础上，不仅提高了灵敏度分布的均匀性，同时改善了Hessian矩阵的病态性，在不影响算法实时性的前提下，有效提高了图像重建质量。

说明书

技术领域

本发明涉及电阻层析成像图像重建，具体是一种基于模型细化的改进牛顿-拉夫逊ERT 图像重建法，属于电学层析成像技术领域。

背景技术

电学层析成像技术是过程层析成像技术的重要分支，包括电阻层析成像(Electrical  Resistance Tomography,ERT)、电容层析成像(Electrical Capacitance Tomography,ECT)、 电阻抗层析成像(Electrical Impedance Tomography,EIT)以及电磁层析成像 (Electromagnetic Tomography,EMT)，该技术具有响应速度快、非侵入、价格低廉等优点， 在医学及工业测量等领域具有广阔的应用前景。

在电学层析成像技术中，有限元模型的拓扑结构对正问题的计算精度与图像重建质量 具有重要的影响。图2中的有限元模型a为传统按等间隔原理剖分的有限元模型，目前已 经通过仿真与实际实验验证了其不合理性；图3中的有限元模型b是以模型a中除最外层 之外的每一层半径为变量，同时将模型中心区域半径变量所在区间剖分成多个子区间，以 敏感场均匀分布时模型均方根值的倒数为适应度函数，并引入三角形最长边与最短边的比 值作为惩罚函数，利用基于区间算法与粒子群算法的改进遗传算法离线优化模型a，得到的 优化后的有元限模型b，该优化后的有限元模型b拓扑结构能有效提高正问题的计算精度与 图像重建质量。但不同拓扑结构的有限元模型不仅影响正问题的计算精度，其对应的敏感 场均匀分布时灵敏度矩阵的条件数与灵敏度分布的均匀性也不同。

此外，电学层析成像技术本身具有的病态性严重影响了图像重建质量，而其病态的根 源是灵敏度矩阵的条件数很大，目前常采用截断奇异值(Truncated Singular Value  Decomposition,TVSD)算法与Tikhonov正则化算法改善其病态性，这两种算法均存在一 定的缺陷，具体如下：

截断奇异值算法通过截断较小的奇异值实现解的稳定性，但解的分辨率损失较严重；

Tikhonov正则化算法通过在目标函数中加入一个罚函数来实现对解的阻尼作用，达到 使解稳定的目的，同时又在一定程度上保证了解的空间分辨率，但图像重建质量依赖于正 则化参数的选取；当正则化参数较小时，很难改善病态性；当正则化参数较大时，虽然可 以减小对误差的敏感性，但其解通常偏离真实值，甚至可能导致所求解无意义。目前常用 的正则化参数确定方法包括：Tikhonov先验估计、Morozov偏差原理、Arcangeli准则、广 义Arcangeli准则、L-曲线准则等。由于电学层析成像技术对实时性要求较高，而按上述 策略选择正则化参数的计算量较大，因此，实际应用中主要通过经验选择正则化参数，具 有很大的局限性。

为了提高灵敏度分布的均匀性与空间分辨率，通常采用有限元细化的方法。有限元细 化的方法大体可分为两类：h-细化与p-细化。其中h-细化可分为单元细分、网格重构、r- 细化，p-细化可分为多项式阶次在所有区域同步增加以及多项式阶次在局部区域逐次升阶 增加。上述细化方法虽然可以提高精度，但同时存在不足之处：单元细分需要在悬空节点 处施加局部约束，网格重构的计算量很大，而且从一种网格划分到另一种网格划分传输数 据时也存在一些问题，r-细化通常达不到理想精度，目前仍停留在理论研究阶段，而p-细 化无法事先预测出多项式的阶次，一般需要多次的再计算，因而导致更多的计算成本。

发明内容

针对上述现存的技术问题，本发明提供一种基于模型细化的改进牛顿-拉夫逊ERT图像 重建法，能有效提高图像重建质量。

本发明为实现上述目的，通过采用在每个三角形有限元的形心位置增加一个节点的方 法，细化有效提高正问题计算精度的有限元模型b，并在算法重建过程中，遵循“计算正问 题时采用细化前有限元模型b，修正电阻率分布时采用细化后的有限元模型c及其对应的灵 敏度矩阵”的原则，其具体步骤是：

①按传统等间隔剖分原理建立初始有限元模型a；

②以初始有限元模型a中除最外层之外的每一层半径为变量，以敏感场均匀分布时模 型均方根值的倒数为适应度函数，并引入三角形最长边与最短边的比值作为惩罚函数，利 用改进遗传算法离线优化初始有限元模型a，得到优化后的有元限模型b；

③细化步骤②中优化后的有限元模型b，得到细化后的有元限模型c，并完成细化后 的有限元模型c的节点与有限元的编号；

④离线计算敏感场均匀分布时细化后的有限元模型c对应的灵敏度矩阵S；

⑤离线计算(S^TS+μ^(k)E)^-1S^T，其中k为算法迭代次数、μ为正则化因子、E为单位 阵、T为转置矩阵；

⑥取边界电压测量值V₀，并利用线性反投影算法重建图像，将其重建结果作为改进牛 顿-拉夫逊算法初始电阻率分布R⁽⁰⁾；

⑦利用步骤②中优化后的有限元模型b计算正问题，得边界电压计算值V^(k)；

⑧计算误差e=1/2(||V^(k)-V₀||₂)²，从算法迭代次数与误差两方面判断是否满足算法结 束条件，若满足，算法结束并显示重建结果；否则利用细化后的有限元模型c及其对应的 灵敏度矩阵S修正电阻率分布R^(k+1)=R^(k)-(S^TS+μ^(k)E)^-1S^T(V^(k)-V₀)，并跳至步骤⑦循 环计算。

其中，所述的细化步骤②中的优化后有限元模型是在可有效提高正问题计算精度的有 限元模型即模型b的基础上，通过采用在每个三角形有限元形心位置增加一个节点的方法 实现的；所述的细化之后的有限元模型，原节点编号不变，新增节点编号N与其所在的原 有限元编号n及原节点总数total的关系满足：N=n+total；所述的细化之后的有限元模型， 遵循“由内而外，逆时针旋转”的原则，根据新增节点所在的原有限元编号n，对包含此新 增节点的三个有限元分别编号3n-2、3n-1、3n。

原有限元编号n以及编号为3n-2、3n-1、3n三个有限元的关系满足： r(n)=R(3n-2)=R(3n-1)=R(3n)，其中，r(n)表示n的电阻率，R(3n-2)、R(3n-1)、R(3n) 分别表示3n-2、3n-1、3n的电阻率。

本发明的有益效果是：充分利用了优化后的有限元模型b与细化后的有限元模型c的 优点，在有效提高了正问题计算精度的基础上，不仅提高了灵敏度分布的均匀性，同时改 善了Hessian矩阵的病态性，在不影响算法实时性的前提下，有效提高了图像重建质量。

附图说明

图1为本发明的基于模型细化的改进牛顿-拉夫逊ERT图像重建法流程图；

图2为有限元模型a即传统按等间隔原理剖分的有限元模型的示意图；

图3为有限元模型b即利用改进遗传算法离线优化有限元模型a得到的有限元模型的 示意图；

图4为有限元模型c即对有限元模型b细化后得到的有限元模型的示意图；

图5为三种不同模型在算法迭代过程中Hessian矩阵条件数的比较示意图；

图6为三种不同模型在敏感场均匀分布时正问题计算精度的比较示意图。

具体实施方式

下面结合附图对本发明的基于模型细化的改进牛顿-拉夫逊ERT图像重建法进一步说 明。

本发明的基于模型细化的改进牛顿-拉夫逊ERT图像重建法，充分利用了优化后有限元 模型与细化之后的有限元模型的优点，在有效提高了正问题计算精度的基础上，不仅提高 了灵敏度分布的均匀性，同时改善了Hessian矩阵的病态性，在不影响算法实时性的前提 下，有效提高了图像重建质量，该算法的具体步骤是：

1、如图2所示，按传统等间隔剖分原理建立初始有限元模型a；

2、如图3所示，以初始有限元模型a中除最外层之外的每一层半径为变量，以敏感场 均匀分布时模型均方根值的倒数为适应度函数，并引入三角形最长边与最短边的比值作为 惩罚函数，利用改进遗传算法离线优化初始有限元模型a，得到优化后的有元限模型b；

3、如图4所示，通过采用在每个三角形有限元形心位置增加一个节点的方法，对优化 后的有限元模型b进行细化，得到细化后的有元限模型c，并完成细化后的有限元模型c的 节点与有限元的编号，编号方式如下：

1）原节点编号不变，新增节点编号N与其所在的原有限元编号n及原节点总数total的 关系满足：N=n+total；

2）遵循“由内而外，逆时针旋转”的原则，根据新增节点所在的原有限元编号n，对 包含此新增节点的三个有限元分别编号3n-2、3n-1、3n；

4、离线计算敏感场均匀分布时细化后的有限元模型c对应的灵敏度矩阵S；

5、离线计算(S^TS+μ^(k)E)^-1S^T其中k为算法迭代次数、μ为正则化因子、E为单位阵、 T为转置矩阵；

6、取边界电压测量值V₀并利用线性反投影算法重建图像，将其重建结果作为改进牛顿 -拉夫逊算法初始电阻率分布R⁽⁰⁾；

7、利用优化后的有限元模型b计算正问题，得边界电压计算值V^(k)；

8、计算误差e=1/2(||V^(k)-V₀||₂)²，从算法迭代次数与误差两方面判断是否满足算法结 束条件，若满足，算法结束并显示重建结果，否则利用细化后的有限元模型c及其对应的 灵敏度矩阵S修正电阻率分布R^(k+1)=R^(k)-(S^TS+μ^(k)E)^-1S^T(V^(k)-V₀)，并跳至步骤7循 环计算，其中计算正问题时采用的有限元模型b中编号为n的有限元电阻率r(n)，与细化后 的有限元模型c中编号为3n-2、3n-1、3n的三个有限元电阻率R(3n-2)、R(3n-1)、R(3n) 的关系满足：r(n)=R(3n-2)=R(3n-1)=R(3n)。

如图5所示，在算法迭代过程中，图3中细化后的有限元模型c对应的Hessian矩阵 条件数最小，与图1中有限元模型a与图2中有限元模型b相比，有效的改善了Hessian 矩阵病态程度。

如图6所示，采用图1中有限元模型a、图2中有限元模型b和图3中有限元模型c计 算正问题时，敏感场均匀分布时RMS值分别为7.2045%、1.4428%、1.4428%，与图1中有限 元模型a相比，图2中有限元模型b和图3中有限元模型c均能够有效提高正问题的计算 精度；而图2中有限元模型b的有限元数目仅为图3中有限元模型c的有限元数目的1/3， 故采用图2中有限元模型b计算正问题时比采用图3中有限元模型c计算正问题时所需的 时间更少。

经验证，图1中有限元模型a、图2中有限元模型b和图3中有限元模型c所对应的灵 敏度矩阵的最大值与最小值的差值分别为0.0722、0.0640、0.0213，与图1中有限元模型 a和图2中有限元模型b相比，图3中的有限元模型c有效提高了灵敏度分布的均匀性。

经过在基于三种不同模型的基础上对相同图像进行图像重建的验证，与基于修正牛顿- 拉夫逊算法在图1中有限元模型a与图2中有限元模型b下的重建结果相比，基于模型细 化的改进牛顿-拉夫逊ERT图像重建法进行图像重建得到的重建图像边缘更清晰、保真度更 高，特别适用于日出流型与十字流型等复合流型。

综上所述，计算正问题时采用细化前的有限元模型b，修正电阻率分布时采用细化后的 有限元模型c及其对应的灵敏度矩阵，能够确保在有效提高了正问题计算精度的基础上， 不仅提高了灵敏度分布的均匀性，同时改善了Hessian矩阵的病态性，在不影响其实时性 的前提下，有效提高了图像重建质量。