A Lagrange-Newton Method for EIT/UT Dual-Modality Image Reconstruction

An image reconstruction method is proposed based on Lagrange-Newton method for electrical impedance tomography (EIT) and ultrasound tomography (UT) dual-modality imaging. Since the change in conductivity distribution is usually accompanied with the change in acoustic impedance distribution, the reconstruction targets of EIT and UT are unified to the conductivity difference using the same mesh model. Some background medium distribution information obtained from ultrasound transmission and reflection measurements can be used to construct a hard constraint about the conductivity difference distribution. Then, the EIT/UT dual-modality inverse problem is constructed by an equality constraint equation, and the Lagrange multiplier method combining Newton-Raphson iteration is used to solve the EIT/UT dual-modality inverse problem. The numerical and experimental results show that the proposed dual-modality image reconstruction method has a better performance than the single-modality EIT method and is more robust to the measurement noise.


Introduction
Electrical impedance tomography (EIT) is a non-invasive technique to reconstruct the distribution of the media within a closed vessel based on the differences of their conductivities. Due to its advantages of high speed, non-radiation, and low cost, EIT has been widely researched and used in many fields, such as, industrial process imaging [1,2], medical imaging [3,4], geological exploration [5,6], structural health monitoring [7,8], and nondestructive evaluation [9,10]. However, as a soft-field imaging technique, EIT usually suffers from the low spatial resolution and imaging accuracy for its inherent nonlinearity and ill-posedness [11]. It has hindered the further application of EIT in many fields, especially for medical imaging, where a higher imaging accuracy is usually expected [12]. Thus, the study of improving the EIT spatial resolution and imaging accuracy can always attract the attention of the researchers.
Many researchers have devoted to the study of improving the spatial resolution and imaging accuracy of EIT and many methods have been developed. For the nonlinearity of EIT, the directly linearization [13], iterative linearization [14], and directly nonlinear methods [15] have been proposed to formulate the EIT inverse problem. For the ill-posedness of EIT, the regularization methods, such as, detected from URT are used to improve the EIT inclusion boundary reconstruction accuracy [34]. In addition, the emerging impediography can improve the EIT image accuracy by coupling the mechanical effects driven from ultrasound waves [35][36][37]. Although the EIT/UT dual-modality imaging has made some progresses, the way of information fusion between EIT and UT is still an open issue in EIT/UT dual-modality imaging researches.
In this paper, a dual-modality image reconstruction method combining EIT and UT is proposed. Some prior information about the background medium distribution is determined by the ultrasound reflection and transmission measurements, and then used as a hard constraint added into the EIT inverse problem through an equality constraint equation. Then, the ultrasound constraint EIT image reconstruction problem is reformed by the Lagrange multiplier method and iteratively solved with Newton-Raphson method. To test the performance of the proposed dual-modality imaging method, the numerical and experimental tests are conducted.

EIT Mathematical Model
EIT mathematical model is derived from Maxwell's electromagnetic theory under the condition of low frequencies. In the quasi-static approximation, the electrical fields can be described in terms of a scalar voltage potential u satisfying the following equation: where σ represents the conductivity in the imaging domain Ω.
Considering the contact impedance between the background medium and electrodes, the complete electrode model (CEM) is used to define the boundary conditions [38]: where ∂Ω is the boundary of the imaging domain Ω, e l ⊂ ∂Ω represents the location covered by the l-th electrode, L is the total number of the electrodes, → n is the outward normal vector on ∂Ω, z l is the contact impedance between the l-th electrode and background medium, I l and U l are the electric current and electric potential at the l-th electrode, respectively. Subject to the law of the conservation of charge, the electric current and voltage in the electrodes satisfy the following equation: Thus, with the definition of EIT mathematical model and the corresponding boundary conditions, the electric potential, conductivity, electric current, and contact impedance can be formulated as the following nonlinear equation: Considering a small change on the conductivity, the corresponding change on the boundary voltage can be linearly represented using Equation (4) as follows: ∆U= u(σ, z, I) − u(σ 0 , z, I) = du(σ, z, I) dσ If the change of conductivity ∆σ is small enough, the high order terms O (∆σ) 2 can be neglected, and Equation (5) will be simplified as follows: further, the Equation (6) is usually discretized as follows: where y ∈ R M×1 represents the boundary voltage difference, x ∈ R N×1 represents the unknown conductivity difference, J ∈ R M×N represents the Jacobian matrix, which can be calculated using the Geselowitz's sensitivity theorem [39]: where J m,n is the element at the m-th row and n column of the Jacobian matrix, u I i is the electric potential for the i-th electrode pair with electric current I i , u I j is the electric potential for the j-th electrode pair with electric current I j , Ω n is the region covered by n-th pixel. For a given conductivity difference x, the boundary voltage y can be determined by the Equations (1)-(3). This process is usually mentioned as the EIT forward problem, where the Boundary Element Method (BEM) solver developed by Ren et al. is used to calculate the EIT forward problem and obtain the sensitivity matrix by the perturbation method [40,41]. In contrast, given the boundary voltage y, the EIT inverse problem is to calculate the conductivity difference x by the Equation (7). Since the sensitivity matrix J is usually singular, the unknown conductivity difference x cannot be directly obtained from the Equation (7) in EIT inverse problem. The regularization methods are commonly used to solve it, which will be introduced in the latter section.

UT Mathematical Model
When the influence of medium heterogeneity on acoustic wave propagation is not considered, the propagation of ultrasonic wave in medium can be approximated as linear propagation model. This model is usually referred as the geometrical/ray acoustic model and widely used in UT image reconstruction. Thus, the mathematical model of UT can be referred to the geometrical optics or ray optics model, and the most commonly considered modes of the ultrasonic wave propagation are transmission and reflection modes, as shown in Figure 1. In addition, the piston transducer is used in the paper, where the energy of the ultrasonic wave is mainly concentrated in a narrow main lobe area [42]. The power of the transmission and reflection wave is mainly affected by the acoustic impedance of the media where the ultrasonic wave pass through, and it can be described as follows: where P i , P r , P t represent incident, reflection, and transmission sound pressures, respectively. α R and α T represent the reflection and transmission coefficients, respectively. Z 1 and Z 2 are the acoustic impedance of the media and equal to the product of medium density ρ and sound velocity c.

Transmission Mode Mathematical Model
Under the assumption that the ultrasonic wave is propagating along the straight line, the mathematical model of transmission mode UT is similar to the x-ray CT. The ultrasonic waves sending from the exciting transducer will pass through the media and arrive at the measuring transducers. Consider that the acoustic attenuation effect of the medium, the transmission mode UT usually utilizes the energy attenuation of the ultrasonic waves to invert the medium distribution in the ultrasonic wave propagation paths. Inspired from the Equation (10), we reformulate the sending acoustic energy e I from the exciting transducer and the receiving acoustic energy m I from the measuring transducer as the following equation: where μ is the medium's attenuation coefficient, Ω k is the region where the k-th ultrasonic wave pass through. Consider a scenario that medium located in the region Ω k changes, the corresponding attenuation coefficient μ in Ω k will be also changing to μ δμ + . Then, the energy attenuation of the ultrasonic wave in the region Ω k will change and the Equation (11) can be rewritten as follows: where δμ is the variation of the attenuation coefficient in Ω k , b m I and a m I represent the acoustic energy before and after the attenuation coefficient changes, respectively.

Reflection Mode Mathematical Model
Similar to the propagation of optics ray in the medium, the mathematical model of the reflection mode UT can be derived from the Fermat's principle, where the ultrasonic waves will reflect at the media interface and continue to propagate to the measuring transducer along the shortest propagation path. Thus, the ultrasonic waves sending from the exciting transducer will be reflected at the interface where the acoustic impedances change and continue to propagate in the background medium to the measuring transducer. When the exciting and measuring transducer are determined, the detected points will fall on an elliptical arc. This elliptical arc can be determined by the coordinates of the exciting and measuring transducer and the time of flight (TOF) of the ultrasonic wave, which can be described as follows:

Transmission Mode Mathematical Model
Under the assumption that the ultrasonic wave is propagating along the straight line, the mathematical model of transmission mode UT is similar to the x-ray CT. The ultrasonic waves sending from the exciting transducer will pass through the media and arrive at the measuring transducers. Consider that the acoustic attenuation effect of the medium, the transmission mode UT usually utilizes the energy attenuation of the ultrasonic waves to invert the medium distribution in the ultrasonic wave propagation paths. Inspired from the Equation (10), we reformulate the sending acoustic energy I e from the exciting transducer and the receiving acoustic energy I m from the measuring transducer as the following equation: where µ is the medium's attenuation coefficient, Ω k is the region where the k-th ultrasonic wave pass through. Consider a scenario that medium located in the region Ω k changes, the corresponding attenuation coefficient µ in Ω k will be also changing to µ + δµ. Then, the energy attenuation of the ultrasonic wave in the region Ω k will change and the Equation (11) can be rewritten as follows: where δµ is the variation of the attenuation coefficient in Ω k , I b m and I a m represent the acoustic energy before and after the attenuation coefficient changes, respectively.

Reflection Mode Mathematical Model
Similar to the propagation of optics ray in the medium, the mathematical model of the reflection mode UT can be derived from the Fermat's principle, where the ultrasonic waves will reflect at the media interface and continue to propagate to the measuring transducer along the shortest propagation path. Thus, the ultrasonic waves sending from the exciting transducer will be reflected at the interface where the acoustic impedances change and continue to propagate in the background medium to the measuring transducer. When the exciting and measuring transducer are determined, the detected points will fall on an elliptical arc. This elliptical arc can be determined by the coordinates of the exciting and measuring transducer and the time of flight (TOF) of the ultrasonic wave, which can be described as follows: where c is the sound speed in the background medium, t f is the TOF of the ultrasonic wave in the background medium, (x d , y d ) is the coordinate of the detected point, which obeys the arc distribution, (x e , y e ) and (x m , y m ) are the coordinates of the exciting and measuring transducer, respectively.

Complementarity of EIT and URT Sensitive Field
The boundary measurements of EIT and UT have different sensitivities about the changes of the sensitive parameters in different spaces of the imaging domain, and higher sensitivity means higher detection resolution and accuracy. In EIT, the electric current is injected into the imaging domain through the electrodes attached at the boundary. The electric field lines near the edge of the imaging domain are dense than the central region of the imaging domain. Therefore, the EIT measurements have a higher sensitivity to the changes in conductivity near the edge of the imaging domain. According to the reciprocity principle, the sensitivity of EIT measurements about the changes of conductivity is proportional to the integration of the electric field intensity [43]. As shown in Figure 2a, EIT has higher sensitivity near the edge of the imaging domain. In UT, the ultrasonic waves emitted by the transducers are mainly concentrated in a narrow main lobe range along the radial of the transducer. Only the changes of acoustic impedance occurring at the propagation path of ultrasonic waves will cause the change on reflection or transmission measurements. In the imaging domain, all the ultrasonic waves propagation paths of the transducers are intersected in the central of the imaging domain. Therefore, the UT measurements have a higher sensitivity about the changes in acoustic impedance near the central of the imaging domain [44]. According to the radon transform theory, the sensitivity of UT measurements about the change of acoustic impedance can be formed by the superimposition of all single transducer ultrasonic waves propagation paths. As shown in Figure 2b, UT has higher sensitivity near the central of the imaging domain. Based on the complementarity of EIT and UT sensitivity space distribution, the combination of EIT and UT is expected to improve the imaging resolution and accuracy. x y are the coordinates of the exciting and measuring transducer, respectively.

Complementarity of EIT and URT Sensitive Field
The boundary measurements of EIT and UT have different sensitivities about the changes of the sensitive parameters in different spaces of the imaging domain, and higher sensitivity means higher detection resolution and accuracy. In EIT, the electric current is injected into the imaging domain through the electrodes attached at the boundary. The electric field lines near the edge of the imaging domain are dense than the central region of the imaging domain. Therefore, the EIT measurements have a higher sensitivity to the changes in conductivity near the edge of the imaging domain. According to the reciprocity principle, the sensitivity of EIT measurements about the changes of conductivity is proportional to the integration of the electric field intensity [43]. As shown in Figure 2a, EIT has higher sensitivity near the edge of the imaging domain. In UT, the ultrasonic waves emitted by the transducers are mainly concentrated in a narrow main lobe range along the radial of the transducer. Only the changes of acoustic impedance occurring at the propagation path of ultrasonic waves will cause the change on reflection or transmission measurements. In the imaging domain, all the ultrasonic waves propagation paths of the transducers are intersected in the central of the imaging domain. Therefore, the UT measurements have a higher sensitivity about the changes in acoustic impedance near the central of the imaging domain [44]. According to the radon transform theory, the sensitivity of UT measurements about the change of acoustic impedance can be formed by the superimposition of all single transducer ultrasonic waves propagation paths. As shown in Figure 2b, UT has higher sensitivity near the central of the imaging domain. Based on the complementarity of EIT and UT sensitivity space distribution, the combination of EIT and UT is expected to improve the imaging resolution and accuracy.

Unification of Mesh Model and Reconstruction Target
In the paper, the purpose of both EIT and UT is to reconstruct the media distribution based on the change of conductivity or acoustic impedance distribution. In some applications, when the media distribution changes, the conductivity and acoustic impedance will also change. For example, in gas/liquid two-phase flow cross-section reconstruction problem, the presence of bubble will change the distribution of conductivity and acoustic impedance, the targets of EIT and UT are all the reconstruction of the bubble distribution. Thus, the reconstruction target of EIT and UT can be

Unification of Mesh Model and Reconstruction Target
In the paper, the purpose of both EIT and UT is to reconstruct the media distribution based on the change of conductivity or acoustic impedance distribution. In some applications, when the media distribution changes, the conductivity and acoustic impedance will also change. For example, in gas/liquid two-phase flow cross-section reconstruction problem, the presence of bubble will change Sensors 2019, 19,1966 7 of 18 the distribution of conductivity and acoustic impedance, the targets of EIT and UT are all the reconstruction of the bubble distribution. Thus, the reconstruction target of EIT and UT can be uniformed to the changes of the media distribution, which means that the background medium has the same representation in EIT and UT. In addition, another factor to consider is the mesh model. Based on this demand of reconstructing the same unknown media distribution using EIT and UT, the mesh model for EIT and UT should be also uniformed in the same imaging region. To achieve the purpose, the same mesh model is used in EIT and UT, as shown in Figure 3. uniformed to the changes of the media distribution, which means that the background medium has the same representation in EIT and UT. In addition, another factor to consider is the mesh model. Based on this demand of reconstructing the same unknown media distribution using EIT and UT, the mesh model for EIT and UT should be also uniformed in the same imaging region. To achieve the purpose, the same mesh model is used in EIT and UT, as shown in Figure 3.

Total Variation Regularization for EIT
It's known that EIT inverse problem is a typical ill-posed inverse problem due to the large condition number of the Jacobian matrix J, which means that it's unrealistic to directly solve the unknown x from the Equation (7). To overcome this problem, the regularization methods are usually adopted. First, we convert the EIT inverse problem to find the minimum from the least-square equation 2 2 Jx y − . Then, a penalty term, also referred as regularization term, is added to stabilize the solution process by constraining the solution space. This process can be formulated in the following equations: where the regularization parameter β is a positive constant and balances the residual term 2 2 Jx y − and the penalty term ( ) R x in the inverse problem solving process. There are many forms of the regularization term based on the difference of the solution space constraints. Here, the Total Variation (TV) regularization term is adopted, and the regularization term can be written as: However, the Equations (14)-(16) cannot be directly solved due to the non-differentiability of x ∇ . To overcome this problem, a smooth approximation is usually adopted [45]. Then, the regularization term ( ) R x can be rewritten as: where 0 θ > is a predetermined small constant and is selected as 6 1 10 − × here.

Constraints Information from UT Measurements
As mentioned above, the purpose of EIT and UT is to reconstruct the media distribution based on the change of conductivity or acoustic impedance from the boundary measurements. Since there

Total Variation Regularization for EIT
It's known that EIT inverse problem is a typical ill-posed inverse problem due to the large condition number of the Jacobian matrix J, which means that it's unrealistic to directly solve the unknown x from the Equation (7). To overcome this problem, the regularization methods are usually adopted. First, we convert the EIT inverse problem to find the minimum from the least-square equation Jx − y 2 2 . Then, a penalty term, also referred as regularization term, is added to stabilize the solution process by constraining the solution space. This process can be formulated in the following equations: where the regularization parameter β is a positive constant and balances the residual term Jx − y 2 2 and the penalty term R(x) in the inverse problem solving process. There are many forms of the regularization term based on the difference of the solution space constraints. Here, the Total Variation (TV) regularization term is adopted, and the regularization term can be written as: However, the Equations (14)-(16) cannot be directly solved due to the non-differentiability of |∇x|. To overcome this problem, a smooth approximation is usually adopted [45]. Then, the regularization term R(x) can be rewritten as: where θ > 0 is a predetermined small constant and is selected as 1 × 10 −6 here.

Constraints Information from UT Measurements
As mentioned above, the purpose of EIT and UT is to reconstruct the media distribution based on the change of conductivity or acoustic impedance from the boundary measurements. Since there is no change in the conductivity and acoustic impedance in the background medium region, the image reconstruction of the background medium satisfies the following relationship in EIT and UT: where Ω b represents the background medium region, Ω b EIT and Ω b UT are the background medium region reconstructed from EIT and UT, respectively.
Thus, the background medium region reconstructed from UT can be used as the prior information in EIT image reconstruction. Consider two scenarios whether there is a target on the ultrasonic wave's propagation path, as shown in Figure 4a,b. If there is no target on the ultrasonic wave propagation path, the ultrasonic waves sending from the exciting transducer will continue to propagate to the measuring transducer and the attenuation of the acoustic energy is only caused by the background medium. And then, we can mark the pixels on the transmission ultrasonic waves propagation path as background medium region Ω b UT . In contrast, if there is a target on the ultrasonic wave propagation path, the ultrasonic waves sending from the exciting transducer will reflect on the target interface, and then continue to propagate to the measuring transducer through the background medium. The time of flight t f can be extracted from the reflection waves. And then, using the Equation (13) and the main lobe area of the transducer, we can calculate the region covered by the reflection ultrasonic waves propagation paths. After that, we can also mark the pixels on the reflection ultrasonic waves propagation path as background medium region Ω b UT .
EIT UT medium region reconstructed from EIT and UT, respectively. Thus, the background medium region reconstructed from UT can be used as the prior information in EIT image reconstruction. Consider two scenarios whether there is a target on the ultrasonic wave's propagation path, as shown in Figure 4a,b. If there is no target on the ultrasonic wave propagation path, the ultrasonic waves sending from the exciting transducer will continue to propagate to the measuring transducer and the attenuation of the acoustic energy is only caused by the background medium. And then, we can mark the pixels on the transmission ultrasonic waves propagation path as background medium region Ω b UT . In contrast, if there is a target on the ultrasonic wave propagation path, the ultrasonic waves sending from the exciting transducer will reflect on the target interface, and then continue to propagate to the measuring transducer through the background medium. The time of flight f t can be extracted from the reflection waves. And then, using the Equation (13) and the main lobe area of the transducer, we can calculate the region covered by the reflection ultrasonic waves propagation paths. After that, we can also mark the pixels on the reflection ultrasonic waves propagation path as background medium region Ω b UT . Since the conductivity of the background medium region does not change, we can derive the following equation: where ( ) h x represents the constraint function of UT measurements with respect to the conductivity difference x, C is the corresponding extend constraint matrix, and it's calculated by the analytical method:   Since the conductivity of the background medium region does not change, we can derive the following equation: where h(x) represents the constraint function of UT measurements with respect to the conductivity difference x, C is the corresponding extend constraint matrix, and it's calculated by the analytical method: where C ∈ R M U ×N , M U is the total number of pixels in region of Ω b UT , n(i) is the serial number of the pixels located at Ω b UT in the set of all pixels.

Lagrange-Newton Solution for Dual-Modality Inverse Problem
In summary, we can convert the media distribution image reconstruction problem into the conductivity change image reconstruction problem for the EIT/UT dual-modality imaging. The EIT/UT dual-modality inverse problem can be described as follows: where f (x) and h(x) are defined by Equations (15) and (19), respectively. Obviously, the problem (21) is a typical equality constraint optimization problem. The Lagrange function combining with Newton method is commonly used to solve it. First, we formulate the following Lagrange function: where λ is the Lagrange multiplier. The stability point of the above Lagrange function corresponds to the optimal solution of the unknown x, and satisfies the following conditions: Then, we solve the problem (22) using the Newton-Raphson method, and the iteration solution can be described as x k+1 = x k + (δx) k , where the iteration step length (δx) k and the Lagrange multiplier (δλ) k satisfy the following equation: where where J T is the transpose of J, L a ∈ R Q×N is a sparse matrix, which can be determined by the topological relationship between the pixels, referring to the literature [45], H ∈ R Q×Q is a diagonal matrix whose diagonal entry can be described as follows: where L a,q ∈ R 1×N is the q-th row element of L a , Q is the total number of adjacent element pairs along the horizontal or vertical direction, referring to the literature [45].
To analyze the iteration convergence and determine the iteration termination condition of the problem (21), we define the following penalty function: In summary, the implement of the proposed EIT/UT dual-modality image reconstruction method can be described as Algorithm 1.

Results
To evaluate the performance of the proposed EIT/UT dual-modality imaging method, a series of numerical and experimental tests are carried out and compared with single-modality EIT method. In the tests, three typical EIT image reconstruction methods, Newton's One-Step Error Reconstructor (NOSER), L1 regularization, and TV regularization, are used in single-modality EIT image reconstruction. Considering that the selection of regularization parameters has a great influence on the imaging results, the regularization parameters in all methods are empirically selected in the range of 10 −6 to 10 −3 , and the one that reached the most similar result to the real model was used in the comparison [45]. In EIT, the adjacent electric current excitation and adjacent voltage measurement model is used to obtain the boundary measurements. In UT, the piston transducers are used in numerical and experimental tests, and the acoustic energy can be better concentrated in a narrow beam range. In transmission mode, three transducers directly opposite to the exciting transducer are used to obtain the transmission measurements. In reflection mode, three transducers on each adjacent side of the exciting transducer are used to obtain the reflection measurements. The ultrasound frequency is 1MHz in the numerical and experimental tests. In addition, the reconstruction results from EIT/UT dual-modality method is simplified as UET in the paper.

Numerical Results and Discussion
The numerical tests are conducted by the commercial multi-physical coupling software COMSOL and MATLAB. The simulation model is designed by simulating the discrete bubbles in a conductive liquid. The conductivity of the bubble is 1 × 10 −6 Sm −1 and the background medium conductivity is 1 × 10 −2 Sm −1 . The density of the bubble is 1.29 kgm −3 and sound velocity is 340 ms −1 . The density of the background medium is 1000 kgm −3 and sound velocity is 1400 ms −1 . The radius of the imaging domain is 75 mm, and 16 electrodes and 16 ultrasound transducers are staggered evenly and arranged on the periphery of the observation domain, as shown in Figure 3. The piston ultrasound transducers sequentially send and receive the three periods sinusoidal pulse signal in the narrow field range.

Quantitative Index
To evaluate the reconstruction results, four quantitative indexes of relative image error (RE), correlation coefficient (CC), position error (PE) and area difference (AD), are used in the paper. The RE and CC are usually used to describe the degree of difference of the true conductivity difference and reconstructed conductivity difference. The definition of RE and CC can be expressed as follows: where x is the calculated conductivity difference vector, x * is the real conductivity difference vector, x n and x * n are the n-th element of x and x * , respectively, x and x * are the mean values of x and x * , respectively. A lower value of RE and a higher value of CC mean a good image reconstruction result.
The PE and AD are respectively used to evaluate the target position and shape reconstructed accuracy by comparing the segmented reconstructed images and the real distribution images. The segmented threshold is selected as 25% of the minimum pixel amplitudes, as follows: where x is the calculated image pixel amplitudes (conductivity difference), x q is the segmented image pixel amplitudes. The PE and AD are defined as follows: where c t is the center of gravity of the target in real conductivity difference distribution image, c q is the center of gravity of the target in reconstructed conductivity difference distribution image, x b is the binary image of the real conductivity difference distribution, x q is the segmented image of the reconstructed conductivity difference distribution image, N is the total number of the pixels.

Image Reconstruction Results Analysis
In numerical tests, five different models are used to compare the proposed EIT/UT dual-modality imaging method (UET) and the single-modality EIT imaging method using NOSER, L1, TV. The image reconstruction results are shown in Figures 5-8 when adding no noise, 40 dB, 25 dB, 20 dB Gaussian noise in EIT measurements, respectively. The real targets boundaries are indicated by the solid black line in the reconstructed images. It can be seen from the Figures 5-8 that the UET has a cleaner background and better imaging performance compared with the other methods, which can be also verified by the quantitative indexes, as shown in Figure 9. For different quantitative indexes, we respectively calculate their mean values in five media distribution models under a certain noise condition. When the noise level increases, the imaging results from single-modality EIT using NOSER, L1 and TV are gradually getting worse, while the proposed UET method can always better reconstruct the media distribution in different noise levels. Thus, the proposed UET method is more robust to the EIT measurement noise. In addition, the expected time of UET is about 0.32 s for per image reconstruction. The expected time of EIT using NOSER, L1, TV are about 0.02 s, 0.2 s, 0.3 s for per image reconstruction, respectively. Thus, the imaging speed of UET is comparable to that of single-modality EIT using TV method, but the imaging accuracy of UET have an obvious improvement. the media distribution in different noise levels. Thus, the proposed UET method is more robust to the EIT measurement noise. In addition, the expected time of UET is about 0.32 s for per image reconstruction. The expected time of EIT using NOSER, L1, TV are about 0.02 s, 0.2 s, 0.3 s for per image reconstruction, respectively. Thus, the imaging speed of UET is comparable to that of single-modality EIT using TV method, but the imaging accuracy of UET have an obvious improvement.   the media distribution in different noise levels. Thus, the proposed UET method is more robust to the EIT measurement noise. In addition, the expected time of UET is about 0.32 s for per image reconstruction. The expected time of EIT using NOSER, L1, TV are about 0.02 s, 0.2 s, 0.3 s for per image reconstruction, respectively. Thus, the imaging speed of UET is comparable to that of single-modality EIT using TV method, but the imaging accuracy of UET have an obvious improvement.

Comparison of Edge Reconstruction Ability
It can be seen from the numerical results that the UET has a cleaner edge than other methods. To verify this conclusion, we further analyze the numerical results of the model with one circle and one square in Figure 5. First, the pixel values of the reconstructed images are normalized to the interval [0, 1]. Then, the contour of the pixel values along the center horizontal line of the reconstructed images is drawn in Figure 10. It's obvious that the contour of UET reconstructed image is closer to the real contour. Thus, the proposed UET method has a better edge reconstruction ability than the single-modality EIT methods.

Comparison of Edge Reconstruction Ability
It can be seen from the numerical results that the UET has a cleaner edge than other methods. To verify this conclusion, we further analyze the numerical results of the model with one circle and one square in Figure 5. First, the pixel values of the reconstructed images are normalized to the interval [0, 1]. Then, the contour of the pixel values along the center horizontal line of the reconstructed images is drawn in Figure 10. It's obvious that the contour of UET reconstructed image is closer to the real contour. Thus, the proposed UET method has a better edge reconstruction ability than the single-modality EIT methods.

Experimental Results and Discussion
To further verify the performance of the proposed method, a set of experimental tests was carried out. The experimental system is shown in Figure 11, which is developed by Tianjin University, China [46]. It mainly contains three parts, the sensor array, the data acquisition and control unit, and the image reconstruction and visualization unit. The radius of the sensor is 75 mm and 16 electrodes and 16 ultrasonic transducers are evenly placed on the periphery of the sensor. The contact resistance between the electrodes and background medium is 37.1 Ω/cm 2 . The background medium consists of Na2SO4 solution with the conductivity of 1.41 mS/cm. The data acquisition speed of EIT is about 1000frames/s, and the data acquisition speed of UT is about 17-frames/s. The expected time of the image reconstruction in experimental tests is similar to numerical tests.
The experimental results are shown in Figure 12, where five different distribution models (A, B, C, D, E) of the nylon rod ( 11 10 − mS/cm) are used. It can be seen from the results that the proposed UET method has a better imaging performance compared with the single-modality EIT methods. The quantitative analysis of the experimental results is shown in Figure 13, which can also verify the conclusion that the proposed UET method shows better reconstruction performance.

Comparison of Edge Reconstruction Ability
It can be seen from the numerical results that the UET has a cleaner edge than other methods. To verify this conclusion, we further analyze the numerical results of the model with one circle and one square in Figure 5. First, the pixel values of the reconstructed images are normalized to the interval [0, 1]. Then, the contour of the pixel values along the center horizontal line of the reconstructed images is drawn in Figure 10. It's obvious that the contour of UET reconstructed image is closer to the real contour. Thus, the proposed UET method has a better edge reconstruction ability than the single-modality EIT methods.

Experimental Results and Discussion
To further verify the performance of the proposed method, a set of experimental tests was carried out. The experimental system is shown in Figure 11, which is developed by Tianjin University, China [46]. It mainly contains three parts, the sensor array, the data acquisition and control unit, and the image reconstruction and visualization unit. The radius of the sensor is 75 mm and 16 electrodes and 16 ultrasonic transducers are evenly placed on the periphery of the sensor. The contact resistance between the electrodes and background medium is 37.1 Ω/cm 2 . The background medium consists of Na2SO4 solution with the conductivity of 1.41 mS/cm. The data acquisition speed of EIT is about 1000frames/s, and the data acquisition speed of UT is about 17-frames/s. The expected time of the image reconstruction in experimental tests is similar to numerical tests.
The experimental results are shown in Figure 12, where five different distribution models (A, B, C, D, E) of the nylon rod ( 11 10 − mS/cm) are used. It can be seen from the results that the proposed UET method has a better imaging performance compared with the single-modality EIT methods. The quantitative analysis of the experimental results is shown in Figure 13, which can also verify the conclusion that the proposed UET method shows better reconstruction performance.

Experimental Results and Discussion
To further verify the performance of the proposed method, a set of experimental tests was carried out. The experimental system is shown in Figure 11, which is developed by Tianjin University, China [46]. It mainly contains three parts, the sensor array, the data acquisition and control unit, and the image reconstruction and visualization unit. The radius of the sensor is 75 mm and 16 electrodes and 16 ultrasonic transducers are evenly placed on the periphery of the sensor. The contact resistance between the electrodes and background medium is 37.1 Ω/cm 2 . The background medium consists of Na 2 SO 4 solution with the conductivity of 1.41 mS/cm. The data acquisition speed of EIT is about 1000-frames/s, and the data acquisition speed of UT is about 17-frames/s. The expected time of the image reconstruction in experimental tests is similar to numerical tests.
The experimental results are shown in Figure 12, where five different distribution models (A, B, C, D, E) of the nylon rod (10 −11 mS/cm) are used. It can be seen from the results that the proposed UET method has a better imaging performance compared with the single-modality EIT methods. The quantitative analysis of the experimental results is shown in Figure 13, which can also verify the conclusion that the proposed UET method shows better reconstruction performance. Sensors 2019, 19, x FOR PEER REVIEW 15 of 18 Figure 11. EIT/UT dual-modality data acquisition system.

Conclusions
A novel EIT/UT dual-modality imaging method is proposed, where the measurement information obtained from UT can be used as prior constraint and improve the EIT image reconstruction performance. The dual-modality inverse problem is constructed by the equality constraint equation, and then, solved by the Lagrange-Newton method. The numerical and experimental tests show that the proposed EIT/UT dual-modality imaging method has a better image reconstruction performance than the traditional single-modality EIT methods, where the reconstructed images of UET have clearer targets edges and less artifacts in the background region. In the experimental tests, the average RE of UET have reduced 29.7%, 31.1%, 14.9% than single-modality EIT using NOSER, L1, TV, respectively. The average CC of UET have improved 7.7%, 16.9%, 7.4% than single-modality EIT using NOSER, L1, TV, respectively. The average PE of UET have reduced 32.9%, 68.2%, 71.3% than single-modality EIT using NOSER, L1, TV, respectively. The average AD of UET have reduced 55.1%, 50.6%, 36.5% than single-modality EIT using NOSER, L1, TV, respectively. In addition, the UET shows better performance of edge reconstruction, where the contour of UET reconstructed image is closer to the real contour than single-modality EIT methods. Furthermore, the proposed UET method is more robust to the EIT measurement noise, where the quantitative indexes of UET have no significant change with the increase of the measurement noise. The future work will focus on the research of enriching the UT measurements by adding the number of the transducers or increasing the diffusion angle of the transducers, and further improve the image reconstruction performance of EIT/UT dual-modality imaging.