Comparison of imaging quality between linear sampling method and time reversal in microwave imaging problems

Abstract

This research investigated imaging quality of two important methods widely used in electromagnetic inverse scattering problems. The algorithms, time reversal (TR) and linear sampling method (LSM), were compared for resolution of point imaging and a correlation indicator to determine image quality. Comparisons were made in single- and multifrequency modes for 2D scenarios in free-space. Comparisons revealed that resolution of TR is much better than LSM. In order to compare the total reconstructed images, several cases were considered to determine a comprehensive conclusion. The simulations were done based on experimental data. In this case, comparisons showed that in the term of correlation indicator, LSM surpasses TR.

Keywords:

• Inverse scattering
• imaging quality
• resolution
• linear sampling method (LSM)
• time reversal (TR)

AMS Subject Classifications:

• 78A25
• 78A46
• 45B05
• 78M05

1. Introduction

Image quality among microwave techniques is open to debate. The term quality is a descriptive expression. It needs to be quantified in order to compare and evaluate different methods. Some prior researches have been done comparing image quality among inverse methods.[1,2] A common technique to quantify this term is considering resolution of method as an index of quality.[3] Unfortunately, resolution is a somewhat ambiguous term and there are several proposed definitions for resolution.[4] The term resolution is an inherent characteristic of each imaging methods and it needs clarification and a mathematical definition. Comparison of imaging resolutions for different imaging methods is commonly applied to determine imaging quality. The intention of this study was to compare resolution between two important inverse methods; linear sampling method (LSM) [5–9] and time reversal (TR).[10–12]

LSM and TR are both methods of qualitative reconstruction. These methods represent different interpretations from a unique physical phenomenon; LSM is a frequency domain method and TR has originated from methods of time domain reconstruction. TR can be considered as an interpretation of distorted wave Born approximation (DWBA).[13] Linear form of inverse scattering theory can be applied in order to implement a matched-filter response in the space-frequency domain, while time domain implementation of TR algorithm with the finite difference time domain (FDTD) method proposes some advantages. Time domain implementations do not need Green function calculations for ambient media, due to inherent characteristics of the FDTD in that inhomogeneity can be inserted directly to FDTD lattice. In addition, FDTD can simulate dispersive media. Furthermore, the FDTD method can simulate full-wave 3D problems. Implementation of TR using FDTD could have practical application in areas such as geophysics, where it is referred as time migration.[14–16]

Several migration techniques have been proposed in recent years. Although these techniques are different but a common characteristic is that most of them can be described using matched-filter concept. The main difference between such methods arises from interpretation of the wave propagation implemented in the filter. A migration method with the focus on reverse time migration is described in [15]. The technique manipulates electromagnetic wave propagation theory and scattering, as well as matched-filtering. The procedure has some suppositions, as described below:

(1)	The profile of the object is divided into several distinct point scatterers. This type of approximation permits evaluation of each single pixel by a separate filter.
(2)	The linear scattering assumption allows each filter to be made independently.
(3)	A priori data or initial information on a background medium is essential for wave propagation (this is a common point of several imaging methods).

Based on these three suppositions, wave migration requires background medium characteristics and weak scatterers. The term ‘weak scatterers’ is ambiguous and there is no standard limit for this definition. This point limits the ability of TR algorithms, however TR methods can be adapted to recover all types of material including perfect electric conductor (PEC) objects.[16] Linearized methods are the first choice for primary reconstruction in distributed radar systems. These groups of methods show very good performance. The main drawback of such methods is failing to incorporate nonlinear nature of multiple-bounce scattering, as such unwanted artefacts are introduced and fine features of the imaging scene could be missed as a result of lowered signal to noise ratio.

Based on described points, LSM reveals advantages in comparison with TR. The LSM is a linear method which is based on solving an integral equation, while needs no linearization. The linear equation in LSM is an exact relation without any approximation. Thus, the LSM can consider the whole aspect of problem including multi-bounces in the data. In the framework defined for LSM, the method acts effective when using a large number of multi-static data. This form could be impractical for several imaging cases. For situations, in which a finite amount of data is available, the quality of LSM image is heavily dependent on the availability of dense multi-static data. In addition, in various practical cases, measurements at high frequencies might be sparse and the distribution of measurement and incident waves could be non-uniform. In these situations, image reconstruction using LSM may fail. It is known in general that the availability of multi-frequency data should compensate for reduced spatial diversity, however it is not directly revealed that how this point could be incorporated for LSM.[17] Thus, one of the major problems for LSM is how to compensate in under-sampled spatial data cases.

This study aimed to compare image quality of the two aforementioned methods in the same framework. The first step was to present a definition of the term resolution. The Rayleigh resolution definition was adopted to make the comparison. Following this definition, resolution was presented for the LSM and TR methods for single- and multifrequency modes.

Resolution describes the image quality for each single pixel. Also important for relating quality of an imaging method is proximity and resemblance between a reconstructed image and the real target; so effectiveness of an algorithm is determined in terms of localization and shape detection of the real target. To the best of our knowledge, a unique agreement in determining this feature of imaging methods has not yet been reached in the literature. An index called image correlation indicator has been applied to compare quality of the reconstructed images between the two methods. This index is defined and used to compare imaging quality accordingly. For the sake of quality comparison between resultant images from these algorithms, reconstructions were made under the Fresnel institute dataset [18] and the different shapes were examined.

2. Mathematical background

A matrix with $N\times N$ dimensions is defined as multistatic data matrix $\mathbf{E}=\{E_{ms}(\mathit{\omega })\}$. Elements of this matrix are determined by the value of the signal received at antenna number $m$, when the excitation is sent as an impulse from the transmitter antenna number $s$. $\mathit{\omega }$ is the signal frequency and $N$ is the number of antennas. In transverse magnetic (TM) mode, and for a two dimensional configuration, it is assumed that the structure is extended across z-axis. Based on this assumption, the vector form of scattering equation could be simplified to scalar equations. Therefore, $E_z$ is used instead of $\overrightarrow{E}$. Figure 1 shows the geometry of the problem in hand with multistatic configuration of filamentary transmitter/receiver antennas which are located on circumference ($\mathrm{\Gamma }$) of the measurement curve. Investigation domain is also denoted with $D$. Assuming time harmonic radiation field ${\mathcal{E}}_z(r,t)=E_z(r)e^{\text{j}\mathit{\omega }t}$, the spatial part $E_z(r)$ of the electric field satisfies the Helmholtz equation:(1) $$\begin{array}{c}\hfill {\mathrm{\nabla }}^2{E}_z(r)+k^2{E}_z(r)=\text{j}\mathit{\omega }{\mathit{\mu }}_0{J}_z(r);r\in D\end{array}$$(1)

where $k=\mathit{\omega }\sqrt{{\mathit{\mu }}_0\mathit{ϵ}}$ represents wavenumber and ${\mathit{\mu }}_0$ and $\mathit{ϵ}(r)={\mathit{ϵ}}_0{\mathit{ϵ}}_r$ stands for permeability and permittivity of the medium, respectively. Electric field can be decomposed into an incident and scattered part, $E_z=E_z^{inc}+E_z^s$, where the incident field is propagated by source antennas, while the object is not present. Scattered wave is the reaction of object to incident wave. The incident field generated by source could be revealed as below:(2) $$\begin{array}{c}\hfill E_z^{inc}(r)=\text{j}\mathit{\omega }{\mathit{\mu }}_{0}I(\mathit{\omega })G_{2D}(r,r^{\prime })\end{array}$$(2)

the $I(\mathit{\omega })$ represents source current spectrum, while $G_{2D}(r,r^{\prime })=\frac{j}{4}{H}_0^{(2)}(k_0\mathit{\rho })$ is free-space scalar Green function and $\mathit{\rho }=|r-r^{\prime }|$. Based on theorem of volume equivalence,[19] $k$ in throughout of the region $D$ can be replaced by $k_0$, provided discontinuity’s effect due to presence of the object is taken into account by an equivalent current source. Thus $J_z(r)=O(r)E_z(r)$ where $O(r)=\text{j}\mathit{\omega }[\mathit{ϵ}(r)-{\mathit{ϵ}}_0]$ denotes the object function. The solution for Equation (1(1) $$\begin{array}{c}\hfill {\mathrm{\nabla }}^2{E}_z(r)+k^2{E}_z(r)=\text{j}\mathit{\omega }{\mathit{\mu }}_0{J}_z(r);r\in D\end{array}$$(1) ) is:(3) $$\begin{array}{c}\hfill E_z^s(r)=E_z(r)-E_z^{inc}(r)=\text{j}\mathit{\omega }{\mathit{\mu }}_0\int O(r^{\prime })E_z(r^{\prime })G_{2D}(r,r^{\prime })\text{d}{r}^{\prime }.\end{array}$$(3)

Figure 1. Configuration of problem for two dimensional structures. Transmitting/receiving antennas for multi-static configuration are placed on measurement curve $\mathrm{\Gamma }$.

3. Methods formulation: LSM & TR

3.1. Linear sampling method

The LSM is based on sampling within the investigation domain. A Fredholm linear integral equation is solved for each point. The subsequent Fredholm linear integral equation of the first kind is the preliminary point for LSM formulation [5](4) $$\begin{array}{c}\hfill {\int }_{{\mathrm{\Gamma }}_s}{E}_z^s(r_{{\mathrm{\Gamma }}_m},r_{{\mathrm{\Gamma }}_s})\mathit{\psi }(p_t,r_{{\mathrm{\Gamma }}_s})\text{d}{r}_{{\mathrm{\Gamma }}_s}=G_{2D}(r_{{\mathrm{\Gamma }}_m},p_t)\end{array}$$(4)

where $p_t\in D$ is a sampling point inside the investigation area and $r_{{\mathrm{\Gamma }}_s}$ and $r_{{\mathrm{\Gamma }}_m}$ are the source and measurement positions on $\mathrm{\Gamma }$, respectively. In the far field region, Green function can be approximated as ${lim}_{\mathit{\rho }\to \infty }{G}_{2D}(\mathit{\rho })\approx \frac{j}{4}\sqrt{\frac{2}{\mathit{\pi }k\mathit{\rho }}}exp(-j(k\mathit{\rho }-\frac{\mathit{\pi }}{4}))$. Parameter $\mathit{\psi }$ is a function that is known from solving the integral equation. Its Euclidian norm is small when $p_t$ is inside the scatterer, as $p_t$ travels close to its exterior contour or outside the object, the norm increases to a very large value. Therefore, the first step is to pixelate the investigated area. In the simulation process, we corrupt the input data with a specified noise in order to overcome inverse crime and to meet the conditions for practical measurement. As a result, $[E_z^s]$ will be substituted by $[E_z^{sn}]$. Then, assessment of noise power is reached when parameter $h$ is defined so that $\Vert [E_z^{sn}]-[E_z^s]\Vert <h$. Generally, $h$ could be originate from several source of uncertainty including that of noise. As the discretized equation is ill-conditioned, a regularization method with an optimum regularization parameter (${\mathit{\alpha }}_{p_t}^{\ast }$) is applied. Using singular value decomposition (SVD) $[E_z^{sn}]=[U][\mathrm{\Sigma }]{[V]}^{\ast }$, ${\mathit{\alpha }}_{p_t}^{\ast }$ is the root of the following equation [20](5) $$\begin{array}{c}\hfill \mathit{\vartheta }[\mathit{\alpha }(p_t)]=\sum _{i=1}^{q^n}\frac{{[\mathit{\alpha }(p_t)]}^2-h^2{({\mathit{\sigma }}_i^n)}^2}{{[{({\mathit{\sigma }}_i^n)}^2+\mathit{\alpha }(p_t)]}^2}{|{〈G_{p_t},u_i^n〉}_N|}^2\end{array}$$(5)

where ${\mathit{\sigma }}_i^n$ the diagonal elements of $[\mathrm{\Sigma }]$ are singular values and ${\{u_i^n\}}_{i=1}^N$ columns of $[U]$ are singular vectors. Parameter $q^n$ is the rank of $[E_z^{sn}]$ and ${〈G_{p_t},u{{}_i}^n〉}_N={\sum }_{j=1}^N{(G_{p_t})}_j{(u_i^n)}_j^{\ast }$. Once ${\mathit{\alpha }}_{p_t}^{\ast }$ is calculated, norm of $\mathit{\psi }$ can be obtained as [20], according to the Tikhonov regularization method:(6) $$\begin{array}{c}\hfill \Vert {\mathit{\psi }}_{p_t,{\mathit{\alpha }}^{\ast }(p_t)}{\Vert }_{{}_N}^2=\sum _{i=1}^{q^n}\frac{{({\mathit{\sigma }}_i^n)}^2}{[{({\mathit{\sigma }}_i^n)}^2+{\mathit{\alpha }}^{\ast }(p_t){]}^2}{|〈G_{p_t},u_i^n〉|}^2.\end{array}$$(6)

Choosing a monotonic indicator function for the norm of $\mathit{\psi }$, $IF=-log[\Vert {\mathit{\psi }}_{p_t,{\mathit{\alpha }}^{\ast }(p_t)}\Vert ]$, provides a tool to identify the scatterer boundary. Specific indicator function has been chosen here, different choices of indicator functions may influence the overall image however.

The basic form of LSM has been proposed for single frequency. However, recently in some studies,[21–24] multi-frequency LSM is proposed to overcome the weakness of LSM under sparse small aperture data. Guzina et al. [22] proposed two candidates for multi-frequency indicator function and concluded that serial indicator $\left(\sqrt{1/{\sum }_{f}I{F}_f^{-2}}\right)$ is not able to distinguish interior boundary of the scatterer from the exterior, when the define frequency bandwidth traverses one of the isolated eigenvalues of interior problem. In general, these eigenvalues are unknown since the contour of obstacle is not known before reconstruction. Contrary to serial form, the parallel indicator $\left(\sqrt{{\sum }_{f}I{F}_f^2}\right)$ is not sensitive to the presence of pertinent interior eigenvalues. Therefore, it seems to be robust in a generic scattering configuration. In order to omit the effect of possible singularities, the result of single frequency reconstruction could be merged to attain a normalized paralleled index as below [21]:(7) $$\begin{array}{c}\hfill \widehat{IF}=\frac{1}{N_F}\sum _{f=1}^{N_F}\frac{I{F}_f^2}{\begin{array}{c}\text{max}(I{F}_f^2(r))\\ r\in D\end{array}}\end{array}$$(7)

where $N_f$ is the number of frequencies.

3.2. Time interval

TR family methods recommend a distinctive prospect to solve electromagnetic focusing techniques. They have been employed in engineering problems such as acoustics for some years and more recently they have been applied to electromagnetics. TR methods can be interpreted as a spatiotemporal matched-filter method.[12] By considering Equation (3(3) $$\begin{array}{c}\hfill E_z^s(r)=E_z(r)-E_z^{inc}(r)=\text{j}\mathit{\omega }{\mathit{\mu }}_0\int O(r^{\prime })E_z(r^{\prime })G_{2D}(r,r^{\prime })\text{d}{r}^{\prime }.\end{array}$$(3) ) and using DWBA, the scattered field in the frequency domain can be represented as below [25](8) $$\begin{array}{c}\hfill E_{ms}^s(r_{{\mathrm{\Gamma }}_m},r_{{\mathrm{\Gamma }}_s})\approx \text{j}\mathit{\omega }{\mathit{\mu }}_0\int G_{2D}(r_{{\mathrm{\Gamma }}_m},r^{\prime })O(r^{\prime })E^{inc}(r^{\prime },r_{{\mathrm{\Gamma }}_s})\text{d}{r}^{\prime }.\end{array}$$(8)

In order to solve the problem, the investigation domain should be discretized into $P$ pixels. After discretizing and applying the reciprocity theorem:(9) $$\begin{array}{c}\hfill E_{ms}^s=\sum _{p=1}^P{\mathit{\tau }}_p{G}_{2D}(p,r_{{\mathrm{\Gamma }}_m})G_{2D}(p,r_{{\mathrm{\Gamma }}_s})\end{array}$$(9)

where ${\mathit{\tau }}_p=-\text{j}{\mathit{\omega }}^3{\mathit{\mu }}_0^{2}I(\mathit{\omega })[\mathit{ϵ}(p)-{\mathit{ϵ}}_0]$ is the scatterer strength of $p\text{th}$ pixel.

In match filter formulation, final migrated data can be interpreted as an image with a finite number of pixels. Suppose that every pixel has a main point scatterer that may be modelled in such a way that an object’s profile can be assumed as an impulse located at $p_t$, the center of the pixel. Mathematically, the transfer function of a match filter can be revealed as a complex conjugate of a presumed received waveform from the object under evaluation to which the filter is going to be matched. In cases when it is acceptable to suppose a linear scattering model, i.e. interactions between scattering centres can be neglected, it is possible to take account of the match filter of the point scatterer directly from the expected received signal waveform, and this makes it essential to conjugate the gathered data prior to application of the filtering process.[16] Therefore, in accordance with (9(9) $$\begin{array}{c}\hfill E_{ms}^s=\sum _{p=1}^P{\mathit{\tau }}_p{G}_{2D}(p,r_{{\mathrm{\Gamma }}_m})G_{2D}(p,r_{{\mathrm{\Gamma }}_s})\end{array}$$(9) ), the output result of the match filter could be presented as below:(10) $$\begin{array}{c}\hfill \mathit{\chi }(r)=\{E^{s\ast }(r_{{\mathrm{\Gamma }}_m},r_{{\mathrm{\Gamma }}_s})G_{2D}(r,r_{{\mathrm{\Gamma }}_m})\}\{G_{2D}(r,r_{{\mathrm{\Gamma }}_s})\};r\in D\end{array}$$(10)

where the highest value of function will occur at or in the neighbourhood of the scatterer point ($r=p_t$). In such a case, this is the image of the point scatterer.

From formula (10(10) $$\begin{array}{c}\hfill \mathit{\chi }(r)=\{E^{s\ast }(r_{{\mathrm{\Gamma }}_m},r_{{\mathrm{\Gamma }}_s})G_{2D}(r,r_{{\mathrm{\Gamma }}_m})\}\{G_{2D}(r,r_{{\mathrm{\Gamma }}_s})\};r\in D\end{array}$$(10) ), the right-hand side of the equation can be divided into two distinct terms and this division corresponds to its physical interpretation. The first is the electrical field of a current source $E^{s\ast }$. If we take into account time dependence of the received wave signal, then the source can be revealed as a TR $E^{s\ast }\Rightarrow E^s(-t)$. The generated field of this source is referred as back propagated electric field $E^{BP}$. The second term of the equation is the incident wave $E^{\mathrm{inc}}$.[16] In conclusion, when the investigation domain is pixelated, and TR procedure is used (i.e. the phase conjugation technique is applied), then complete expression for migrated data could be defined as suggested by Liao and Dogaru [26]:(11) $$\begin{array}{cc}\hfill \mathit{\chi }(r)=& \sum _{\mathit{\omega }}\sum _s\sum _m{E}_{ms}^{BP}(r;\mathit{\omega })E_s^{inc}(r;\mathit{\omega })\hfill \\ \hfill =& \sum _{\mathit{\omega }}\sum _s\sum _m{E}_z^{s\ast }(r_{{\mathrm{\Gamma }}_m},r_{{\mathrm{\Gamma }}_s},\mathit{\omega })G_{2D}(r,r_{{\mathrm{\Gamma }}_m},\mathit{\omega })G_{2D}(r,r_{{\mathrm{\Gamma }}_s},\mathit{\omega }).\hfill \end{array}$$(11)

4. Comparison of algorithms

4.1. PSF and resolution

Resolution can be considered as an appropriate indicator for quality of an imaging system or method; however, resolution needs mathematical definition. For many years, optic researches have used several definitions of resolution. In classical definitions, resolution means the ability of an imaging system or method to resolve a point in the vicinity of another point whereby these two points generate a scattered field of equal intensity. A common definition for resolution limit is proposed by Rayleigh’s criterion. In Rayleigh’s criterion, two neighbouring points are resolved if the maxima of the scattered wave function of one point occurs at the first minima of the scattered wave function of other points, or at a point far away from the minima. If resolution can be achieved beyond these limits then it is termed a super resolution.[3]

The definition of resolution used in this study is the Rayleigh criterion. According to the Rayleigh criterion, when a single point is imaged by a system or method, due to diffraction, the image will not be a single point and the image has non-zero width, while the mathematical width of a point is zero. The resulting distribution is termed point spreading function (PSF). In PSF, the image index has a peak at the location of an imaged point and it falls for the surrounding points. In an ideal imaging system, the image has its maximum in the location of a point and for any other points in the vicinity, the value of an index is zero, while in practical imaging methods, the value of an index for the neighborhood is much lower than that of the point, however the value is not zero.[4]

In order to compare resolutions of the two methods according to Rayleigh resolution criterion, a single point needs to be reconstructed by both methods. Mathematical definition of a point is as follows; the point to be recovered is a dielectric cylinder with a circular cross section located at the centre of a coordinate. Radius of the circle is 7.5 mm. The dielectric relative permittivity is ${\mathit{ϵ}}_r=1.8$ and it is lossless. Source and measurement locations on a circle are $24$ equally spaced in 15 m of the object. All simulations are done in TM polarization.

Figure 2 shows results of the reconstruction from single frequency synthetic data obtained by method of moment forward solver at 3 GHz. Based on this figure and the Rayleigh criterion, resolution of the two methods can be calculated as shown below:(12) $$\begin{array}{cc}\hfill LS{M}_{SF}& :x_3-x_1=0.8011\mathit{\lambda }\hfill \end{array}$$(12) (13) $$\begin{array}{cc}\hfill T{R}_{SF}& :x_5-x_1=0.26703\mathit{\lambda }.\hfill \end{array}$$(13)

Figure 2. Normalized PSF of (a) LSM and (b) TR in single-frequency free-space reconstruction. (c) Comparison of the profiles of the single-frequency PSFs for LSM and TR. Extremums of main lobe for each method have been denoted with appropriate numbers.

Figure 2 and equations show that the TR method had finer and thus better resolution. Results demonstrate that TR image quality surpassed that of LSM for imaging a single point. However, analysis of the figure determines that the null of LSM was obviously deeper than that of TR. This cannot be considered an advantage but it represents some kind of better contrast for the LSM method in comparison with TR.

To investigate the effect of frequency diversity on resolution of algorithms, the above-mentioned point scatterer is imaged by the multifrequency signal excitation that consists of a combination of 11 single tone signal from 2.5 to 3.5 GHz, with equal step increments of 100 MHz. The resulting images are depicted in Figure 3 ($\mathit{\lambda }$ is the free-space wavelength at the centre frequency 3 GHz)(14) $$\begin{array}{cc}\hfill LS{M}_{MF}& :x_3-x_1=0.76427\mathit{\lambda }\hfill \end{array}$$(14) (15) $$\begin{array}{cc}\hfill T{R}_{MF}& :x_5-x_1=0.26703\mathit{\lambda }.\hfill \end{array}$$(15)

Figure 3. Normalized PSF of (a) LSM and (b) TR in multi-frequency free-space reconstruction. (c) Comparison of the profiles of the multi-frequency PSFs for LSM and TR. Extremums of main lobe for each method have been denoted with appropriate numbers.

Results of multi frequency demonstrate some similarity to those of single frequency. It means that this attribute of methods almost independent from number of applied frequencies though the overall quality of the resultant images in this case have been improved.

According to Figures 2 and 3, it can be said that the maximum peak for the LSM is sharper in comparison to the TR but next values of energy do not drop off as fast. So it is difficult to assess resolution in a clear way. This can go towards the point that a resolution measure in isolation is a hard measure to ascertain quality. Therefore, in order to clarify about this point, separation resolution of algorithms between two canonical objects [3,27] is reported in Figure 4 by sweeping the distance ($\mathit{\delta }$) between two identical dielectric cylinders with radius of 25 mm and ${\mathit{ϵ}}_r=1.8$. Simulations for this configuration are carried out at 3 GHz with $24$ antennas uniformly located at a circle at distance of 1.5 m.

Figure 4. Separation resolution study: (a) resolution ratio for various distances between two dielectric cylinders (b) representation of profiles corresponding to $\mathit{\delta }=25\text{mm}$ (c), (d) and (e) reconstructed images obtained from LSM of cylinders scenario with $\mathit{\delta }=5,15,25\text{mm}$, respectively, and (f), (g) and (h) reconstructed images obtained from TR under similar conditions (The solid blue lines in each reconstruction plot indicate the true target boundary).

Semenov et al. [27] used a cross-sectional profile across the centres of two cylinders and considered the ratio of the minimum value between them ($U_{\text{min}}$) to the first maximum closest to this minimum ($U_{\text{max}}$) as a metric to measure resolvability of experimental imaging systems. Following the same procedure and from Figure 4(a) and (b), it can be stated that the separation resolution for TR obviously better than that for LSM. This result is consistent with the previous discussion based on the Rayleigh criterion.

4.2. Correlation index

Another common standard for comparing methods is similarity between the reconstructed image and the real target. An indicator is needed to compare performance of the two methods. The Rayleigh resolution index was used to compare the quality of single point reconstruction between the two methods. Although all reconstructed images by different imaging methods consist of single points, the Rayleigh index is not appropriate for making comparisons of total image quality. Accordingly, it was necessary to propose a new indicator to determine image quality.

The focus of this study was to investigate qualitative methods to recover the shape and location of objects under evaluation. Based on this, the results of recovery related to whether or not a sample point belonged to an investigated object. This problem resembles a binary form: for a point on the object $1$ was assigned and for others, $0$. The indicator was then introduced based on this description. A $N\times N$ matrix was considered to image an area, which is discretized to $N\times N$ cells. In This matrix, elements that correspond to the surrounding free-space are equal to $0$ while a point corresponding to the object is 1. This matrix is called the real target matrix.

For each method of imaging, when an image is constructed, a definition for the threshold and binarized images were defined accordingly. Two-dimentsional gradients of images demonstrated that maximum variation occurred next to the boundary between the scatterer and its host medium. Thus, normalized values of pixels with the maximum gradient could be chosen as a proper threshold for binarization. Based on this threshold, values below the threshold were assigned to $0$ and values higher than the threshold were assigned to 1 (object). Then, the correlation indicator was defined as described below; the indicator equals correlation of the two matrices; the real target matrix and the image matrix.

The quality of the whole image was then determined based on the proposed indicator. Definition of correlation for the two binary matrices $A$ and $B$ is described below:(16) $$\begin{array}{c}\hfill corr2(A,B)=\frac{{\sum }_m{\sum }_n(A_{mn}-\overline{A})(B_{mn}-\overline{B})}{\sqrt{(\sum _m\sum _n{(A_{mn}-\overline{A})}^2)(\sum _m\sum _n{(B_{mn}-\overline{B})}^2)}}\end{array}$$(16)

where $\overline{A}$ and $\overline{B}$ are the mean values (average of matrix elements) of their corresponding matrices.

In this section, objects were reconstructed based on measurement data i.e. Fresnel institute data-set. This data-set is related to an aspect-limited configuration in which the source is transferred along the circumference of a circle, while the spacing of the steps is 10${}^{\circ }$ per step and for each wave illumination , the measurement receiver antenna was moved in $5^{\circ }$ steps along a $240$ arc which excluded the $120$ angular sector centred around the incidence direction. Four configurations were considered for imaging. These configurations are depicted in Figure 5.

Figure 5. Configuration of test structures: (a) an off-centred dielectric cylinder (b) a rectangular conductor (c) two dielectric, symmetrical located on axis, cylinder (d) a u-shaped conductor. All dielectric cylinders are with circular cross section of radius 15 mm and the real part of the relative permittivity ${\mathit{ϵ}}_r=3\pm 0.3$.[18]

Imaging was done with frequency diversity in a set of frequencies: 2, 4, 6 & 8 GHz. For the u-shaped structure, due to its complicated shape, more resolution was needed. Hence for this scenario, reconstruction was repeated at higher frequency bands, i.e. 10, 12, 14 & 16 GHz and named the fifth case study. Results of the reconstructed images are depicted in Figures 6 and 7. In these figures, $3\mathit{\lambda }\times 3\mathit{\lambda }$ investigation domain has been pixelated into $91\times 91$ pixels in all reconstructions and $\mathit{\lambda }$ denotes the free-space wavelength at the frequency 2 GHz.

Figure 6. LSM Recovery.

Figure 7. TR recovery.

As aforementioned, binarization was necessary to evaluate quality. The resultant binary images from LSM and TR are shown in Figures 8 and 9, respectively. Also, binary images of the test scenarios under evaluation, which are necessary tools for obtaining a correlation index in each case, should be calculated with the same number of pixels.

Figure 8. Binary images obtained from Figure 6 with Gradient-based threshold set to: (a) 0.6595, (b) 0.6377, (c) 0.7976, (d) 0.6966 and (e) 0.5813.

Figure 9. Binary images obtained from Figure 7 with Gradient-based threshold set to: (a) 0.3833, (b) 0.3909, (c) 0.4531, (d) 0.7232 and (e) 0.4728.

Based on these results, the two methods were compared by the suggested correlation index. The correlation index is calculated for five case studies. The results are illustrated in Figure 10.

Figure 10. Correlation index for five case studies.

Based on the proposed threshold, the quality of total image for LSM was much better than that of TR. This is a fascinating result because when quality is compared between resolutions, TR was much better than LSM. This demonstrates that resolution parameter is not sufficient enough to determine quality of an imaging algorithm. By comparing the indices of fourth and fifth case studies, it can be stated that increasing the imaging frequency evidently increased the quality of TR, while it had no effect on the LSM image.

At this point, it should be noted that the value of the correlation index was highly dependent on the value of the threshold. As discussed in several references, the decision threshold criterion is an open dispute in imaging methods. Usually, if another threshold criterion is adopted then correlations will show different results. In order to present the importance of threshold value setting and evaluating the maximum capability of algorithms in retrieving known scenarios, an F-measure score can be exploited. The F-measure score is a quite common indicator in imaging and vision.[17,28] This score takes into account the precision as well as recall of a detector with following definition which is a function of threshold ($\mathcal{T}$) and seeks the argument that achieves the maximum F-measure:(17) $$\begin{array}{c}\hfill Q=\begin{array}{c}\max \\ 0\le \mathcal{T}\le 1\end{array}2\frac{precision(\mathcal{T})recall(\mathcal{T})}{precision(\mathcal{T})+recall(\mathcal{T})}\end{array}$$(17)

where(18) $$\begin{array}{c}\hfill precision(\mathcal{T})=\frac{|\{I=1\}\cap \{\mathcal{X}=1\}|}{|\{I=1\}|}\end{array}$$(18) (19) $$\begin{array}{c}\hfill recall(\mathcal{T})=\frac{|\{I=1\}\cap \{\mathcal{X}=1\}|}{|\{\mathcal{X}=1\}|}\end{array}$$(19) $I$ and $\mathcal{X}$ denote the reconstructed image after binarization and true target for corresponding scenario, respectively. Table 1 presents a summary of quantified assessment of algorithms based on Fresnel institute experiments. Although ideal threshold determination is not possible in practice, in order to study F-measure score, here the known real target is considered, which provides the maximum accessible performance of each method under different conditions which in turn provides an even-handed quantitative measure to assess imaging quality.

Table 1. Quantitative assessment of image reconstruction using of the F-measure. Higher recall-precision performance leads to a higher score.

Case-study	Method	Ideal threshold ($\widehat{\mathcal{T}}$)	Quality score ($Q$)
#1	LSM	0.875	0.8276
	TR	0.875	0.8667
#2	LSM	0.825	1.0
	TR	0.975	0.8333
#3	LSM	0.85	0.6825
	TR	0.8	0.5906
#4	LSM	0.875	0.6275
	TR	0.55	0.4958
#5	LSM	0.875	0.5806
	TR	0.675	0.6984

Table 1 reveals that in second, third and fourth case studies LSM shows a better upper limit of performance than TR. Since these scenarios contain PEC or multiple dielectric targets, these results are due to Born approximation usage in TR formulation which assumes weak scattering nature model and neglect multiple scattering events. Indeed these results confirm the claim that LSM can overcome Born-approximated methods, when multi-bounce events are considered. However, in the first case, TR acts quite similar to LSM. Even TR performance is slightly better than LSM in retrieving single low-contrast dielectric cylinder. In fifth case, TR has better score than LSM in exploiting high frequency raw data.

5. Conclusion

This comprehensive study was done to compare image quality of TR and LSM. Two indexes were introduced to compare imaging quality, the Rayleigh resolution for point imaging quality and a proposed index called the correlation index for total image quality. In terms of the resolution parameter, results showed that quality was better in TR than LSM. Resolutions of both methods were similar for single frequency and multifrequency imaging. These results demonstrate that the TR method could produce finer points than LSM. In order to compare total image between the reconstruction methods, a new indicator called the correlation index was recommended whereby imaging was made using real measurements. Simulation showed that image quality for LSM was much better than that of TR. In addition, the recovery quality for LSM was somehow constant and not dependent on an object’s structure or material, while this was not the case for TR. Increasing imaging frequency improved the quality of a TR image, but had no effect on quality of the LSM image.

Acknowledgements

The authors would like to thank Dr M. Dehmollaian at University of Tehran, Tehran, Iran, for helpful consultations on various steps of this study. The authors would also like to acknowledge the anonymous reviewers for their constructive comments and suggestions to improve this work.

Notes

No potential conflict of interest was reported by the authors.