Differential Electric Field Tomography

Abstract—The article considers electric field tomography, which has some unique properties. A solution is proposed to the problem of detecting the useful signal when a studied object interacts with an electromagnetic field, which is due to the high measurement error associated with variation in the capacitances of the electrodes and the object with respect to the ground. Numerical simulation is used to demonstrate the feasibility of differential excitations of measurements, which can significantly reduce the error and increase the noise immunity of the system. It is shown that a system of eight electrodes is sufficient for reconstructing the image by convolution and back projection. In this case, phase measurements can be used, similarly to unipolar electric field tomography, and, in some cases, amplitude measurements.

INTRODUCTION

Electric field tomography (EFT) is a type of quasistatic electromagnetic tomography that allows contactless visualization of the spatial distribution of an object’s electrical properties. It has several advantages over other types of quasistatic electromagnetic tomography. In addition to electrical capacitance and magnetic induction tomography, it does not require contact with the object of study, unlike electrical impedance tomography [1], which eliminates problems of accurate electrode positioning and spurious contact phenomena. In addition, the electrical impedance method cannot yield a tomogram of the skull, since weakly conducting bone tissue prevents penetration of electric current. Magnetic induction tomography (MIT) copes with a similar problem [2, 3]; however, the signal-to-noise ratio suffers from interference from the electric field. Despite the constructive similarity of EFT and electric capacitance tomography [4], the latter is unsuitable for visualizing media with an electrical conductivity greater than 0.1 S/m and, therefore, is of no interest for biomedicine.

EFT involves probing an object with an external electromagnetic field. When the object and the electric component of the electromagnetic field interact, changes in the phase and amplitude of the latter can be recorded. The authors of [5, 6] presented a variant that measures the phase shift of the probe field due to the Maxwell–Wagner relaxation phenomenon, the essence of which is the delay in the redistribution of free carriers under the action of a probe field due to the finite conductivity of the studied object.

Experimental implementation of the described method was demonstrated in [7, 8]. In a multielectrode system, a unipolar method of exciting and recording the electric field parameters was used, in which all voltages were measured with respect to a common wire—a grounded metal screen. This method has a certain drawback—a significant error in measuring the useful signal due to variations in the capacitances of the electrodes and the object itself with respect to the ground.

The possibility of differential excitations and measurements for leveling the described errors was demonstrated in [9]. In a differential circuit, the voltage source and recorder are each connected symmetrically to two electrodes of the measuring system. In this case, not only phase shifts, but also relative changes in the electric field amplitude can be used for image reconstruction, which has been confirmed by numerical simulation for a system of four electrodes. The sensitivity zone for differential measurements can be approximated by equipotential lines of the probe field. This allows the back projection method to be used for image reconstruction.

A full demonstration of the feasibility of a differential circuit for switching on the electrodes requires simulation of a multielectrode tomographic system. This is the main objective of the study.

1 MATERIALS AND METHODS

Historically, the finite element method (FEM) has occupied the dominant position when numerically simulating quasistatic tomographic systems. Due to the relatively high operating frequency in MIT and EFT, it is preferable to use the finite difference time domain (FDTD) method, based on finite-difference approximation of Maxwell’s differential equations. FDTD makes it possible to take into account the wave propagation of the field, in contrast to the FEM method. In addition, the latter may be numerically unstable, in particular, when simulating the boundaries of a conducting and nonconducting medium. FDTD has no such drawback and ensures numerically stability and reliability of the simulation results.

The initial FDTD equations are Maxwell’s equations (in differential form):

$$\nabla \times \vec {E} = - ~\frac{{\partial \vec {B}}}{{\partial t}},\,\,\,\,~\nabla \times \vec {H} = \,\,\,~\frac{{\partial \vec {D}}}{{\partial t}} + \vec {J},~$$

(1)

where $\vec {E}$ is the electric field strength, $\vec {H}$ is the magnetic field strength, $\vec {D}$ is the electric field induction, $\vec {B}$ is the magnetic field induction, and $\vec {J}$ is the current density vector. If the medium is linear and has no dispersion, then $\vec {D}$ = $\varepsilon {{\varepsilon }_{0}}\vec {E}$, $\vec {B}$ = $\mu {{\mu }_{0}}\vec {H}$, $\vec {J}$ = σ$\vec {E}$, where ε and μ are the dielectric constant and magnetic permeability, ε₀ and µ₀ are the electric and magnetic constants, and σ is the electrical conductivity of the medium. The remaining two Maxwell equations for divergences are automatically satisfied if the correct boundary and initial conditions are imposed. The simplest boundary conditions correspond to a perfect electrical conductor (PEC): the tangential electric field at the boundary is zero; i.e., electromagnetic waves are 100% returned to the computed volume. Zero fields are the typical initial conditions. In FDTD, the derivatives in Eqs. (1) are approximated by second-order finite differences using a technique developed by K. Yi [10]. It should be mentioned that explicit finite-difference schemes require special conditions for stable operation. For FDTD, this condition has the form

$$\Delta t \leqslant {1 \mathord{\left/ {\vphantom {1 c}} \right. \kern-0em} c}\sqrt {{{{(\Delta x)}}^{{ - 2}}} + {{{\left( {\Delta y} \right)}}^{{ - 2}}} + {{{\left( {\Delta z} \right)}}^{{ - 2}}}} ~,$$

(2)

where c is the speed of light in vacuum; ∆x, ∆y, ∆z are the parameters of the calculated volume. Ensuring condition (2) is a significant problem when using the method in simulating quasistatic systems, because its implementation requires an extremely large number of time steps if a high spatial resolution is required, and the finite simulation time should overlap the slowest relaxation processes in the system.

The simulation was performed in the program¹ FDTDpro. The measuring system is a hexagonal chamber made of a conductor (conductivity 10⁶ S/m), in which eight (No. 1–No. 8) electrodes are placed with a rectangular cross section of the same material. The system uses dipole electric field excitation, for which the exciting signal is supplied to a pair of adjacent electrodes and the useful signal is recorded by all remaining free pairs of adjacent electrodes. The dimensions of the electrodes, as well as the dimensions of the gap between the wall and the object and the gaps between the electrodes, were chosen so as to minimize cross capacitance coupling. Two antiphase voltage sources are connected to the exciting electrodes; the recorders of potential difference of the electric field are attached to the rest. The studied object was put in the space between the electrodes. The geometry of the measuring system is shown in Fig. 1.

Fig. 1.

Geometry of measuring system.

The system was excited by a Gaussian voltage pulse, the shape and spectrum of which are shown in Fig. 2. Such an excitation signal is convenient in that in just one run of the program, results in a wide frequency range are obtained by Fourier transform during data processing. At frequencies near the relaxation frequency of the chosen medium (approximately 21 MHz), the pulse has a sufficiently large spectral amplitude.

Fig. 2.

Shape of exciting voltage pulse (a) and its spectral modulus (b).

For image reconstruction, this measurement method can use a filtering algorithm and weighted back projection along the lines of maximum sensitivity [5]. Let there be a reference dataset corresponding to measurements taken at a preceding time instant. Then the conductivity value S (in arbitrary units) assigned to the point of the reconstructed image with polar coordinates r and θ can be calculated by the formula describing the back projection procedure along the equipotentials of the electric field:

$$S\left( {r,\theta } \right) = \sum\limits_{i = 0}^{N - 1} W \left( {r,\theta - \frac{{2\pi }}{N}i} \right){{\left. {{{{\tilde {\lambda }}}_{i}}\left( x \right)} \right|}_{{\Phi \left( x \right){\kern 1pt} = {\kern 1pt} \varphi \left( {r,\theta } \right)}}}.$$

Here i is the number of a pair of active electrodes (profile number), ${{W}_{i}}$ is the geometric weight factor ensuring independence of the sensitivity of the tomograph from the coordinates of the points inside the investigated space, ${{\tilde {\lambda }}_{i}}$ is the result of linear interpolation of the discrete argument function ${{\tilde {\lambda }}_{i}} = {{u_{r}^{i}\left( j \right)} \mathord{\left/ {\vphantom {{u_{r}^{i}\left( j \right)} {u_{m}^{i}\left( j \right)}}} \right. \kern-0em} {u_{m}^{i}\left( j \right)}} - 1$, where $u_{m}^{i}\left( j \right)$ is the voltage measured by the jth pair of electrodes (j = 2, …, N – 2, counting takes place from the active pair in this measurement, which corresponds to j = 0), $u_{r}^{i}\left( j \right)$ is the corresponding voltage from the reference dataset; the value of argument x (0 < x <N) is determined by the intersection of the equipotential of the electric field passing through the point r, θ, with the external boundary of the object, and ${{\Phi }_{i}}\left( x \right)$ is the distribution of the potential along the boundary of an object with uniform conductivity. The weight factor ${{W}_{i}}$ is determined by the approximate analytical solution to the direct problem for a given tomograph configuration assuming a uniform conductivity distribution:

$$\begin{gathered} W\left( {r,\theta } \right) \\ = \frac{{8r_{1}^{4}r_{2}^{4}\sqrt {r_{3}^{2} + 3{{{\sin }}^{2}}{{\theta }_{x}}} }}{{r_{3}^{4}\sqrt {\left( {r_{1}^{2} + 3{{r}^{2}}{{{\sin }}^{2}}\theta } \right)\left( {r_{2}^{2} + 3{{r}^{2}}{{{\sin }}^{2}}\left( {\theta - {{\theta }_{x}}} \right)} \right)} }}. \\ \end{gathered} $$

It can also be determined experimentally.

2 DISCUSSION

In the numerical simulation, the dependence of the phase shift introduced by the object on the signal frequency was studied when the object was located at the center of the system. For this, the recorded signals (reference, without an embedded object, and with an embedded object) were subjected to Fourier transform and the obtained images were subtracted. Figure 3 shows the dependence of the phase shift on pair of electrodes nos. 6 and 7. The maximum phase shift is near the relaxation frequency of the medium.

Fig. 3.

Frequency dependence of phase shift introduced by object.

First of all, a model with an object at the center of the system was used to reconstruct the image. Figure 4 shows an image of the object based on the measured relative amplitudes and phases of the useful signal at a frequency of 20.5 MHz, where numerals 1 and 8 show the location of electrode nos. 1 and 8.

Fig. 4.

Image of object at center of system: (a) based on signal amplitude values, (b) based on phase shift values.

Then, before proceeding to visualization of inhomogeneous media, it was decided to obtain an image of a compact object shifted with respect to the center. Figure 5a shows the model itself, and Figs. 5b and 5c, the images reconstructed from the signal amplitudes and phases, respectively.

Fig. 5.

Image of object shifted with respect to center: (a) measurement model, (b) image reconstructed from relative amplitude values, (c) image reconstructed from phase shift values.

To obtain a tomographic image of an extended object with an inhomogeneous electrical conductivity distribution, we used the model shown in Fig. 6.

Fig. 6.

Model with conductive inclusion object against less conductive cylinder.

An object with parameters σ = 0.1 S/m and ε = 81 is placed inside a weakly conducting cylinder with parameters σ = 0.001 S/m and ε = 81 with an offset with respect to the center of about 1/3 the chamber radius. As reference data, we used a model with a weakly conducting cylinder without an inclusion object, which in electrical impedance tomography terminology is called dynamic visualization. The reconstruction result is shown in Fig. 7. If the measurement results in an empty system (static visualization) were used as reference data, then a satisfactory result is obtained only for phase measurements (Fig. 8). Figure 8 shows that the inclusion object is shifted toward the center; therefore, it is necessary to refine the reconstruction algorithm and, when conducting a physical experiment, increase the number of electrodes.

Fig. 7.

Dynamic tomography: (a) based on signal amplitude values, (b) based on phase shift values.

Fig. 8.

Static tomography: (a) based on signal amplitude values, (b) based on phase shift values.

CONCLUSIONS

By simulating a differential measuring system of eight electrodes, it was possible to obtain a tomographic image of a conductive inclusion object against a less conductive cylinder. In this case, only phase measurements are suitable for static tomography. The perturbations created by the object are large for amplitude measurements. Therefore, the linear model is not applicable and the internal structure is lost.

The next step in the study should be a physical experiment with a setup for electric field tomography. The total number of electrodes should be expanded to 16; this will increase the amount of data for reconstruction and, as a result, increase the image resolution.