CN110501587A - 一种基于电力线分布的Radon逆变换电容层析成像方法 - Google Patents

Abstract

本发明涉及一种基于电力线分布的Radon逆变换电容层析成像方法，主要包括以下步骤：对于N电极电容层析成像系统，首先循环激励测量得到传统激励模式下空场和有物体条件下的N(N‑1)/2个独立边界测量值，并等效计算得到多电极激励模式下边界测量值，进而得到N×(N‑k)的电容比值矩阵；然后计算场域内电力线曲线方程及不同电力线交点坐标，并建立其与平行线交点的映射；利用Radon逆变换重建得到不同激励电极数目的等间距网格划分场域的重建值，累加不同激励模式重建值，将其映射到对应电力线交点，输出重建结果。该方法将硬场经典重建方法与多电极激励策略相结合，用于电学成像领域物质分布重建，能有效重建出被测区域中介质分布，在电学成像领域具有重要应用价值。

Description

一种基于电力线分布的Radon逆变换电容层析成像方法

技术领域

本发明涉及一种基于电力线分布的Radon逆变换电容层析成像方法，属于电学成像领域。将硬场重建经典算法Radon逆变换应用于软场重建，并与多电极激励策略相结合，增加了测量数据量，实现了感兴趣区域内介电常数分布重建。

背景技术

于20世纪80年代发展起来的电容层析成像(electrical capacitancetomography，ECT)技术以其快速、安全、非接触测量等优点广泛应用于等石油、化工、电力、冶金、建材等工业过程中的参数检测和过程监测。通过扫描空间敏感阵列获取被测物场在激励下的电信号，调制解调后的数据采用图像重建算法反演出被测场域内物质分布信息。电容层析成像系统由三部分构成：(1)信息获取单元，即空间敏感阵列，根据物场实际分布得到边界测量值；(2)信息处理单元：将得到的电信号进行调制解调和滤波；(3)信息恢复单元：通过对逆问题求解反演被测区域物质的空间分布(赵玉磊等.电容层析成像技术的研究进展与分析[J].仪器仪表学报,2012,33(8).)。

ECT重建是一个典型的欠定、非线性、病态性的逆问题求解过程，主要表现在：(1)空间敏感阵列测得的信息量远小于待求解的未知量，故该问题的解不唯一；(2)传感器内部被测区域的敏感场分布受介质分布的变化而变化，具有典型的“软场”特性；(3)反问题的求解本质上是微商求解问题，故边界测量值的微小数据波动会使求解出的物质分布有很大误差。故逆问题的求解具有严重的不适定性(王化祥等.电学层析成像[M].2013.)。常用的ECT图像重建算法大致可分为基于灵敏度矩阵的图像重建算法和直接重建算法两类。

1972年Lehr发表在《IEEE生物医学工程会刊》(IEEE Transactions onBiomedical Engineering)19卷2期，156-157页上题为《一种用于阻抗体积场的向量推导方法》(A vector derivation useful in impedance plethysmographic fieldcalculations)的文章中，在四端口网络的阻抗灵敏度表达基础上，推导了一种更简洁、更具有普适性的阻抗灵敏度表达计算式。典型基于灵敏度矩阵的图像重建算法主要包括线性反投影算法(linear back projection,简称LBP)、Tikhonov正则化法、兰德韦伯算法(Landweber)等。LBP是最早提出的图像重建算法，最早由Barber和Brown于1983年发表在《电子通讯》(Electronic Letters)19卷22期，933-935页《应用电位层析成像技术重建电阻率空间分布》(Imaging spatial distributions of resistivity using appliedpotential tomography—APT)一文中提出并加以验证，是一种基于等位线反投影的定性反投影重建算法，该算法忽略了敏感场分布的变化，沿着投影域将所有的测量结果反投影回去，将经过某点的所有投影射线进行累加进而得到该点的密度估计值，该方法对简单分布的重建结果较好，对复杂分布的中间质量一般。2000年Peng等人发表在《粒子与粒子系统特性》(Particle&Particle Systems Characterization)第17卷第3期题为《采用正则化方法用于电容层析成像图像重建》(Using regularization methods for imagereconstruction of electrical capacitance tomography)一文中提出了标准Tikhonov正则化法，通过引入正则化参数来对逆问题进行求解，但由于该方法过度光滑，因不光滑信息的缺失导致重建图中细节信息丢失、空间分辨率不高。Landweber方法最早是由Landweber于1951年发表在《美国数学杂志》(American journal of mathematics)第73卷，615-624页，题为《一类自由保有积分方程的迭代公式》(An iteration formula forfredholm integral equations of the first kind)提出；Yang等人于1999年发表在《测量科学与技术》(Measurement science and technology)10卷11期，第1065-1069页，题为《一种基于Landweber迭代方法的电容层析成像图像重建算法》(An image-reconstructionalgorithm based on Landweber's iteration method for electrical-capacitancetomography)一文中，通过最速下降法求解目标函数的极小值进而对逆问题进行求解，但该方法重建速度较慢。

直接图像重建算法因其无需计算灵敏度矩阵、可以独立直接地反演出被测区域内任一点上的介电常数分布，能显著提高图像重建的灵活性和实时性，近年来得到了广泛应用，常用的有Calderon算法、D-bar方法等。1980年Calderon在《计算应用和数学》(Computational and applied mathematics)发表的题为《关于逆边界问题》(On aninverse boundary value problem)的文章中提出了Calderon算法，通过所提出的线性化准则成功重建出低对比度的电导率分布。2013年Cao等人发表在《仪表学报》(Journal ofInstrumentation)8卷3期，题为《利用Calderon方法和相邻电流模式对人体胸部进行二维图像重建》(2D image reconstruction of a human chest by using Calderon’s methodand the adjacent current pattern)文章中，利用保角变换将单位圆映射到非圆区域，并利用Caledron方法成功重建出非圆区域内介质分布。D-bar方法是采用非线性傅里叶变换重建被测区域内介质分布，Zhao等人于2018年《IEEE仪器与测量汇刊》(IEEE Transactionson instrumentation and measurement)，题为《一种快速D-bar的直接电容层析成像重建算法》(Direct image reconstruction for electrical capacitance tomography usingshortcut D-bar method)提出了改进的快速D-bar方法，该方法可以直接独立地恢复感兴趣区域内各点的介电常数，具有较好的计算性能。

Radon逆变换是层析成像的基础，指出了二维或三维的物体可以通过其无限或连续的投影反演得到。Gregory于1987年在《IEEE声学语音信号处理学报》(IEEETransactions on acoustics speech&signal processing)第35卷，第162-172页，题为《离散Radon变换》(Discrete radon transform)中推导了将逆问题归结为一个线性代数问题，推导了Radon逆变换的离散实现形式，并将其应用于多变量情况。Clackdoyle等人在2004年发表在《逆问题》(Inverse problems)第20卷第4期，第1281-1291页题为《二维Radon变换的反演公式》(A large class of inversion formulae for the 2D Radon transform offunctions of compact support)中详细介绍了二维Radon变换的几种替代公式，不同于原始Radon逆变换公式，所提出的公式都需要目标函数的信息。

考虑到Radon逆变换重建图像时要求数据量足够多，而电容层析成像传感器在实际制作时会受到测量电路、测量精度、传感器信噪比、加工工艺等多方面的限制，无法布置足够多的电容。基于以上背景，本发明提出一种基于电力线分布的Radon逆变换电容层析成像方法。将层析成像算法的基础Radon逆变换与多电极激励策略相结合，能够获取更多边界电容信息，大大增多了有效数据量，重建速度快且重建精度高。

发明内容

本发明的目的在于提出一种基于电力线分布的Radon逆变换电容层析成像方法，将经典硬场重建方法Radon逆变换应用于软场重建，所提出的方法重建速度快、精度高，能有效重建被测区域内介电常数分布，是一种有效新颖的电容层析成像图像重建方法。

本发明的技术方案是：

步骤一、对于N个电极均匀分布排列的圆形电容层析成像传感器，按照顺时针顺序分别将N个电极标号为r(1≤r≤N)。采用传统激励模式，即在N个均匀分布排列的电极中选择一个电极并施加幅值为1V的激励电压，其余N-1个电极接地，测量激励电极与接地电极之间的电容值。在一次循环激励测量过程中，具体实现方式为：(1)在电极1上施加幅值为1V的激励电压，其余N-1个电极均接地或与地同电位，可分别测量得到激励电极1与电极2到电极N共N-1个电容值；(2)在电极2上施加幅值为1V的激励电压，其余N-1个电极均接地或与地同电位，可测量得到激励电极2与电极3到电极N共N-2个电容值；(3)在电极3上施加幅值为1V的激励电压，其余N-1个电极均接地或与地同电位，可测量得到激励电极3与电极4到电极共N-3个电容值，以此类推；(N)在电极N-1上施加幅值为1V的激励电压，其余N-1个电极均接地或与地同电位，可测量得到激励电极N-1与电极N共1个电容值。传统激励模式得到的N(N-1)/2个独立测量值可以表示成如下形式的电容矩阵：

式中：C_s,t为电极s和电极t之间的电容(s≠t)；C_s,s为第s个电极上的自电容，其大小为该电极与其他电极之间的电容总和，即：

采用多电极激励模式，即改变激励电极数，在k(2≤k≤N/2)个相邻的电极上施加激励，其余N-k个电极接地或与地同电位。其具体实现方法为：(1)电极1到电极k共k个电极同步施加电压激励，其余N-k个电极接地或与地同电位，可分别测量得到激励电极组与电极k+1到电极N共N-k个电容值；(2)电极2到电极k+1共k个电极同步施加激励，其余N-k个电极接地或与地同电位，可测量得到激励电极组与电极k+2到电极N、电极1共N-k个电容值；(3)电极3到电极k+2共k个电极同步施加激励，其余N-k个电极接地或与地同电位，可测量得到激励电极组与电极k+3到电极N、电极1、电极2共N-k个电容值，以此类推；(N)电极N到电极k-1共k个电极同步施加激励，其余N-k个电极接地或与地同电位，可测量得到激励电极组与电极k到电极N-1共N-k个电容值。该种激励模式共能得到N(N-k)/2个独立测量值

多电极激励模式(k≥2)的电容矩阵可表示为式中，列向量为k个电极同时激励时电极s作为激励电极对起始电极所测得电容值，可通过传统激励模式下电容矩阵变换等效成多个电极激励等效得到，即：

其中，为0、1分布的电压激励向量；为k-by-k激励模式下第s个电极激励状态，若第s个电极上施加电压激励，则否则，

步骤二、电容比值矩阵的获取。采用循环电极对激励模式，根据式(1)、(3)分别测量得到k(2≤k≤N/2)个电极激励时空场状态的N(N-k)/2个独立测量值和有物体状态的N(N-k)/2个独立测量值。其中，空场状态指待测场域内充满均匀介质，有物体状态指场域内部存在扰动，利用有物体测量值与空场测量值之比可得到N×(N-k)的电容比值矩阵，如下所示：

其中，Q_{s～s+k-1,s+k+j}(j∈[0 N-k-1])表示激励电压施加在电极s到电极s+k-1共k个电极，距离电极激励对终止电极第j个电极上测得绝对电容值。电容比值矩阵每一行表示在该电极对上施加电压激励时测量电极上得到的N-k个绝对电容测量值。

步骤三、电力线交点计算并建立与平行线交点之间的映射。顺时针均匀排列N个电极的圆形ECT传感器，对于不同激励电极数，场域内的电力线均从激励电极对发出并终止于测量电极，电力线是电场中的“等间距平行线”，圆周上均匀分布的测量电极等同于等间距均匀分布的测量点。循环测量电极坐标便可得到某一特定激励电极对所对应的电力线曲线函数，循环激励电极对便可得到不同激励模式下所有电力线曲线函数表达式，并计算得到不同激励电极对电力线的交点坐标。利用等间距均匀分布测量点生成等间距平行线并计算得到平行线交点坐标，建立实际电力线交点坐标与等间距平行线交点坐标之间的映射关系。

步骤四、Radon逆变换对电容比值矩阵重建。对每一个电极激励对得到的绝对测量值可利用Radon逆变换重建得到被测区域内介电常数分布，其解析重建公式如下所示：

其中，(R,θ)为待测场域内重建点极坐标参数，ε(R,θ)表示极坐标(R,θ)处重建得到的介电常数，p表示测量电极的坐标，φ表示激励电极的角度。N电极圆形电容层析成像传感器的电极数及电极坐标排布范围有限，因此可将式(5)离散为如下形式：

其中，φ_r表示电极r(1≤r≤N)的角度，即L1和L2分别为测量电极坐标分布的上下界限。

式中对离散的绝对电容值求导可通过先将离散绝对电容值拟合为关于测量电极坐标p的连续多项式，再进行一阶求导实现，现以三次拟合为例，则式(6)可写为：

其中a_j为不同阶次的拟合系数。

测量电极等间距排列在区间[L1,L2]范围内，根据测量电极坐标采用分段积分形式可将式(7)分解为：

对不同的激励模式得到N*(N-k)阶的电容比值矩阵可根据式(8)得到被测区域重建值ε_k(R,θ)。

步骤五、判断激励电极数k是否达到激励电极数上限N/2，如果是，则执行步骤六；如果否，则激励电极数k加1并执行步骤二。

步骤六、停止迭代，计算N/2种典型激励模式重建值总和ε_sum(R,θ)，计算公式如下所示：

根据步骤三所建立等间距平行束交点与电力线交点之间的映射关系，将重建值ε_sum(R,θ)映射到对应电力线交点处，并输出重建图像。

本发明的有益效果为：以典型的Radon逆变换为基础，将多电极激励策略与典型硬场重建算法Radon逆变换相结合，实现了被测区域内介电常数的快速且精确重建。

附图说明

图1是本方法的具体实施图。

图2是16电极圆形电容层析成像传感器。

图3具体实施原始模型。

图4是传统激励模式下电力线分布图。

图5基于电力线分布的Radon逆变换重建结果图。

具体实施方式

下面结合如图2所示的16电极圆形电容层析成像传感器和如图3所示的具体实施模型对本发明作进一步说明。

步骤一、圆形电容层析成像传感器内壁均匀分布有16个电极，如图2所示，将16个电极按照顺时针分别标号为r(1≤r≤16)。采用16电极圆形电容层析成像传感器对如图3所示的原始分布模型采用传统激励模式进行测量。其具体实现方式为：(1)在电极1上施加幅值为1V的激励电压，其余15个电极均接地或与地同电位，可分别测量得到电极1与电极2到电极16共15个电容值；(2)在电极2上施加幅值为1V的激励电压，其余15个电极均接地或与地同电位，可分别测量得到电极2与电极3到电极16共14个电容值；(3)在电极3上施加幅值为1V的激励电压，其余15个电极均接地或与地同电位，可分别测量得到电极3与电极4到电极16共13个电容值，以此类推；(16)在电极15上施加幅值为1V的激励电压，其余15个电极均接地或与地同电位，测量得到电极15与电极16共1个电容值。

传统激励模式，即激励电极数k＝1，测量得到的120个独立测量值可记为：

其中：C_s,t为电极s和t之间的电容(s≠t)；C_s,s为第s个电极上的自电容，其值为该电极与其他电极之间的电容总和，即：

改变激励电极数，其中多电极激励模式(k≥2)的电容矩阵可通过单电极激励模式电容矩阵等效得到，列向量为k个电极激励时第s个电极作为激励电极对起始电极所测得电容值，计算式如下：

式中：电压激励向量为0、1分布向量；为k个电极同时激励时第s个电极激励状态，若第s个电极上存在激励电压，则否则，

每种激励模式对应的独立测量值如下表所示。

表1不同激励模式下的独立测量值个数

激励模式	独立测量值个数	激励模式	独立测量值个数
				1-by-1	16×(16-1)/2＝120	5-by-5	16×(16-5)/2＝88
2-by-2	16×(16-2)/2＝112	6-by-6	16×(16-6)/2＝80
				3-by-3	16×(16-3)/2＝104	7-by-7	16×(16-7)/2＝72
4-by-4	16×(16-4)/2＝96	8-by-8	16×(16-1)/2＝64

步骤二、获取电容比值矩阵。采用k电极激励模式分别测量空场条件下的16(16-k)/2个独立测量值和有物体时的16(16-k)/2个独立测量值。利用有物体分布电容矩阵和空场电容矩阵的比值可得到16*(16-k)阶的电容比值矩阵，如下所示：

其中，Q_{s～s+k-1,s+k+j}(j∈[0 N-k-1])为电极s及其顺时针方向的k个相邻电极施加电压激励时，距离电极激励对终止电极第j个电极测量得到的绝对电容测量值。

步骤三、计算电力线交点并与平行线交点建立映射。顺时针均匀分布16个电极的圆形ECT传感器，整个被测区域内，电力线从激励电极对发出并终止于测量电极，认为是电场中的“等间距平行线”，圆周上均匀分布的测量电极等同于等间距分布测量点。计算得到不同数目激励电极所生成的电力线曲线方程，并计算得到电力线交点坐标。根据等间距分布测量点生成等间距平行线并计算平行线交点坐标，建立实际电力线交点坐标与等间距平行线交点坐标之间的映射关系。

步骤四、利用Radon逆变换反演场域介电常数分布。对激励模式下得到的绝对电容测量值利用Radon逆变换重建得到坐标点(R,θ)处的介电常数，解析重建公式如下所示：

其中，待重建点用极坐标参数表示为(R,θ)，在该点重建得到的介电常数记为ε_k(R,θ)，p表示测量电极坐标分布，φ_r表示激励电极r的角度分布，可写为

式中：L1和L2分别为测量电极分布的界限。积分项中对离散绝对电容值求导可通过先将离散电容值拟合为关于p的连续多项式，再对p求导实现，现以三次拟合为例，根据测量电极坐标采用分段积分形式可将式(14)简化为：

对不同的激励模式得到16×(16-k)的电容比值矩阵可根据式(15)重建得到被测区域重建值ε_k(R,θ)。

步骤五、判断激励电极数k是否达到激励电极数上限8，如果是，则执行步骤六；如果否，则执行步骤二。

步骤六、停止迭代，根据下式计算8种典型激励模式重建值总和ε_sum(R,θ)：

通过步骤三所建立平行束交点与电力线交点之间的映射关系，将重建值ε_sum(R,θ)映射到对应电力线交点处，并输出重建图像。

Claims (5)

1.一种基于电力线分布的Radon逆变换电容层析成像方法，其特征在于将电容层析成像传感器的独立测量值转化得到不同激励模式下的测量值，利用Radon逆变换对不同激励模式下测量值进行被测区域内介电常数分布重建，实现感兴趣区域内介电常数分布重建。

2.根据权利要求1所述的一种基于电力线分布的Radon逆变换电容层析成像方法，其特征在于：对于N个电极均匀分布排列的圆形电容层析成像传感器，按照顺时针顺序分别将N个电极标号为r(1≤r≤N)，采用传统激励测量模式循环改变激励电极，其具体实施方案为：(1)在电极1上施加激励，其余N-1个电极接地或与地同电位，可分别测量得到电极1与电极2到电极N共N-1个电容值；(2)在电极2上施加激励，其余N-1个电极接地或与地同电位，可分别测量得到电极2与电极3到电极N共N-2个电容值；(3)在电极3上施加激励，其余N-1个电极接地或与地同电位，可分别测量得到电极3与电极4到电极N共N-3个电容值；(N)以此类推，在电极N-1上施加激励，其余N-1个电极接地或与地同电位，测量得到电极N-1与电极N共1个电容值，通过传统激励测量模式测量得到N(N-1)/2个独立测量值，可表示成如式(1)所示的电容矩阵：

其中：C_s,t为电极s和t之间的电容(s≠t)；C_s,s为第s个电极上的自电容，其值为该电极与其他电极之间的电容总和，即：

采用多电极激励模式(2≤k≤N/2)，即在k个相邻的电极上施加激励，其余N-k个电极接地或与地同电位。其中，k电极激励模式的具体实现方法为：(1)电极1到电极k共k个电极同步施加激励，其余N-k个电极接地或与地同电位；(2)电极2到电极k+1共k个电极同步施加激励，其余N-k个电极接地或与地同电位；(3)电极3到电极k+2共k个电极同步施加激励，其余N-k个电极接地或与地同电位，以此类推；(N)电极N到电极k-1共k个电极同步施加激励，其余N-k个电极接地或与地同电位。该种激励模式共能得到N(N-k)/2个独立测量值，循环激励测量得到多电极激励模式下的电容矩阵式中，列向量为k个激励模式下第s个电极作为激励电极对起始电极所测得电容值，多电极激励模式的电容矩阵可由传统激励模式的电容矩阵通过变换等效得到，计算公式如下：

式中：为0、1分布的电压激励向量；为k个电极激励模式下第s个电极激励状态，若第s个电极上施加激励，则否则，

3.根据权利要求1所述的一种基于电力线分布的Radon逆变换电容层析成像方法，其特征在于：计算电力线交点并与平行束交点之间建立映射，对于N电极的圆形ECT(electricalcapacitance tomography)系统，对于不同激励模式，电力线均从激励电极对发出终止于测量电极，电力线是电场中的“等间距平行线”，圆周上均匀分布的测量电极等同于等间距分布测量点，循环测量电极得到某固定激励电极对时所有电力线曲线函数，采用循环激励策略进而获得不同激励模式下所有电力线曲线函数表达式，并计算所有电力线交点坐标。根据等间距分布测量点生成等间距平行线并计算平行线交点坐标，建立实际电力线交点坐标与等间距平行线交点坐标之间的映射关系。

4.根据权利要求1所述的一种基于电力线分布的Radon逆变换电容层析成像方法，其特征在于：获取不同激励模式下的电容比值矩阵，采用循环激励测量模式，根据式(1)、(3)分别计算k个电极激励时空场条件的N(N-k)/2个独立测量值和有物体条件的N(N-k)/2个独立测量值；其中，空场指待测场域内部充满均匀分布的介质，有物体指待测场域内部存在扰动，利用有物体测量值与空场测量值的比值可得到N*(N-k)阶的电容比值矩阵，可表示为：

其中，Q_{s～s+k-1,s+k+j}(j∈[0 N-k-1])为在电极s到电极s+k-1共k个电极同步施加激励，距离电极激励对终止电极第j个电极上测得绝对电容值，电容比值矩阵的每一行表示该电极对作为激励电极对时测量电极上的绝对电容值。

5.根据权利要求1所述的一种基于电力线分布的Radon逆变换电容层析成像方法，其特征在于：利用Radon逆变换方法根据电容比值矩阵进行重建，对k个电极激励模式下的电容比值矩阵Q_N×(N-k)利用Radon逆变换方法重建被测场域内介电常数分布，计算式为：

其中，p表示测量电极的坐标排列分布，φ表示激励电极的角度分布，Q_N×(N-k)(p,φ)为电容比值向量，(R,θ) 为待测场域内重建点极坐标参数，ε(R,θ)表示极坐标为(R,θ)点处的重建介电常数，实际上，激励电极个数及测量电极坐标排布范围均有限，故可将式(5)写成如下形式：

式中：L1和L2分别为测量电极分布的上下界限。积分项中对离散绝对电容值求导可通过先将离散电容值拟合为关于p的连续多项式，再对p求导实现，现以三次拟合为例，则式(6)可写为：

式中：a_j为不同阶次的拟合系数；测量电极等间距排列在区间[L1,L2]范围内，根据测量电极坐标采用分段积分形式可得：

对不同的激励模式得到N×(N-k)的电容比值矩阵可根据式(8)得到被测区域重建值ε_k(R,θ)。

直到激励电极数k达到激励电极数上限N/2，停止迭代，根据式(9)计算N/2种典型激励模式重建值总和ε_sum(R,θ)：

根据权利要求3建立的电力线交点并与平行束交点之间的映射，将重建值ε_sum(R,θ)映射到电力线交点，并输出重建图像，完成感兴趣区域重建。