# A New Iterative Algorithm Based on Correction of Sensitivity Matrix for Electrical Resistance Tomography

Authors: Yutong Chen; Yan Han; Wuqiang Yang; Kun Li
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

Electrical resistance tomography (ERT) is used to reconstruct the resistance/conductivity distribution. Usually, a uniform distribution is assumed as the initial condition to obtain a generic sensitivity matrix, which may be very different from a theoretical sensitivity matrix, resulting in a large error. The aim of this study is to analyse the difference between a generalized sensitivity matrix and a theoretical sensitivity matrix and to improve image reconstruction. The effect of the generic sensitivity matrix and theoretical sensitivity matrix on image reconstruction is analyzed. The error caused by the use of the generic sensitivity matrix is estimated. To reduce the error, an improved iterative image reconstruction algorithm is proposed, which is based on calculation of the error between the generic sensitivity matrix and the theoretical sensitivity matrix, and a correction coefficient with a penalty. During the iterative process, the resistivity distribution and sensitivity matrix are alternatively corrected. Simulation and experimental results show that the proposed algorithm can improve the quality of images, e.g., of two-phase distributions.


---

## Body

## 1. Introduction

As a fast, nonradiation, nonintrusive, and low-cost imaging technique, electrical resistance tomography (ERT) [1] has been developed for visualizing the internal conductivity distribution of industrial processes, by providing 2D/3D images. ERT has been widely used for chemical [2], food processing [3], biology [4], geology [5], and other industrial applications [6]. In ERT, usually a current is injected to a pair of electrodes and the boundary voltages are measured from other pairs of electrodes [7]. A particular difficulty for ERT is that its “soft field” [8] makes the inverse problem be nonlinear and ill posed, and hence the spatial resolution is limited. Various image reconstruction algorithms have been developed for ERT, such as linear back-projection (LBP) [9], Landweber iteration [10], and Tikhonov regularization [11]. Usually, a uniform distribution is assumed as the initial condition to obtain a generic sensitivity matrix.Recently, several modified algorithms have been proposed. For example, Zhang [12] and Xiao et al. [13], respectively, used the results obtained by Landweber iteration and Newton–Raphson algorithms as the initial condition for updating a sensitivity matrix. These modified algorithms can reduce the reconstruction error caused by a generic sensitivity matrix and improve the image quality.In this paper, the effect of the generic sensitivity matrix on image quality is analyzed, and a new image reconstruction algorithm is proposed, which is based on calculating the error introduced by the generic sensitivity matrix according to the variation in the boundary voltage. This algorithm includes a penalty factor, which can correct the sensitivity matrix during the iterative image reconstruction process.


## 2. Principle of ERT

### 2.1. Mathematical Model

An ERT sensor is shown in Figure1 [14]. It shows that a current is injected into a pair of adjacent electrodes and the voltages between the remaining pairs of adjacent electrodes are measured. The injection and measurement are sequentially conducted until each adjacent electrode pair is used for current injection. For a 16-electrode ERT sensor, 104 independent measurements can be obtained.Figure 1
16-electrode ERT sensor model. (a) Uniform conductivity. (b) Perturbed conductivity.
(a)
(b)When a perturbationΔσ≪σ0 (uniform conductivity) is small, according to the sensitivity coefficient principle, a linearised model for current injection and voltage measurement (see Figure 1) is given by [15](1)S⋅Δσ=ΔZ,where ΔZ=ΔZ1,ΔZ1,…,ΔZMZT is a MZ×1 vector of change in measured boundary voltage, Δσ=Δσ1,Δσ2,…,ΔσMσT is a Mσ×1 vector of change in conductivity, and S=SijMZ×Mσ is a calculated sensitivity matrix, where Sij is the sensitivity coefficient of the jth pixel relative to the ith measurement.Once the change in boundary voltageΔZ and the sensitivity matrix S have been obtained, a vector of change in conductivity Δσ can be calculated, e.g., using Landweber iteration algorithm [10](2)Δσk+1=Δσk+α⋅STΔZ−S⋅Δσk,Δσ1=STΔZ,where Δσk is a vector of change in conductivity at the kth iteration and α is the step length.FromΔσ, the conductivity distribution vector σ can be obtained by(3)σ=Qσ0+Δσ,where Q is the projection operator [16], so that the conductivity distribution is limited between 0 and 1.


### 2.2. Effect of Sensitivity Matrix on Image Reconstruction

In equation (2), the sensitivity matrix S changes with the conductivity distribution σ, causing serious ill posedness and nonlinearity, as shown in Figure 2. σ can be solved by Fσ=Z [17], where Z is the boundary voltage vector. It is common that the sensitivity matrix is used to estimate the boundary voltage based on the conductivity distribution, which is called forward projection:(4)Sσ=Z.Figure 2
Different conductivity distributions and corresponding sensitivity matrix. (a) Conductivity distribution of D1. (b) Mean of sensitivity matrix of D1. (c) Uniform distribution. (d) Mean of sensitivity matrix of uniform distribution.
(a)
(b)
(c)
(d)In practice, the sensitivity matrixS cannot be used to solve the problem because S changes with the conductivity distribution σ. Therefore, only the uniform field sensitivity matrix S0 can be used to find the approximate solution to the conductivity distribution by equation (1). Figure 3 shows a reconstructed image using a theoretical sensitivity matrix and another image using a generic sensitivity matrix based on a uniform field distribution.Figure 3
Reconstructed image of model D1. (a) Using theoretical sensitivity matrix of D1. (b) Using sensitivity matrix of uniform field.
(a)
(b)Comparing Figure3(a) with 3(b), it can be seen that the area of discrete-phase is enlarged and the edge of two phases is not distinct with artifacts near the boundary, i.e., the image is significantly degraded. The main reason for this is that the sensitivity matrix in equation (1) is only suitable for the case, where the additional perturbation Δσ of the sensing field is small [17] and the error is larger with a larger Δσ.


## 2.1. Mathematical Model

An ERT sensor is shown in Figure1 [14]. It shows that a current is injected into a pair of adjacent electrodes and the voltages between the remaining pairs of adjacent electrodes are measured. The injection and measurement are sequentially conducted until each adjacent electrode pair is used for current injection. For a 16-electrode ERT sensor, 104 independent measurements can be obtained.Figure 1
16-electrode ERT sensor model. (a) Uniform conductivity. (b) Perturbed conductivity.
(a)
(b)When a perturbationΔσ≪σ0 (uniform conductivity) is small, according to the sensitivity coefficient principle, a linearised model for current injection and voltage measurement (see Figure 1) is given by [15](1)S⋅Δσ=ΔZ,where ΔZ=ΔZ1,ΔZ1,…,ΔZMZT is a MZ×1 vector of change in measured boundary voltage, Δσ=Δσ1,Δσ2,…,ΔσMσT is a Mσ×1 vector of change in conductivity, and S=SijMZ×Mσ is a calculated sensitivity matrix, where Sij is the sensitivity coefficient of the jth pixel relative to the ith measurement.Once the change in boundary voltageΔZ and the sensitivity matrix S have been obtained, a vector of change in conductivity Δσ can be calculated, e.g., using Landweber iteration algorithm [10](2)Δσk+1=Δσk+α⋅STΔZ−S⋅Δσk,Δσ1=STΔZ,where Δσk is a vector of change in conductivity at the kth iteration and α is the step length.FromΔσ, the conductivity distribution vector σ can be obtained by(3)σ=Qσ0+Δσ,where Q is the projection operator [16], so that the conductivity distribution is limited between 0 and 1.


## 2.2. Effect of Sensitivity Matrix on Image Reconstruction

In equation (2), the sensitivity matrix S changes with the conductivity distribution σ, causing serious ill posedness and nonlinearity, as shown in Figure 2. σ can be solved by Fσ=Z [17], where Z is the boundary voltage vector. It is common that the sensitivity matrix is used to estimate the boundary voltage based on the conductivity distribution, which is called forward projection:(4)Sσ=Z.Figure 2
Different conductivity distributions and corresponding sensitivity matrix. (a) Conductivity distribution of D1. (b) Mean of sensitivity matrix of D1. (c) Uniform distribution. (d) Mean of sensitivity matrix of uniform distribution.
(a)
(b)
(c)
(d)In practice, the sensitivity matrixS cannot be used to solve the problem because S changes with the conductivity distribution σ. Therefore, only the uniform field sensitivity matrix S0 can be used to find the approximate solution to the conductivity distribution by equation (1). Figure 3 shows a reconstructed image using a theoretical sensitivity matrix and another image using a generic sensitivity matrix based on a uniform field distribution.Figure 3
Reconstructed image of model D1. (a) Using theoretical sensitivity matrix of D1. (b) Using sensitivity matrix of uniform field.
(a)
(b)Comparing Figure3(a) with 3(b), it can be seen that the area of discrete-phase is enlarged and the edge of two phases is not distinct with artifacts near the boundary, i.e., the image is significantly degraded. The main reason for this is that the sensitivity matrix in equation (1) is only suitable for the case, where the additional perturbation Δσ of the sensing field is small [17] and the error is larger with a larger Δσ.


## 3. Error Penalty Reconstruction Based on Landweber Iteration

### 3.1. Error Introduced by Generic Sensitivity Matrix

Assuming that the theoretical sensitivity matrixS is known, the following equation can be obtained from equation (4):(5)S⋅σ−S0⋅σ0=Z−Z0,where Z0 is the boundary voltage vector of uniform field σ0.BecauseS cannot be obtained when σ is unknown, the inverse problem is usually solved by equation (1) using a generic sensitivity matrix S0. Equation (1) can be rewritten as(6)S0⋅Δσ=ΔZ,where ΔZ=Ζ−Ζ0. Then, equation (6) can be expanded to(7)S0σ−S0σ0=Z−Z0.Comparing with equation (5), the theoretical sensitivity matrix S is replaced by the generic sensitivity matrix S0. To analyse the error, equation (7) is rewritten as(8)Sσ−S0σ0+S0−Sσ=Z−Z0,and hence(9)η=S0−Sσ.It is the error vector introduced when the theoretical sensitivity matrix is replaced. For simplicity,η is called “matrix substitution error” in this paper, which can be ignored when Δσ=σ−σ0 is small. However, when the conductivity distribution changes significantly compared with the uniform field, η also becomes larger and would result in a larger reconstruction error.If the conductivity distributionσ˜ can be found, which is close enough to the real distribution σ as the initial distribution, and its corresponding sensitivity matrix S˜ can be used as the initial sensitivity matrix, then(10)S˜σ−S˜σ˜=Z−Z˜.As a result, the errorη=S˜−Sσ can be effectively reduced, improving the image quality.


### 3.2. Sensitivity Matrix Correction with Penalty Factor

Considering Landweber’s solution as the initial distributionσ˜ in equation (10), the distribution is obtained by the generic sensitivity matrix and the size of the error η is different from each other for different distribution. To reduce the error, the change in distribution calculated by conventional Landweber iterative algorithm needs to be constrained differently according to the distribution. In equation (9), it can be seen that if the difference between the generic sensitivity matrix S0 and the theoretical sensitivity matrix S is small, the error η is small and the reconstruction result Δσ˜ is accurate. The greater difference between them would result in the larger error of reconstruction. Therefore, a penalty factor μ associated with η is chosen to penalize the error of Δσ˜. However, the theoretical error η cannot be obtained when the true distribution σ is unknown. Therefore, it is necessary to use other related known quantities to determine the penalty factor μ.In equation (4), the measured boundary voltage vector Z is related to the theoretical sensitivity matrix S and the true distribution σ. Therefore, a method is proposed to estimate the error η using the change in boundary voltage ΔZ. The relative changes in the boundary voltage is(11)μ=ΔZZ0.It is set as the penalty factor. From equations (5) and (9), it can be seen that μ would increase as η increases. A correction coefficient is defined as(12)β=B−μ.It is used to correctΔσ˜ obtained by Landweber algorithm, where B is a constant. Considering the convergence rate and stability, B can be selected between 0.5 and 1.5. Then, equation (3) can be rewritten as(13)σ˜=Qσ0+β⋅Δσ˜.After the conductivity distribution is corrected, a corrected sensitivity matrixS˜ can be obtained. Then, the corrected distribution σ˜ and sensitivity matrix S˜ can be used as the initial value in equation (10), to further reduce the image error.


### 3.3. Improved Landweber Algorithm

The following is about approximating the true field distribution by multiple corrections. During this process, the error changes as(14)ηp=Sp−S,where p=0,1,…,P is the correction time, Sp is the sensitivity matrix corresponding to the pth distribution corrected. When p = 0, the initial value is corresponding to the uniform field.The correction factor for thep + 1st correction is adjusted as(15)βp+1=B−μp+1,where μp+1=ΔZp/Zp and ΔZp=Z−Zp, Zp is the calculated boundary voltage corresponding to σp after the pth correction. Because μp+1 is related to the error ηp, with the approximation, whether or not further correction of the sensitivity matrix is needed can be determined by μp+1:(16)μp+1<ε,No need to correct the sensitivity matrix,≥ε,Need to correct the sensitivity matrix.The improved reconstruction algorithm is shown in Figure4 with the following steps:Step 1: input the boundary voltage measurementZ and the uniform field conductivity distribution σ0.Step 2: calculate the sensitivity matrixSp and change in boundary voltage value ΔZp using σp, and then calculate the correction coefficient βp+1.Step 3: calculate the approximate increment of conductivity distributionΔσp after M iterations using the conventional Landweber iterate algorithm.Step 4: correctΔσp using βp+1.Step 5: determine whether or not it is necessary to correct the sensitivity matrix. If necessary, return to Step 2. If not, output the conductivity distributionσp+1.Figure 4
Block diagram of improved Landweber iterate algorithm.


## 3.1. Error Introduced by Generic Sensitivity Matrix

Assuming that the theoretical sensitivity matrixS is known, the following equation can be obtained from equation (4):(5)S⋅σ−S0⋅σ0=Z−Z0,where Z0 is the boundary voltage vector of uniform field σ0.BecauseS cannot be obtained when σ is unknown, the inverse problem is usually solved by equation (1) using a generic sensitivity matrix S0. Equation (1) can be rewritten as(6)S0⋅Δσ=ΔZ,where ΔZ=Ζ−Ζ0. Then, equation (6) can be expanded to(7)S0σ−S0σ0=Z−Z0.Comparing with equation (5), the theoretical sensitivity matrix S is replaced by the generic sensitivity matrix S0. To analyse the error, equation (7) is rewritten as(8)Sσ−S0σ0+S0−Sσ=Z−Z0,and hence(9)η=S0−Sσ.It is the error vector introduced when the theoretical sensitivity matrix is replaced. For simplicity,η is called “matrix substitution error” in this paper, which can be ignored when Δσ=σ−σ0 is small. However, when the conductivity distribution changes significantly compared with the uniform field, η also becomes larger and would result in a larger reconstruction error.If the conductivity distributionσ˜ can be found, which is close enough to the real distribution σ as the initial distribution, and its corresponding sensitivity matrix S˜ can be used as the initial sensitivity matrix, then(10)S˜σ−S˜σ˜=Z−Z˜.As a result, the errorη=S˜−Sσ can be effectively reduced, improving the image quality.


## 3.2. Sensitivity Matrix Correction with Penalty Factor

Considering Landweber’s solution as the initial distributionσ˜ in equation (10), the distribution is obtained by the generic sensitivity matrix and the size of the error η is different from each other for different distribution. To reduce the error, the change in distribution calculated by conventional Landweber iterative algorithm needs to be constrained differently according to the distribution. In equation (9), it can be seen that if the difference between the generic sensitivity matrix S0 and the theoretical sensitivity matrix S is small, the error η is small and the reconstruction result Δσ˜ is accurate. The greater difference between them would result in the larger error of reconstruction. Therefore, a penalty factor μ associated with η is chosen to penalize the error of Δσ˜. However, the theoretical error η cannot be obtained when the true distribution σ is unknown. Therefore, it is necessary to use other related known quantities to determine the penalty factor μ.In equation (4), the measured boundary voltage vector Z is related to the theoretical sensitivity matrix S and the true distribution σ. Therefore, a method is proposed to estimate the error η using the change in boundary voltage ΔZ. The relative changes in the boundary voltage is(11)μ=ΔZZ0.It is set as the penalty factor. From equations (5) and (9), it can be seen that μ would increase as η increases. A correction coefficient is defined as(12)β=B−μ.It is used to correctΔσ˜ obtained by Landweber algorithm, where B is a constant. Considering the convergence rate and stability, B can be selected between 0.5 and 1.5. Then, equation (3) can be rewritten as(13)σ˜=Qσ0+β⋅Δσ˜.After the conductivity distribution is corrected, a corrected sensitivity matrixS˜ can be obtained. Then, the corrected distribution σ˜ and sensitivity matrix S˜ can be used as the initial value in equation (10), to further reduce the image error.


## 3.3. Improved Landweber Algorithm

The following is about approximating the true field distribution by multiple corrections. During this process, the error changes as(14)ηp=Sp−S,where p=0,1,…,P is the correction time, Sp is the sensitivity matrix corresponding to the pth distribution corrected. When p = 0, the initial value is corresponding to the uniform field.The correction factor for thep + 1st correction is adjusted as(15)βp+1=B−μp+1,where μp+1=ΔZp/Zp and ΔZp=Z−Zp, Zp is the calculated boundary voltage corresponding to σp after the pth correction. Because μp+1 is related to the error ηp, with the approximation, whether or not further correction of the sensitivity matrix is needed can be determined by μp+1:(16)μp+1<ε,No need to correct the sensitivity matrix,≥ε,Need to correct the sensitivity matrix.The improved reconstruction algorithm is shown in Figure4 with the following steps:Step 1: input the boundary voltage measurementZ and the uniform field conductivity distribution σ0.Step 2: calculate the sensitivity matrixSp and change in boundary voltage value ΔZp using σp, and then calculate the correction coefficient βp+1.Step 3: calculate the approximate increment of conductivity distributionΔσp after M iterations using the conventional Landweber iterate algorithm.Step 4: correctΔσp using βp+1.Step 5: determine whether or not it is necessary to correct the sensitivity matrix. If necessary, return to Step 2. If not, output the conductivity distributionσp+1.Figure 4
Block diagram of improved Landweber iterate algorithm.


## 4. Simulation and Results

### 4.1. Experimental Conditions

MatLab is used for simulation. The configuration of the simulated ERT sensor is listed below:(i)
Circular inner diameter: 50 mm(ii)
Liquid conductivity: 1 S/m(iii)
Gas conductivity: 0 S/m(iv)
Number of electrodes: 16(v)
Excitation and measurement mode: current injection to adjacent electrode pair and voltage measurement from the other adjacent electrode pairs (see Figure1)(vi)
Number of independent measurement: 104(vii)
Number of finite elements: 1536 (see Figure5)To verify the proposed algorithm, three algorithms are used to reconstruct images for the four models, as shown in Figure6. A uniform field is used as the initial distribution. The three algorithms are shown in Table 1.Figure 5
Mesh grid.Figure 6
Four conductivity distribution models for simulation. (a) Model D1 (b) Model D2. (c) Model D3. (d) Model D4.
(a)
(b)
(c)
(d)Table 1
Algorithm table.
Algorithm 1 [10, 16]
Δ
σ
1
=
S
T
Δ
Z
Δ
σ
k
+
1
=
Δ
σ
k
+
α
⋅
S
T
Δ
Z
−
S
⋅
Δ
σ
k
σ
=
Q
σ
0
+
Δ
σ
M
Algorithm 2 [12]
Δ
σ
p
=
Δ
σ
M
S
P
,
Δ
Z
P
Conventional Landweber
σ
p
+
1
=
Q
σ
p
+
β
p
+
1
Δ
σ
p
β
p
+
1
=
1
Algorithm 3
Δ
σ
p
=
Δ
σ
M
S
P
,
Δ
Z
P
Conventional Landweber
σ
p
+
1
=
Q
σ
p
+
β
p
+
1
Δ
σ
p
β
p
+
1
=
B
−
Δ
Z
p
/
Z
pIn Algorithm 1 (uncorrected), the conventional Landweber iterative algorithm [10] is used and the sensitivity matrix is not corrected during the iterative process. In Algorithm 2 [12] (βp+1=1), βp+1=1. M=500, when p=1, and M=100, when p>1. p=0,1,…,P,P=10 and ε=0. In Algorithm 3 (βp+1=1.2−ΔZp/Zp), the algorithm proposed in this paper is used, βp+1=B−ΔZp/Zp,B=1.2. M=500, when p=1, and M=100, when p>1. p=0,1…P,P=10 and ε=0.A reconstruction error and a correlation coefficient are used to evaluate the reconstruction results. The smaller the reconstruction error and the greater the reconstruction correlation means higher image quality [10]:(17)Reconstruction error=σ^−σσ×100%,Correlation coefficient=∑i=1Kσi−σ¯σ^i−σ^¯∑i=1Kσi−σ¯2∑i=1Kσi−σ^¯2,where σ^ is the calculated conductivity distribution, σ is the true conductivity distribution, and K = 1536 is the number of elements in the sensing area.Quantitative analysis is performed by calculating the square sum of the substitution error of the matrix:(18)Matrix substitution errorηp=Sp−S⋅σ.


### 4.2. Results and Analysis

#### 4.2.1. Results with Noise Free Data

Without considering noise, the above three algorithms are used for reconstruction and the results are shown in Figure7. With Algorithm 1, the reconstructed images with 500 iterations, 600 iterations, and 1400 iterations are given. For Algorithms 2 and 3, the reconstructed images after 1, 2, and 10 corrections are given, and the images correspond to iterations 500 times, 600 times, and 1400 times. The third column (βp+1=1.2−ΔZp/Zp) is obtained using Algorithm 3.Figure 7
Reconstruction of three algorithms without noise. (a) Reconstructed image of model D1. (b) Reconstructed image of model D2. (c) Reconstructed image of model D3. (d) Reconstructed image of model D4.
(a)
(b)
(c)
(d)In Figure7, the first column is obtained by Algorithm 1, the second column by Algorithm 2, and the third column by Algorithm 3. As can be seen, after sensitivity matrix correction (2nd column and 3rd column), the artifacts of the reconstructed images are reduced and the boundary between two phases is clear. Comparing with Algorithm 2 (the 2nd column in Figure 7), Algorithm 3 reduces the artifacts of reconstruction images, where the contour of the discrete phase shrinks obviously and is close to the actual model distribution, as shown in Figure 6.Quantitative analyses of the reconstructed images are showed in Figures8 and 9. Note that the sensitivity matrix is not corrected in Algorithm 1, and the number of corrections represents the number of iterations corresponding to the other two algorithms.Figure 8
Reconstruction errors of three algorithms without noise. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)Figure 9
Correlation coefficients of three algorithms without noise. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)From Figures8 and 9, it can be seen that the quality of the images obtained by Algorithm 3 are obviously improved. For models D2 and D3, the differences between the model boundary voltage and the uniform field boundary voltage are relatively small (see Figure 10(a)) and the initial correction coefficient is larger (see Figure 10(b)). After the first time distribution correction, the distributions reconstructed by Algorithm 3 have larger image errors (as shown in Figures 8(b) and 8(c)) and lower correlations (as shown in Figures 9(b) and 9(c)), but the quality of reconstructed images obtained after the subsequent correction are improved, compared with Algorithms 1 and 2. For the models D1 and D4, the initial changes in voltage are greater (see Figure 10(a)). After punitive restraint on the reconstruction distribution, the reconstruction errors are decreased (as shown in Figures 8(a) and 8(d)) and reconstruction correlation coefficients are obviously increased (as shown in Figures 9(a) and 9(d)).Figure 10
Penalty factors and correction coefficients of four models. (a) Penalty factors/relative changes in boundary voltage. (b) Correction Coefficients.
(a)
(b)Quantitative analysis of matrix substitution error is shown in Figure11. It can be seen that for the models D1, D3, and D4 that the square sum of the matrix substitution error ηp is significantly reduced after the sensitivity matrix is modified by the proposed algorithm. For the model D2, ηp increases after the first correction, the reason for this is the constant term B=1.2, which makes the β>1 (see Figure 10(b)), resulting in over relaxation. However, the error after the second correction is markedly reduced. The error in the subsequent correction process has also been in a downward trend. Compared with Algorithm 2 (βp+1=1) and Algorithm 3, if the initial matrix substitution error η0 is small (see Figures 11(b) and 11(c)), both algorithms can effectively reduce the matrix substitution error ηp, and the effect of the penalty factor is not obvious. However, when the initial error η0 is large (see Figures 11(a) and 11(d)), the reduction of the error ηp after the correction of the sensitivity matrix is more obvious with Algorithm 3. It can also be seen from Figures 8 and 9 that for the models D1 and D4, the decrease in error and the improvement of the correlation coefficient is more significantly with Algorithm 3.Figure 11
Matrix substitution errorηp. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)Because the correction requires additional calculation of the sensitivity matrix, the amount of computation increases significantly. Therefore, the number of matrix corrections must be limited. It can be seen in Figures8 and 9 that after twice corrections, the quality of the reconstructed images has been greatly improved and the sensitivity matrix has only been recalculated once at this time. Therefore, the distribution after correcting twice can be used as the final result. Figure 12 shows data comparison with different algorithms in terms of the reconstruction error and correlation coefficient after corrected twice. It can be seen that errors are obviously reduced and correlations improved with Algorithm 3.Figure 12
(a) Reconstruction errors and (b) image correlation coefficients with different algorithms after 600 iterations.
(a)
(b)From the quantitative analyses, it can be found that the reconstruction errors are significantly reduced with the proposed algorithm and the correlations improved.


#### 4.2.2. Results with Noisy Data

Because noise is inevitable and affects image reconstruction, it is necessary to analyse the antinoise ability of the algorithm. By adding random noise to the boundary voltage vector, the signal-to-noise ratio (SNR) becomes 30 dB. The three algorithms are used to perform 600 iterations for the four models. The obtained images are shown in Figure13. With Algorithms 2 and 3, the sensitivity matrices are corrected once and reconstructed distributions are corrected twice.Figure 13
Reconstructed images after 600 iterations with noise. (a) Reconstructed images of D1. (b) Reconstructed images of D2. (c) Reconstructed images of D3. (d) Reconstructed images of D4.
(a)
(b)
(c)
(d)It can be seen in Figure13 that reconstructed images by Algorithm 3 are less distorted with fewer artifacts than those obtained by the other two algorithms. Quantitative analyses of the reconstructed images with noise data are showed in Figure 14. It can be seen that the errors of reconstructed distributions obtained by Algorithm 3 are lower (see Figure 14(a)) and correlation coefficients are greater (see Figure 14(b)).Figure 14
(a) Reconstruction errors and (b) image correlation coefficients of different algorithms after 600 iterations with noise.
(a)
(b)From the above results, it is concluded that the reconstruction error can be effectively reduced and the correlation can be improved by the proposed algorithm.


#### 4.2.3. Experimental Results

Experiment was conducted to validate the proposed algorithm. The experimental setup is shown in Figure15, which was developed at the University of Leeds. Configurations of the experimental system are listed below.(i)
Injected current: 15.10 mA(ii)
Current Frequency: 9600 Hz(iii)
Diameter of test filed: 50 mm(iv)
ERT sensor: single plane 16-electrodes.Figure 15
Experimental setup of the ERT test system.Six distribution models are selected, as shown in Figure16. The three algorithms are used to reconstruct six models. The reconstructed images after 600 iterations with real data are shown in Figure 17.Figure 16
Conductivity distribution models for experimentation. (a) Model 1. (b) Model 2. (c) Model 3. (d) Model 4. (e) Model 5. (f) Model 6.
(a)
(b)
(c)
(d)
(e)
(f)Figure 17
Reconstructed images after 600 iterations with real data. (a) Reconstructed images of Model 1. (b) Reconstructed images of Model 2. (c) Reconstructed images of Model 3. (d) Reconstructed images of Model 4. (e) Reconstructed images of Model 5. (f) Reconstructed images of Model 6.
(a)
(b)
(c)
(d)
(e)
(f)It can be seen from Figure17 that the proposed algorithm effectively distinguishes targets. For both models, the artifacts with the proposed algorithm are less and the two-phase boundary of the distribution is clear.


## 4.1. Experimental Conditions

MatLab is used for simulation. The configuration of the simulated ERT sensor is listed below:(i)
Circular inner diameter: 50 mm(ii)
Liquid conductivity: 1 S/m(iii)
Gas conductivity: 0 S/m(iv)
Number of electrodes: 16(v)
Excitation and measurement mode: current injection to adjacent electrode pair and voltage measurement from the other adjacent electrode pairs (see Figure1)(vi)
Number of independent measurement: 104(vii)
Number of finite elements: 1536 (see Figure5)To verify the proposed algorithm, three algorithms are used to reconstruct images for the four models, as shown in Figure6. A uniform field is used as the initial distribution. The three algorithms are shown in Table 1.Figure 5
Mesh grid.Figure 6
Four conductivity distribution models for simulation. (a) Model D1 (b) Model D2. (c) Model D3. (d) Model D4.
(a)
(b)
(c)
(d)Table 1
Algorithm table.
Algorithm 1 [10, 16]
Δ
σ
1
=
S
T
Δ
Z
Δ
σ
k
+
1
=
Δ
σ
k
+
α
⋅
S
T
Δ
Z
−
S
⋅
Δ
σ
k
σ
=
Q
σ
0
+
Δ
σ
M
Algorithm 2 [12]
Δ
σ
p
=
Δ
σ
M
S
P
,
Δ
Z
P
Conventional Landweber
σ
p
+
1
=
Q
σ
p
+
β
p
+
1
Δ
σ
p
β
p
+
1
=
1
Algorithm 3
Δ
σ
p
=
Δ
σ
M
S
P
,
Δ
Z
P
Conventional Landweber
σ
p
+
1
=
Q
σ
p
+
β
p
+
1
Δ
σ
p
β
p
+
1
=
B
−
Δ
Z
p
/
Z
pIn Algorithm 1 (uncorrected), the conventional Landweber iterative algorithm [10] is used and the sensitivity matrix is not corrected during the iterative process. In Algorithm 2 [12] (βp+1=1), βp+1=1. M=500, when p=1, and M=100, when p>1. p=0,1,…,P,P=10 and ε=0. In Algorithm 3 (βp+1=1.2−ΔZp/Zp), the algorithm proposed in this paper is used, βp+1=B−ΔZp/Zp,B=1.2. M=500, when p=1, and M=100, when p>1. p=0,1…P,P=10 and ε=0.A reconstruction error and a correlation coefficient are used to evaluate the reconstruction results. The smaller the reconstruction error and the greater the reconstruction correlation means higher image quality [10]:(17)Reconstruction error=σ^−σσ×100%,Correlation coefficient=∑i=1Kσi−σ¯σ^i−σ^¯∑i=1Kσi−σ¯2∑i=1Kσi−σ^¯2,where σ^ is the calculated conductivity distribution, σ is the true conductivity distribution, and K = 1536 is the number of elements in the sensing area.Quantitative analysis is performed by calculating the square sum of the substitution error of the matrix:(18)Matrix substitution errorηp=Sp−S⋅σ.


## 4.2. Results and Analysis

### 4.2.1. Results with Noise Free Data

Without considering noise, the above three algorithms are used for reconstruction and the results are shown in Figure7. With Algorithm 1, the reconstructed images with 500 iterations, 600 iterations, and 1400 iterations are given. For Algorithms 2 and 3, the reconstructed images after 1, 2, and 10 corrections are given, and the images correspond to iterations 500 times, 600 times, and 1400 times. The third column (βp+1=1.2−ΔZp/Zp) is obtained using Algorithm 3.Figure 7
Reconstruction of three algorithms without noise. (a) Reconstructed image of model D1. (b) Reconstructed image of model D2. (c) Reconstructed image of model D3. (d) Reconstructed image of model D4.
(a)
(b)
(c)
(d)In Figure7, the first column is obtained by Algorithm 1, the second column by Algorithm 2, and the third column by Algorithm 3. As can be seen, after sensitivity matrix correction (2nd column and 3rd column), the artifacts of the reconstructed images are reduced and the boundary between two phases is clear. Comparing with Algorithm 2 (the 2nd column in Figure 7), Algorithm 3 reduces the artifacts of reconstruction images, where the contour of the discrete phase shrinks obviously and is close to the actual model distribution, as shown in Figure 6.Quantitative analyses of the reconstructed images are showed in Figures8 and 9. Note that the sensitivity matrix is not corrected in Algorithm 1, and the number of corrections represents the number of iterations corresponding to the other two algorithms.Figure 8
Reconstruction errors of three algorithms without noise. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)Figure 9
Correlation coefficients of three algorithms without noise. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)From Figures8 and 9, it can be seen that the quality of the images obtained by Algorithm 3 are obviously improved. For models D2 and D3, the differences between the model boundary voltage and the uniform field boundary voltage are relatively small (see Figure 10(a)) and the initial correction coefficient is larger (see Figure 10(b)). After the first time distribution correction, the distributions reconstructed by Algorithm 3 have larger image errors (as shown in Figures 8(b) and 8(c)) and lower correlations (as shown in Figures 9(b) and 9(c)), but the quality of reconstructed images obtained after the subsequent correction are improved, compared with Algorithms 1 and 2. For the models D1 and D4, the initial changes in voltage are greater (see Figure 10(a)). After punitive restraint on the reconstruction distribution, the reconstruction errors are decreased (as shown in Figures 8(a) and 8(d)) and reconstruction correlation coefficients are obviously increased (as shown in Figures 9(a) and 9(d)).Figure 10
Penalty factors and correction coefficients of four models. (a) Penalty factors/relative changes in boundary voltage. (b) Correction Coefficients.
(a)
(b)Quantitative analysis of matrix substitution error is shown in Figure11. It can be seen that for the models D1, D3, and D4 that the square sum of the matrix substitution error ηp is significantly reduced after the sensitivity matrix is modified by the proposed algorithm. For the model D2, ηp increases after the first correction, the reason for this is the constant term B=1.2, which makes the β>1 (see Figure 10(b)), resulting in over relaxation. However, the error after the second correction is markedly reduced. The error in the subsequent correction process has also been in a downward trend. Compared with Algorithm 2 (βp+1=1) and Algorithm 3, if the initial matrix substitution error η0 is small (see Figures 11(b) and 11(c)), both algorithms can effectively reduce the matrix substitution error ηp, and the effect of the penalty factor is not obvious. However, when the initial error η0 is large (see Figures 11(a) and 11(d)), the reduction of the error ηp after the correction of the sensitivity matrix is more obvious with Algorithm 3. It can also be seen from Figures 8 and 9 that for the models D1 and D4, the decrease in error and the improvement of the correlation coefficient is more significantly with Algorithm 3.Figure 11
Matrix substitution errorηp. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)Because the correction requires additional calculation of the sensitivity matrix, the amount of computation increases significantly. Therefore, the number of matrix corrections must be limited. It can be seen in Figures8 and 9 that after twice corrections, the quality of the reconstructed images has been greatly improved and the sensitivity matrix has only been recalculated once at this time. Therefore, the distribution after correcting twice can be used as the final result. Figure 12 shows data comparison with different algorithms in terms of the reconstruction error and correlation coefficient after corrected twice. It can be seen that errors are obviously reduced and correlations improved with Algorithm 3.Figure 12
(a) Reconstruction errors and (b) image correlation coefficients with different algorithms after 600 iterations.
(a)
(b)From the quantitative analyses, it can be found that the reconstruction errors are significantly reduced with the proposed algorithm and the correlations improved.


### 4.2.2. Results with Noisy Data

Because noise is inevitable and affects image reconstruction, it is necessary to analyse the antinoise ability of the algorithm. By adding random noise to the boundary voltage vector, the signal-to-noise ratio (SNR) becomes 30 dB. The three algorithms are used to perform 600 iterations for the four models. The obtained images are shown in Figure13. With Algorithms 2 and 3, the sensitivity matrices are corrected once and reconstructed distributions are corrected twice.Figure 13
Reconstructed images after 600 iterations with noise. (a) Reconstructed images of D1. (b) Reconstructed images of D2. (c) Reconstructed images of D3. (d) Reconstructed images of D4.
(a)
(b)
(c)
(d)It can be seen in Figure13 that reconstructed images by Algorithm 3 are less distorted with fewer artifacts than those obtained by the other two algorithms. Quantitative analyses of the reconstructed images with noise data are showed in Figure 14. It can be seen that the errors of reconstructed distributions obtained by Algorithm 3 are lower (see Figure 14(a)) and correlation coefficients are greater (see Figure 14(b)).Figure 14
(a) Reconstruction errors and (b) image correlation coefficients of different algorithms after 600 iterations with noise.
(a)
(b)From the above results, it is concluded that the reconstruction error can be effectively reduced and the correlation can be improved by the proposed algorithm.


### 4.2.3. Experimental Results

Experiment was conducted to validate the proposed algorithm. The experimental setup is shown in Figure15, which was developed at the University of Leeds. Configurations of the experimental system are listed below.(i)
Injected current: 15.10 mA(ii)
Current Frequency: 9600 Hz(iii)
Diameter of test filed: 50 mm(iv)
ERT sensor: single plane 16-electrodes.Figure 15
Experimental setup of the ERT test system.Six distribution models are selected, as shown in Figure16. The three algorithms are used to reconstruct six models. The reconstructed images after 600 iterations with real data are shown in Figure 17.Figure 16
Conductivity distribution models for experimentation. (a) Model 1. (b) Model 2. (c) Model 3. (d) Model 4. (e) Model 5. (f) Model 6.
(a)
(b)
(c)
(d)
(e)
(f)Figure 17
Reconstructed images after 600 iterations with real data. (a) Reconstructed images of Model 1. (b) Reconstructed images of Model 2. (c) Reconstructed images of Model 3. (d) Reconstructed images of Model 4. (e) Reconstructed images of Model 5. (f) Reconstructed images of Model 6.
(a)
(b)
(c)
(d)
(e)
(f)It can be seen from Figure17 that the proposed algorithm effectively distinguishes targets. For both models, the artifacts with the proposed algorithm are less and the two-phase boundary of the distribution is clear.


## 4.2.1. Results with Noise Free Data

Without considering noise, the above three algorithms are used for reconstruction and the results are shown in Figure7. With Algorithm 1, the reconstructed images with 500 iterations, 600 iterations, and 1400 iterations are given. For Algorithms 2 and 3, the reconstructed images after 1, 2, and 10 corrections are given, and the images correspond to iterations 500 times, 600 times, and 1400 times. The third column (βp+1=1.2−ΔZp/Zp) is obtained using Algorithm 3.Figure 7
Reconstruction of three algorithms without noise. (a) Reconstructed image of model D1. (b) Reconstructed image of model D2. (c) Reconstructed image of model D3. (d) Reconstructed image of model D4.
(a)
(b)
(c)
(d)In Figure7, the first column is obtained by Algorithm 1, the second column by Algorithm 2, and the third column by Algorithm 3. As can be seen, after sensitivity matrix correction (2nd column and 3rd column), the artifacts of the reconstructed images are reduced and the boundary between two phases is clear. Comparing with Algorithm 2 (the 2nd column in Figure 7), Algorithm 3 reduces the artifacts of reconstruction images, where the contour of the discrete phase shrinks obviously and is close to the actual model distribution, as shown in Figure 6.Quantitative analyses of the reconstructed images are showed in Figures8 and 9. Note that the sensitivity matrix is not corrected in Algorithm 1, and the number of corrections represents the number of iterations corresponding to the other two algorithms.Figure 8
Reconstruction errors of three algorithms without noise. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)Figure 9
Correlation coefficients of three algorithms without noise. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)From Figures8 and 9, it can be seen that the quality of the images obtained by Algorithm 3 are obviously improved. For models D2 and D3, the differences between the model boundary voltage and the uniform field boundary voltage are relatively small (see Figure 10(a)) and the initial correction coefficient is larger (see Figure 10(b)). After the first time distribution correction, the distributions reconstructed by Algorithm 3 have larger image errors (as shown in Figures 8(b) and 8(c)) and lower correlations (as shown in Figures 9(b) and 9(c)), but the quality of reconstructed images obtained after the subsequent correction are improved, compared with Algorithms 1 and 2. For the models D1 and D4, the initial changes in voltage are greater (see Figure 10(a)). After punitive restraint on the reconstruction distribution, the reconstruction errors are decreased (as shown in Figures 8(a) and 8(d)) and reconstruction correlation coefficients are obviously increased (as shown in Figures 9(a) and 9(d)).Figure 10
Penalty factors and correction coefficients of four models. (a) Penalty factors/relative changes in boundary voltage. (b) Correction Coefficients.
(a)
(b)Quantitative analysis of matrix substitution error is shown in Figure11. It can be seen that for the models D1, D3, and D4 that the square sum of the matrix substitution error ηp is significantly reduced after the sensitivity matrix is modified by the proposed algorithm. For the model D2, ηp increases after the first correction, the reason for this is the constant term B=1.2, which makes the β>1 (see Figure 10(b)), resulting in over relaxation. However, the error after the second correction is markedly reduced. The error in the subsequent correction process has also been in a downward trend. Compared with Algorithm 2 (βp+1=1) and Algorithm 3, if the initial matrix substitution error η0 is small (see Figures 11(b) and 11(c)), both algorithms can effectively reduce the matrix substitution error ηp, and the effect of the penalty factor is not obvious. However, when the initial error η0 is large (see Figures 11(a) and 11(d)), the reduction of the error ηp after the correction of the sensitivity matrix is more obvious with Algorithm 3. It can also be seen from Figures 8 and 9 that for the models D1 and D4, the decrease in error and the improvement of the correlation coefficient is more significantly with Algorithm 3.Figure 11
Matrix substitution errorηp. (a) D1. (b) D2. (c) D3. (d) D4.
(a)
(b)
(c)
(d)Because the correction requires additional calculation of the sensitivity matrix, the amount of computation increases significantly. Therefore, the number of matrix corrections must be limited. It can be seen in Figures8 and 9 that after twice corrections, the quality of the reconstructed images has been greatly improved and the sensitivity matrix has only been recalculated once at this time. Therefore, the distribution after correcting twice can be used as the final result. Figure 12 shows data comparison with different algorithms in terms of the reconstruction error and correlation coefficient after corrected twice. It can be seen that errors are obviously reduced and correlations improved with Algorithm 3.Figure 12
(a) Reconstruction errors and (b) image correlation coefficients with different algorithms after 600 iterations.
(a)
(b)From the quantitative analyses, it can be found that the reconstruction errors are significantly reduced with the proposed algorithm and the correlations improved.


## 4.2.2. Results with Noisy Data

Because noise is inevitable and affects image reconstruction, it is necessary to analyse the antinoise ability of the algorithm. By adding random noise to the boundary voltage vector, the signal-to-noise ratio (SNR) becomes 30 dB. The three algorithms are used to perform 600 iterations for the four models. The obtained images are shown in Figure13. With Algorithms 2 and 3, the sensitivity matrices are corrected once and reconstructed distributions are corrected twice.Figure 13
Reconstructed images after 600 iterations with noise. (a) Reconstructed images of D1. (b) Reconstructed images of D2. (c) Reconstructed images of D3. (d) Reconstructed images of D4.
(a)
(b)
(c)
(d)It can be seen in Figure13 that reconstructed images by Algorithm 3 are less distorted with fewer artifacts than those obtained by the other two algorithms. Quantitative analyses of the reconstructed images with noise data are showed in Figure 14. It can be seen that the errors of reconstructed distributions obtained by Algorithm 3 are lower (see Figure 14(a)) and correlation coefficients are greater (see Figure 14(b)).Figure 14
(a) Reconstruction errors and (b) image correlation coefficients of different algorithms after 600 iterations with noise.
(a)
(b)From the above results, it is concluded that the reconstruction error can be effectively reduced and the correlation can be improved by the proposed algorithm.


## 4.2.3. Experimental Results

Experiment was conducted to validate the proposed algorithm. The experimental setup is shown in Figure15, which was developed at the University of Leeds. Configurations of the experimental system are listed below.(i)
Injected current: 15.10 mA(ii)
Current Frequency: 9600 Hz(iii)
Diameter of test filed: 50 mm(iv)
ERT sensor: single plane 16-electrodes.Figure 15
Experimental setup of the ERT test system.Six distribution models are selected, as shown in Figure16. The three algorithms are used to reconstruct six models. The reconstructed images after 600 iterations with real data are shown in Figure 17.Figure 16
Conductivity distribution models for experimentation. (a) Model 1. (b) Model 2. (c) Model 3. (d) Model 4. (e) Model 5. (f) Model 6.
(a)
(b)
(c)
(d)
(e)
(f)Figure 17
Reconstructed images after 600 iterations with real data. (a) Reconstructed images of Model 1. (b) Reconstructed images of Model 2. (c) Reconstructed images of Model 3. (d) Reconstructed images of Model 4. (e) Reconstructed images of Model 5. (f) Reconstructed images of Model 6.
(a)
(b)
(c)
(d)
(e)
(f)It can be seen from Figure17 that the proposed algorithm effectively distinguishes targets. For both models, the artifacts with the proposed algorithm are less and the two-phase boundary of the distribution is clear.


## 5. Conclusions

In this paper, an image reconstruction algorithm for ERT is proposed. It is based on correction of the sensitivity matrix and Landweber iteration. The error introduced by the initial sensitivity matrix is estimated using the relative variation in the boundary voltage. According to the estimated value, the correction coefficient with a penalty factor is obtained and the distribution reconstructed from Landweber iteration and the initial sensitivity matrix are corrected. To improve the reconstruction accuracy, the distribution and sensitivity matrix are alternatively corrected.Both simulation and experimental results show that the boundary between two phases of the distribution can be reconstructed by the proposed algorithm clearly. The artifacts are less and the distribution is more accurate than the traditional Landweber algorithm and the algorithm without the penalty factor. Quantitative analysis shows that the proposed algorithm is effective in reducing errors in particular when the initial error is large.The correlation between the variation in the boundary voltage and the matrix substitution error is analyzed. In the future, research will be carried out on correlation to further reduce the image error and improve the reconstruction accuracy.



---

*Source: 6384132-2019-12-28.xml*