Open AccessArticle

An Accurate Millimeter-Wave Imaging Algorithm for Close-Range Monostatic System

Xinyi Nie

^1,2

Chuan Lin

¹,

Yang Meng

Anyong Qing

^1,*

Jan K. Sykulski

⁴

and

Ian D. Robertson

School of Electrical Engineering, Southwest Jiaotong University, Chengdu 611756, China

School of Electronic and Electrical Engineering, University of Leeds, Leeds LS2 9JT, UK

School of Optoelectronic Engineering, Chongqing University of Posts and Telecommunications, Chongqing 400065, China

⁴

School of Electronics and Computer Science, University of Southampton, Southampton SO17 1BJ, UK

Author to whom correspondence should be addressed.

Sensors 2023, 23(10), 4577; https://doi.org/10.3390/s23104577

Submission received: 7 March 2023 / Revised: 10 April 2023 / Accepted: 24 April 2023 / Published: 9 May 2023

(This article belongs to the Section Electronic Sensors)

Download

Browse Figures

Figure 1
The geometry of imaging system. "> Figure 2
The reconstruction results of point scatterer in scale normalized to the maximum. (a) The reflectivity contour projected to xz plane. (b) The 3-D view of the reflectivity projected to the xz plane. "> Figure 3
Profiles along x(z) in dB scale normalized to the maximum. "> Figure 4
Profiles along x = z in dB scale normalized to the maximum. "> Figure 5
Profiles along y in dB scale normalized to the maximum. "> Figure 6
The 2-D object under test. "> Figure 7
The reconstruction results. (a) Sheen’s method. (b) Meng’s method. (c) The proposed method. "> Figure 8
The reconstructed and true normalized magnitude along the line z = 0 m. "> Figure 9
The CAD model of three blocks within different ranges. "> Figure 10
The front view of the reconstruction results. (a) Sheen’s method. (b) Meng’s method. (c) The proposed method. "> Figure 11
The side view of the reconstruction results. (a) Sheen’s method. (b) Meng’s method. (c) The proposed method. "> Figure 12
The top view of the reconstruction results. (a) Sheen’s method. (b) Meng’s method. (c) The proposed method. "> Figure 13
Profiles along z = 0 m in <a href="#sensors-23-04577-f010" class="html-fig">Figure 10</a>. "> Figure 14
Profiles along z = 0 m in <a href="#sensors-23-04577-f011" class="html-fig">Figure 11</a>. "> Figure 15
Photograph of the target under test. "> Figure 16
Zoomed-in reconstructed result. "> Figure 17
Profiles along the red horizontal line in <a href="#sensors-23-04577-f016" class="html-fig">Figure 16</a>. "> Figure 18
Profiles along the red vertical line in <a href="#sensors-23-04577-f016" class="html-fig">Figure 16</a>. "> Figure 19
The scenario with a volunteer scanned by prototype. "> Figure 20
The reconstruction results (a) by Meng’s method and (b) by the proposed method. ">

Versions Notes

Abstract

An efficient and more accurate millimeter-wave imaging algorithm, applied to a close-range monostatic personnel screening system, with consideration of dual path propagation loss, is presented in this paper. The algorithm is developed in accordance with a more rigorous physical model for the monostatic system. The physical model treats incident waves and scattered waves as spherical waves with a more rigorous amplitude term as per electromagnetic theory. As a result, the proposed method can achieve a better focusing effect for multiple targets in different range planes. Since the mathematical methods in classical algorithms, such as spherical wave decomposition and Weyl identity, cannot handle the corresponding mathematical model, the proposed algorithm is derived through the method of stationary phase (MSP). The algorithm has been validated by numerical simulations and laboratory experiments. Good performance in terms of computational efficiency and accuracy has been observed. The synthetic reconstruction results show that the proposed algorithm has significant advantages compared with the classical algorithms, and the reconstruction by using full-wave data generated by FEKO further verifies the validity of the proposed algorithm. Finally, the proposed algorithm performs as expected over real data acquired by our laboratory prototype.

Keywords:

concealed weapon detection; Fourier transform technique; method of stationary phase; microwave imaging; national security

1. Introduction

Terrorist attacks around the world are still at a high level because regional and religious conflicts have been intensifying [1]. Crowded places such as transportation hubs (airports, railway stations, metro stations, etc.), entertainment venues, government agencies, and schools are primary targets. These threats can result in not only loss of life and property, but also social anxiety and panic. Therefore, it is extremely urgent and challenging to ensure public security. It is critical that prohibited items such as weapons and explosives are detected before they enter those crowded public areas.

Conventional approaches, including metal detectors, infrared detectors, and X-ray systems [2], already exist. To some extent, these approaches are effective but suffer some drawbacks in practical applications. Metal detectors can only detect metal targets. Meanwhile, many modern hazardous articles are now made using advanced technologies, such as plastic, ceramics (e.g., knives), and liquid explosives, that cannot be detected by using metal detectors. Infrared detectors can realize imaging but are seriously affected by the environment, and the image quality is inherently not high [3]. X-ray systems have remarkable performance in terms of penetrability and image quality, so they are a good detection method for carry-on luggage. However, even low-dose X-ray cannot be accepted as a means of personnel security inspection, since X-rays are ionizing radiation and have a cumulative effect.

Millimeter waves are high-frequency electromagnetic waves with relatively short wavelengths in the 1–10 mm range and offer good penetrability. They are capable of providing good spatial resolution and penetrating common clothing and packaging materials. Moreover, millimeter waves are nonionizing. Thus, a millimeter-wave system is well suited for personnel security screening.

In terms of imaging techniques, the back-projection (BP) algorithm is one of the most accurate methods in the spatial domain [2,4,5,6,7,8,9]. Boag [8,9] attempted to speed up the BP algorithm by using multilevel domain decomposition. Unfortunately, the demonstrated acceleration was achieved at the cost of accuracy. As a matter of fact, at acceptable accuracy for engineering practice, the demonstrated efficiency may even be worse than that of the original BP algorithm. Therefore, its heavy computational burden is still a severe issue, and it is impractical for real-time security checks. It is noted that the scattering data are the convolution of reflectivity function and Green’s function in spatial domain, which can be more efficiently handled as the multiplication of their Fourier transforms in the wavenumber domain and converted back through FFT technique. Therefore, many efficient wavenumber-domain algorithms have been developed, such as the holographic algorithm [10,11,12,13,14,15,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30], range Doppler (RD) algorithm [31,32], chirp scaling (CS) algorithm [33,34], and range migration algorithm (RMA) [35,36,37], among which the holographic algorithm has been being heavily studied for its proven commercial potential. In particular, monostatic configuration is often assumed for simplicity and efficiency. In this paper, this practice is also adopted.

Close-range imaging systems based on a monostatic configuration have been studied in many works [7,11,13,16,17,18,19,20,21,22,23,24,25,26,27,28,29,32,36,38]. Most of the studies neglected the amplitude terms of both incident and scattered spherical waves for simplicity [11,16,17,18,19,20,21,22,23,32,36]. Such negligence might result in some severe issues because the propagation loss is very important while the undulate scale of objects in security screening is comparable to the range of distances in close-range millimeter-wave imaging [10,11,35,39,40].

A recent remedy partly took propagation loss into account [25,26]. The scattered signal is treated as a spherical wave, which is closer to the actual physical model. Therefore, the quality of the reconstruction image can be improved to some extent. Nevertheless, the physical model is still not consistent with the actual physical model.

To eradicate the problem, the dual path propagation loss has to be adequately taken into account [4,7,38]. It is straightforward for BP algorithm to compensate the dual path propagation loss by directly including the factor R² in its formulation [4,7]. However, the computation burden is too heavy. A similar computationally inefficient approach is observed in the range stacking algorithm, where the scattering data were transformed into the slant range spatial domain and then multiplied by an amplitude factor R² [38]. Surprisingly, a factor R⁴, inconsistent with electromagnetic theory, appeared in [27].

In research to date, the dual path propagation loss,

1 / R^{2}

, has hardly been considered in monostatic holographic algorithms. A close look at the mathematical formulation reveals that the mathematical methods in holographic algorithms, such as spherical wave decomposition [11] and Weyl identity [25], stop working if the term

1 / R^{2}

is involved. The method of stationary phase (MSP) is therefore chosen to overcome this difficulty. A holographic reconstruction algorithm for a monostatic system with full consideration of dual path propagation loss for close-range millimeter-wave imaging is accordingly proposed in this paper.

The proposed method treats incident waves and scattered waves as spherical waves with more rigorous amplitude terms as per electromagnetic theory. That is to say, the propagation loss is more accurately taken into account. Thus, the physical method is more consistent with the actual scenario. Consequently, it can achieve a better focusing effect for multiple targets in different range planes. The proposed algorithm is derived through MSP. Meanwhile, FFT is also employed to transform spatial domain scattering data to the wavenumber domain for efficient imaging. The proposed algorithm has been validated by improved reconstructed images.

It has to be pointed out that the dual path propagation loss in MIMO and bistatic systems has already been addressed [35,41,42,43]. However, as far as the authors know, our imaging formula is by no means a simple degeneration but a brand-new formulation. More importantly, our algorithm significantly outperforms those degenerated algorithms when the focusing effect plays a role. Due to page limitations, details will be published elsewhere separately [44].

2. Theory

The geometry of the imaging system is shown in Figure 1. The transmitting and receiving antennas are assumed to be positioned at (x′, y′, z′) on the plane of y = y′. An arbitrary point on the target is represented as (x, y, z) with reflectivity function f(x, y, z).

In essence, a transmitter emanates a spherical wave, which will illuminate the object and be reflected by the object. Each point on the object can be regarded as a secondary source transmitting spherical waves. Consequently, the echo signal measured by the corresponding receiver would be the superposition of all the spherical waves emanating from all the points on the object.

Supposing that the scattering process satisfies the Born approximation according to electromagnetic theory [13,45], the measured backscattered data is

s (x^{'}, z^{'}, k) = \int \int \int_{V} f (x, y, z) \frac{\exp (- j 2 k R)}{{(4 π R)}^{2}} d x d y d z

(1)

where

j = \sqrt{- 1}

k = 2 π / λ

is the spatial wave number, and

λ

is the wavelength. The term

\exp (- j 2 k R) / R^{2}

stands for the strict spherical wave expression containing the incident wave and scattered wave. It is noted that the constant term,

1 / {(4 π)}^{2}

, will be omitted for convenience since the constant term has no effect on imaging results. The distance from the transceiver to the object is

R = \sqrt{{(x - x^{'})}^{2} + {(y - y^{'})}^{2} + {(z - z^{'})}^{2}}

(2)

The Fourier transform of the measured backscatter data

s (x^{'}, z^{'}, k)

can be expressed as

S (k_{x}, k_{z}, k) = \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} s (x^{'}, z^{'}, k) e^{- j k_{x} x^{'}} e^{- j k_{z} z^{'}} d x^{'} d z^{'}

(3)

Substituting (1) into (3) yields

\begin{array}{l} S (k_{x}, k_{z}, k) = & \int \int \int_{V} f (x, y, z) \\ \cdot \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} \frac{\exp (- j 2 k R)}{R^{2}} e^{- j k_{x} x^{'}} e^{- j k_{z} z^{'}} d x d y d z d x^{'} d z^{'} \end{array}

(4)

For the sake of simplicity, the double integral about

x^{'}

and

z^{'}

in (4) can be represented as

E (k_{x}, k_{z}, k) = \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} \frac{\exp (- j 2 k R)}{R^{2}} e^{- j k_{x} x^{'}} e^{- j k_{z} z^{'}} d x^{'} d z^{'}

(5)

The asymptotic expansion of the integral of (5) can be obtained by the MSP as [36,41,45]

E (k_{x}, k_{z}, k) = \frac{j π}{k} \frac{1}{y - y^{'}} e^{- j (\sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}} (y - y^{'}) + k_{x} x + k_{z} z)}

(6)

Please see the Appendix A for more details. Substituting (6) into (4), we have

\begin{array}{l} S (k_{x}, k_{z}, k) = & \int \int \int_{V} \frac{1}{y - y^{'}} f (x, y, z) \\ \cdot \frac{j π}{k} e^{- j (\sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}} (y - y^{'}) + k_{x} x + k_{z} z)} d x d y d z \end{array}

(7)

Let

k_{y} = \sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}}

(8)

Then, (7) can be re-formulated as

\begin{array}{l} \frac{k}{j π} S (k_{x}, k_{z}, k) e^{- j k_{y} y^{'}} \\ = \int \int \int_{V} \frac{f (x, y, z)}{y - y^{'}} e^{- j (k_{x} x + k_{y} y + k_{z} z)} d x d y d z ≜ F (k_{x}, k_{y}, k_{z}) \end{array}

(9)

where

F (k_{x}, k_{y}, k_{z})

is the Fourier transform of

f (x, y, z) / (y - y^{'})

Thus, the object can be reconstructed as

\begin{array}{l} f (x, y, z) = \frac{y - y^{'}}{j π} \\ \cdot {FT}_{3 D}^{- 1} \{{[k \cdot {FT}_{2 D} {[s (x^{'}, z^{'}, k)]|}_{x^{'}, z^{'}} e^{- j k_{y} y^{'}}]}_{Interpolation}\} \end{array}

(10)

where FT_2D represents the 2-D Fourier transform and

{FT}_{3 D}^{- 1}

represents the 3-D inverse Fourier transform. It should be noted that

y^{'}

on the right-hand side of (10) is the transceiver plane. Therefore, it is accurately known in both measurement of backscattered data and device under test (DUT) profile reconstruction once the geometry of the imaging system is fixed.

The proposed method takes into account the propagation loss of the spherical wave in free space. Therefore, the proposed algorithm is based on a more precise physical model, which will bring an improvement in the reconstruction results. It should also be pointed out that (10) will be close to the traditional method if the DUT plane is parallel to the image plane. However, a parallel DUT plane is ideal, so it is hardly possible for practical personnel screening. Additionally, a DUT is intermediate rather than far in practical personnel screening. Our algorithm will keep its advantage in practical personnel screening when the wavy DUT distance plays an important role.

3. Results and Discussions

Numerical simulations and experiments are used to test the performance of the proposed algorithm. Different kinds of test scenarios are purposely employed. A point target simulation is used to test resolution. Three rectangle objects are used to test the performance of the algorithm. The capability of 3-D imaging is presented by 3-D full-wave data. Finally, reconstruction with real data also validates the performance.

In order to verify the improvement of our proposed method, Sheen’s method [11] and Meng’s method [25] are chosen as competitors. All the algorithms are realized with self-developed MATLAB codes. The computational platform is a desktop computer with a regular Intel 64-bit 3.19-GHz CPU and 32 GB RAM.

In the numerical experiments, the simulated scattered data were obtained through (1). The system parameters used to generate the synthetic scattered data are presented in Table 1

3.1. Point Spread Function

The point spread function (PSF) quantifies the capability of a system to image an arbitrary point scatterer. Meanwhile, it can also be considered as an essential approach to evaluate imaging algorithms, since point scatterer simulations provide distinctive insight into other aspects of image quality unavailable by realistic target.

A point scatterer is placed at (0 m, 0.4 m, 0 m) in front of the transceiver plane. The transceiver plane is located at y = 0. In order to better address the point object, the equidistant sample is 2 mm instead of 5 mm, which strictly satisfies the Nyquist sampling theorems. Figure 2 shows the reconstructed PSF projected to the xz plane. It can be seen that the energy of the scattered field is mainly focused on (0 m, 0 m) and that the spreading phenomenon occurs.

The normalized PSF along the x and z direction in the y = 0.4 m plane is shown in Figure 3. The width of both of the amplitude lines at −4 dB [46,47] is 5.94 mm. The normalized PSF along x = z in the y = 0.4 m plane is shown in Figure 4. The width of the amplitude line at −4 dB is 4.35 mm. The PSF along y is shown in Figure 5. The width of the amplitude line at −4 dB is 24.3 mm. The resolution of the PSF is thus 5.94 mm along x and z, 4.35 mm along the diagonal at 45° relative to x, and 24.3 mm in the range direction. Theoretically, the lateral resolution and range resolutions are about

λ_{c} / 2 \approx

5 mm and

c / (2 B) \approx

25.8 mm [35,36,48]. Generally speaking, the numerical resolutions from PSF are consistent with the theoretical ones.

3.2. The Synthetic Data

Three identical rectangular objects, which are considered the objects under test, are illustrated in Figure 6. The length and width of the rectangular objects are 24 mm and 48 mm, respectively. The centers of rectangles are placed at (−0.09 m, 0.3 m, 0 m), (0 m, 0.4 m, 0 m), and (0.09 m, 0.5 m, 0 m) from left to right, respectively. The reflection coefficient of the object under test is 1.

The rectangles are reconstructed by Sheen’s method, Meng’s method, and the proposed method. The normalized results are projected to xz plane as shown in Figure 7. All three rectangles are well reconstructed by the three methods. From a visual point of view, Figure 7a,b are quite similar, while Figure 7c is better than the other two. More specifically, the reconstructed rectangles centered at y = 0.3 m are similar across all three methods, but the reconstructed rectangles centered at y = 0.4 m and y = 0.5 m are better for the proposed method than for the other two methods. The magnitude along the line z = 0 m is plotted in Figure 8. It is evident that the reconstructed profiles of the two rectangles on the right by the proposed method are closer to the true profiles, which means that the proposed method can better focus when there are multiple objects within different ranges.

In order to evaluate the algorithms’ performance fairly, the same quantitative performance indicators, namely, the root mean square error (RMSE) and the structural similarity index measure (SSIM), are selected as evaluation indexes. The RMSE is defined as

RMSE = \sqrt{\frac{1}{m} \sum_{i = 1}^{m} {[f_{r e a l} (i) - f_{c a l} (i)]}^{2}}

(11)

where

f_{r e a l} (i)

and

f_{c a l} (i)

are the real and calculated values at the i-th pixel in the true and recovered image, respectively.

And the SSIM is evaluated as

SSIM = {[l (i)]}^{α} \cdot {[c (i)]}^{β} \cdot {[s (i)]}^{γ}

(12)

where

l (i)

c (i)

, and

s (i)

are the luminance comparison function, the contrast comparison function, and the structure comparison function, respectively [49]. In addition,

α

β

, and

γ

are 1 in this paper. When the recovered image is closer to the reference image, the SSIM is closer to 1 and RMSE is closer to 0.

The RMSE and SSIM are computed as shown in Table 2. Hence, it can be concluded that the more accurate physical model does yield a better performance of the proposed algorithm in terms of imaging quality, especially RMSE. Meanwhile, computational efficiency of the proposed algorithm drops only slightly as shown in Table 3. It is noted that there is a significant improvement in Figure 7 and Figure 8 but only very modest improvements in SSIM and RMSE. That is because the area of the object is a small proportion of the total imaging area.

3.3. The Full-Wave Data

In this section, the scattered electromagnetic field is computed by the commercial electromagnetic computation software FEKO. The CAD model of the target is presented in Figure 9. The target is assumed to be perfect electric conductors (PECs). The target is composed of three metal blocks of thickness 20 mm. The length and width of the block objects are 24 mm and 48 mm. The frontal center of the three block objects (namely A, B, and C) are located at (−0.09 m, 0.3 m, 0 m), (0 m, 0.4 m, 0 m), and (0.09 m, 0.5 m, 0 m), respectively.

The front, side, and top views of the reconstruction results are shown in Figure 10, Figure 11 and Figure 12, respectively. Overall, the quality of reconstruction results by the proposed algorithm is better than those by the other two methods. Specific to each block, block A is the best reconstructed of the three blocks by all three methods since it is closest to the transceiver plane. The reconstruction results of block A are almost the same across all three methods. For blocks B and C, it is evident that the reconstruction results by the competing methods become worse with the increasing distance between the blocks and the transceiver plane. However, the proposed method is less sensitive to distance, so its reconstruction results for blocks B and C are much better than those of the other two methods. Furthermore, the normalized magnitude along z = 0 in Figure 10 and Figure 11 is plotted in Figure 13 and Figure 14. Likewise insensitivity to distance, and accordingly better reconstruction results as shown in Figure 13 and Figure 14, are also very easy to observe. The proposed method better reconstructed the farther parts than did the other two methods. The SSIM and RMSE values given in Table 4 further confirm the superior performance of the proposed algorithm against its competitors. Once again, the proposed algorithm only suffers very little deterioration in computational efficiency as shown in Table 5.

3.4. Experimental Results

A millimeter-wave imaging prototype was built in the laboratory, which can collect real data to test in-house-developed imaging algorithms. There are 157 equivalent sampling points with the same spacing of 5 mm along the horizontal direction (x-direction) and 397 equivalent sampling points with 5 mm steps along the vertical direction (z-direction) in the imaging area of 0.8 m × 2 m. The operation frequency ranges from 27 GHz to 32.8 GHz with 220 equidistant sampling points. The acquisition time is only 2 s. The full details of the prototype in the laboratory can be found in [25].

The target contains three groups of perfectly conducting strips of length 100 mm as shown in Figure 15 Both the horizontal and diagonal groups have five sets of strips. The space between neighboring sets is 50 mm. The widths of strips in each set are 10 mm, 7 mm, 6 mm, 5 mm, and 4 mm, and the space between any two neighboring strips in the same set is the same as the strip width. However, the vertical group involves one more set of strips of 20 mm width and spacing. The test target is positioned in front of the transceiver plane center at a distance of 0.4 m.

The reconstruction result is given in Figure 16, which shows only the object under test for better visualization. The profiles along the red lines in Figure 16 are plotted in Figure 17 and Figure 18. It is evident that the strips of 5 mm and greater can be easily distinguished. Hence, the actual resolution of the proposed algorithm for the prototype is about 5 mm.

Finally, another more complex real scenario is presented. A person who carried a cleaver, a fruit knife, and a cellphone, all hidden under a down jacket, was scanned by the prototype as shown in Figure 19. The distance between the volunteer and the antenna plane is also 0.4 m. The reconstruction results are given in Figure 20, which are marked from the circles and squares. All the concealed objects are detected and identifiably presented in the reconstructed images. The SSIM between the reconstruction results by Meng’s method and the proposed method is 0.8852.

4. Conclusions

In this paper, an electromagnetic imaging algorithm for a close-range monostatic system with dual path propagation loss has been presented for use in personnel screening for security applications. In this application, the scale of objects under scrutiny is comparable to the range of distances. The physical model treats incident waves and scattered waves as spherical waves with more rigorous amplitude term as per electromagnetic theory. As a consequence, the mathematical method in classical algorithms, such as spherical wave decomposition and Weyl identity, cannot handle the mathematical model derived from the physical model. The proposed algorithm is derived through MSP and is more accurate than the algorithms without dual path propagation loss. Consequently, the proposed method can achieve a better focusing effect for multiple targets in different range planes.

To verify the proposed algorithm, simulations and experiments have been carried out. PSF quantifies the algorithm resolution to be approximately 5.94 mm in the x and z direction and 24.3 mm in the range direction. Although the proposed algorithm shows slightly worse computing efficiency compared with Meng’s method and Sheen’s method, it outperforms its competitors in terms of quantitative indicators. Good performance has been observed by FEKO-simulated full-wave data from a three-dimensional target. The farther parts of the object were better reconstructed by the proposed method than by its competitors. Finally, millimeter-wave imaging using real data collected from our in-house-developed prototype demonstrated that the proposed algorithm performs well, as expected.

Author Contributions

Conceptualization, X.N., C.L., A.Q. and J.K.S.; methodology, X.N., C.L., A.Q. and I.D.R.; validation, X.N. and Y.M.; writing—original draft preparation, X.N.; writing—review and editing, X.N., C.L., A.Q., J.K.S. and I.D.R. All authors have read and agreed to the published version of the manuscript.

Funding

This work was supported in part by the Sichuan Science and Technology Program under Grant 2023NSFSC0463; National Young Thousand Talent under Grant A0920502051826, Grant YH199911041801, and Grant YX1199912371901; Foreign Talent in Culture and Education under Grant 110000207520190055; National Key Research and Development Plan under Grant 2018YFC0809500; Stable-Support Scientific Project of China Research Institute of Radio Wave Propagation under Grant A132003W02; and Fundamental Research Funds for the Central Universities under Grant 2682018CX20.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

Not applicable.

Acknowledgments

Xinyi Nie and Anyong Qing would like to thank Timothy Amsdon, School of Electronic and Electrical Engineering, University of Leeds, UK, for reviewing the manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

Appendix A

Here, the derivation of (6) by using MSP [36,41,45] is presented, which cannot be solved by the mathematical methods in classical algorithms, such as spherical wave decomposition and Weyl identity [11,25]. Let

\begin{array}{l} E (k_{x}, k_{z}, k) & = \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} \frac{\exp (- j 2 k R)}{R^{2}} e^{- j k_{x} x^{'}} e^{- j k_{z} z^{'}} d x^{'} d z^{'} \\ = \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} D (x^{'}, y^{'}) \cdot e^{j 2 k Φ (x^{'}, z^{'})} d x^{'} d z^{'} \end{array}

(A1)

where

D (x^{'}, z^{'}) = \frac{1}{{(x^{'} - x)}^{2} + {(y^{'} - y)}^{2} + {(z^{'} - z)}^{2}}

(A2)

\begin{array}{l} Φ (x^{'}, z^{'}) = & - \sqrt{{(x^{'} - x)}^{2} + {(y^{'} - y)}^{2} + {(z^{'} - z)}^{2}} \\ - \frac{k_{x}}{2 k} x^{'} - \frac{k_{z}}{2 k} z^{'} \end{array}

(A3)

When k is large and

Φ (x^{'}, z^{'})

reaches extreme value at a stationary point

(x_{s}, z_{s})

in the integrating range, the major contribution to the integral in (A1) comes from the integral value at the neighborhood near stationary point

(x_{s}, z_{s})

At the stationary point, we have

{\frac{\partial Φ}{\partial x^{'}}|}_{x_{s}, z_{s}} = \frac{- (x_{s} - x)}{\sqrt{{(x_{s} - x)}^{2} + {(y^{'} - y)}^{2} + {(z_{s} - z)}^{2}}} - \frac{k_{x}}{2 k} = 0

(A4)

{\frac{\partial Φ}{\partial z^{'}}|}_{x_{s}, z_{s}} = \frac{- (z_{s} - z)}{\sqrt{{(x_{s} - x)}^{2} + {(y^{'} - y)}^{2} + {(z_{s} - z)}^{2}}} - \frac{k_{z}}{2 k} = 0

(A5)

The stationary point can be solved as

x_{s} = x - \frac{k_{x}}{\sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}}} (y^{'} - y)

(A6)

z_{s} = z - \frac{k_{z}}{\sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}}} (y^{'} - y)

(A7)

Here, the negative roots are extracted for each case, and it has been assumed that

y^{'} - y < 0

The asymptotic approximation to (A1) is therefore

\begin{array}{l} E (k_{x}, k_{z}, k) & = D (x_{s}, z_{s}) e^{j 2 k Φ (x_{s}, z_{s})} \\ \cdot \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} e^{\frac{1}{2} j 2 k [α {(x_{s} - x)}^{2} + β {(z_{s} - z)}^{2} + 2 γ (x_{s} - x) (z_{s} - z)]} d x d z \\ = \frac{2 π j σ}{\sqrt{|α β - γ^{2}|}} D (x_{s}, z_{s}) \frac{e^{j 2 k Φ (x_{s}, z_{s})}}{2 k} \end{array}

(A8)

where

D (x_{s}, z_{s}) = \frac{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}}{{(2 k)}^{2} {(y - y^{'})}^{2}}

(A9)

α = {\frac{\partial^{2} Φ}{\partial {x^{'}}^{2}}|}_{x_{s}, z_{s}} = - \frac{{(2 k)}^{2} - k_{x}^{2}}{y^{'} - y} \cdot \frac{\sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}}}{{(2 k)}^{3}} > 0

(A10)

{β = \frac{\partial^{2} Φ}{\partial {z^{'}}^{2}}|}_{x_{s}, z_{s}} = - \frac{{(2 k)}^{2} - k_{z}^{2}}{y^{'} - y} \cdot \frac{\sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}}}{{(2 k)}^{3}}

(A11)

{γ = \frac{\partial^{2} Φ}{\partial x^{'} \partial z^{'}}|}_{x_{s}, z_{s}} = \frac{k_{x} k_{z} \sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}}}{(y^{'} - y) \cdot {(2 k)}^{3}}

(A12)

σ = \{\begin{cases} + 1 for α β > γ^{2}, α > 0 \\ - 1 for α β > γ^{2}, α < 0 \\ - j for α β < γ^{2} \end{cases}

(A13)

Φ (x_{s}, z_{s}) = - \frac{\sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}}}{2 k} (y - y^{'}) - \frac{k_{x}}{2 k} x - \frac{k_{z}}{2 k} z

(A14)

Eventually, substituting (A9)–(A14) into (A8), the asymptotic approximation expression is

E (k_{x}, k_{z}, k) = \frac{j π}{k} \cdot \frac{1}{y - y^{'}} \cdot e^{- j (\sqrt{{(2 k)}^{2} - k_{x}^{2} - k_{z}^{2}} (y - y^{'}) + k_{x} x + k_{z} z)}

(A15)

References

Global Terrorism Database. Available online: https://www.start.umd.edu/gtd/ (accessed on 3 February 2023).
Hounsfield, G.N. Method and Apparatus for Measuring X-or γ-Radiation Absorption or Transmission at Plural Angles and Analyzing the Data. US Patent 3778614, 11 December 1973. [Google Scholar]
Vollmer, M.; Klaus-Peter, M. Infrared Thermal Imaging: Fundamentals, Research and Applications; John Wiley & Sons: Weinheim, Germany, 2018. [Google Scholar]
Demirci, S.; Cetinkaya, H.; Yigit, E.; Ozdemir, C.; Vertiy, A.A. A study on millimeter-wave imaging of concealed objects: Application using back-projection algorithm. Prog. Electromagn. Res. 2012, 128, 457–477. [Google Scholar] [CrossRef]
Moll, J.; Schops, P.; Krozer, V. Towards three-dimensional millimeter-wave radar with the bistatic fast-factorized back-projection algorithm—Potential and limitations. IEEE Trans. Terahertz Sci. Technol. 2012, 2, 432–440. [Google Scholar] [CrossRef]
Abbasi, M.; Shayei, A.; Shabany, M.; Kavehvash, Z. Fast Fourier-based implementation of synthetic aperture radar algorithm for multistatic imaging system. IEEE Trans. Instrum. Meas. 2018, 68, 3339–3349. [Google Scholar] [CrossRef]
Broquetas, A.; Palau, J.; Jofre, L.; Cardama, A. Spherical wave near-field imaging and radar cross-section measurement. IEEE Trans. Antennas Propag. 1998, 46, 730–735. [Google Scholar] [CrossRef]
Boag, A. A fast multilevel domain decomposition algorithm for radar imaging. In Proceedings of the IEEE Antennas and Propagation Society International Symposium, Orlando, FL, USA, 11–16 July 1999; Volume 3, pp. 1788–1791. [Google Scholar]
Boag, A. A fast multilevel domain decomposition algorithm for radar imaging. IEEE Trans. Antennas Propag. 2001, 49, 666–671. [Google Scholar] [CrossRef]
Tan, K.; Chen, X. Precise near-range 3-D image reconstruction based on MIMO circular synthetic aperture radar. IEEE Trans. Microwave Theory Tech. 2021, 69, 2651–2661. [Google Scholar] [CrossRef]
Sheen, D.M.; McMakin, D.L.; Hall, T.E. Three-dimensional millimeter-wave imaging for concealed weapon detection. IEEE Trans. Microwave Theory Tech. 2001, 49, 1581–1592. [Google Scholar] [CrossRef]
Gabor, D. A new microscopic principle. Nature 1948, 161, 777–778. [Google Scholar] [CrossRef]
Wolf, E. Three-dimensional structure determination of semi-transparent objects from holographic data. Opt. Commun. 1969, 1, 153–156. [Google Scholar] [CrossRef]
Qiao, L.; Wang, Y.; Shen, Z.; Zhao, Z.; Chen, Z. Compressive sensing for direct millimeter-wave holographic imaging. Appl. Opt. 2015, 54, 3280–3289. [Google Scholar] [CrossRef]
Qiao, L.; Wang, Y.; Zhao, Z.; Chen, Z. Exact reconstruction for near-field threedimensional planar millimeter-wave holographic imaging. J. Infrared Millim Terahertz Waves. 2015, 36, 1221–1236. [Google Scholar] [CrossRef]
McMakin, D.L.; Sheen, D.M.; Collins, H.D.; Hall, T.E.; Severtsen, R.H. Wideband millimeter-wave holographic weapons surveillance systems. In Proceedings of the European Symposium on Optics for Environmental and Public Safety, Munich, Germany, 19–23 June 1995; Volume 2511, pp. 131–141. [Google Scholar]
Sheen, D.M.; McMakin, D.L.; Barber, J.; Hall, T.E.; Severtsen, R.H. Active imaging at 350 GHz for security applications. In Proceedings of the SPIE Defense and Security Symposium, Orlando, FL, USA, 16–20 March 2008; Volume 6948, p. 69480M. [Google Scholar]
Sheen, D.; McMakin, D.; Hall, T. Near-field three-dimensional radar imaging techniques and applications. Appl. Opt. 2010, 49, E83–E93. [Google Scholar] [CrossRef]
Sheen, D.M.; Clark, R.T.; Tedeschi, J.; Jones, A.M.; Hall, T.E. High-resolution 3D microwave imaging of a moving target using optical motion capture. In Proceedings of the SPIE Defense + Commercial Sensing, Baltimore, MD, USA, 14–18 April 2019; Volume 10994, p. 109940D. [Google Scholar]
Smith, D.; Yurduseven, O.; Livingstone, B.; Schejbal, V. Microwave imaging using indirect holographic techniques. IEEE Antennas Propag. Mag. 2014, 56, 104–117. [Google Scholar] [CrossRef]
Cheng, H.; Gu, J.; Cheng, Y.; Li, S.; Sun, H. Multi-section auto-focus millimeter-wave holography. In Proceedings of the 2016 IEEE MTT-S International Microwave Symposium (IMS), San Francisco, CA, USA, 22–27 May 2016; pp. 1–3. [Google Scholar]
Cheng, H.; Li, S.; Zheng, H.; Jing, H.; Sun, H. A W-band auto-focus holographic imaging system for security screening. IEICE Electron. Express. 2017, 14, 20170347. [Google Scholar] [CrossRef]
Moulder, W.F.; Krieger, J.D.; Majewski, J.J.; Coldwell, C.M.; Nguyen, H.T.; Maurais-Galejs, D.T.; Anderson, T.L.; Dufilie, P.; Herd, J.S. Development of a high-throughput microwave imaging system for concealed weapons detection. In Proceedings of the 2016 IEEE International Symposium on Phased Array Systems and Technology (PAST), Waltham, MA, USA, 18–21 October 2016. [Google Scholar]
Gao, J.; Deng, B.; Qin, Y.; Wang, H.; Li, X. An efficient algorithm for MIMO cylindrical millimeter-wave holographic 3-D imaging. IEEE Trans. Microwave Theory Tech. 2018, 66, 5065–5074. [Google Scholar] [CrossRef]
Meng, Y.; Lin, C.; Zang, J.; Qing, A.; Nikolova, N.K. General theory of holographic inversion with linear frequency modulation radar and its application to whole-body security scanning. IEEE Trans. Microwave Theory Tech. 2020, 68, 4694–4705. [Google Scholar] [CrossRef]
Meng, Y.; Lin, C.; Qing, A.; Nikolova, N.K. Accelerated holographic imaging with range stacking for linear frequency modulation radar. IEEE Trans. Microwave Theory Tech. 2022, 70, 1630–1638. [Google Scholar] [CrossRef]
Zhu, Y.-K.; Yang, M.-H.; Wu, L.; Sun, Y.; Sun, X.-W. Millimeter-wave holographic imaging algorithm with amplitude corrections. Prog. Electromagn. Res. M. 2016, 49, 33–39. [Google Scholar] [CrossRef]
Gao, J.; Deng, B.; Qin, Y.; Wang, H.; Li, X. Efficient terahertz wide-angle NUFFT-based inverse synthetic aperture imaging considering spherical wavefront. Sensors 2016, 16, 2120. [Google Scholar] [CrossRef]
Detlefsen, J.; Dallinger, A.; Huber, S.; Schelkshorn, S. Effective reconstruction approaches to millimeter-wave imaging of humans. In Proceedings of the 28th General Assembly of International Union of Radio Science, New Delhi, India, 23–29 October 2005; pp. 23–29. [Google Scholar]
Yanik, M.E.; Wang, D.; Torlak, M. Development and demonstration of MIMO-SAR mmWave imaging testbeds. IEEE Access 2020, 8, 126019–126038. [Google Scholar] [CrossRef]
Hughes, W.; Gault, K.; Princz, G. A comparison of the range-Doppler and chirp scaling algorithms with reference to RADARSAT. In Proceedings of the IGARSS ’96. 1996 International Geoscience and Remote Sensing Symposium, Lincoln, NE, USA, 31 May 1996; Volume 2, pp. 1221–1223. [Google Scholar]
Soumekh, M. Bistatic synthetic aperture radar inversion with application in dynamic object imaging. IEEE Trans. Signal Process. 1991, 39, 2044–2055. [Google Scholar] [CrossRef] [PubMed]
Bamler, R. A comparison of range-Doppler and wavenumber domain SAR focusing algorithms. IEEE Trans. Geosci. Remote Sens. 1992, 30, 706–713. [Google Scholar] [CrossRef]
Raney, R.K.; Runge, H.; Bamler, R.; Cumming, I.G.; Wong, F.H. Precision SAR processing using chirp scaling. IEEE Trans. Geosci. Remote Sens. 1994, 32, 786–799. [Google Scholar] [CrossRef]
Zhuge, X.; Yarovoy, A.G. Three-dimensional near-field MIMO array imaging using range migration techniques. IEEE Trans. Image Process. 2012, 21, 3026–3033. [Google Scholar] [CrossRef] [PubMed]
Lopez-Sanchez, J.M.; Fortuny-Guasch, J. 3-D radar imaging using range migration techniques. IEEE Trans. Antennas Propag. 2000, 48, 728–737. [Google Scholar] [CrossRef]
Tan, K.; Wu, S.; Liu, X.; Fang, G. Omega-K algorithm for near-field 3-D image reconstruction based on planar SIMO/MIMO array. IEEE Trans. Geosci. Remote Sens. 2019, 57, 2381–2394. [Google Scholar] [CrossRef]
Tan, W.; Huang, P.; Huang, Z.; Qi, Y.; Wang, W. Three-dimensional microwave imaging for concealed weapon detection using range stacking technique. Int. J. Antennas Propag. 2017, 2017, 1480623. [Google Scholar] [CrossRef]
Zhuge, X.; Yarovoy, A.G.; Savelyev, T.; Ligthart, L. Modified Kirchhoff migration for UWB MIMO array-based radar imaging. IEEE Trans. Geosci. Remote Sens. 2010, 48, 2692–2703. [Google Scholar] [CrossRef]
Wang, Z.; Guo, Q.; Tian, X.; Chang, T.; Cui, H. Millimeter-wave image reconstruction algorithm for one-stationary bistatic SAR. IEEE Trans. Microwave Theory Tech. 2020, 68, 1185–1194. [Google Scholar] [CrossRef]
Fromenteze, T.; Yurduseven, O.; Berland, F.; Decroze, C.; Smith, D.R.; Yarovoy, A.G. A transverse spectrum deconvolution technique for MIMO short-range Fourier imaging. IEEE Trans. Geosci. Remote Sens. 2019, 57, 6311–6324. [Google Scholar] [CrossRef]
Wang, J.; Aubry, P.; Yarovoy, A. 3-D short-range imaging with irregular MIMO arrays using NUFFT-based range migration algorithm. IEEE Trans. Geosci. Remote Sens. 2020, 58, 4730–4742. [Google Scholar] [CrossRef]
Fromentèze, T.; Yurduseven, O.; Berland, F.; Cyril, D.; Smith, D.R.; Yarovoy, A.G. Accelerating short range MIMO imaging with optimized Fourier processing. In Proceedings of the SPIE Defense + Commercial Sensing, Online, 27 April–9 May 2020; Volume 11411D, pp. 56–61. [Google Scholar]
Nie, X.; Cheng, L.; Wang, S.; Amsdon, T.; Roberston, I.D.; Qing, A. The effect of range on close-Range millimeter-wave imaging systems. In Proceedings of the 2023 IEEE International Symposium on Antennas and Propagation and North American Radio Science Meeting, Portland, OR, USA, July 2023. (accepted). [Google Scholar]
Born, M.; Wolf, E. Principles of Optics: Electromagnetic Theory of Propagation, Interference and Diffraction of Light, 15th ed.; Cambridge University Press: Cambridge, UK, 2019. [Google Scholar]
Nikolova, N.K. Introduction to Microwave Imaging; Cambridge University Press: New York, NY, USA, 2017. [Google Scholar]
Kettler, D.; Redding, N.; Gray, D. Exact point spread functions for SAR imaging. In Proceedings of the Committee Earth Observation Satellites, Synthetic Aperture Radar Workshop, September 2005; pp. 1–8. Available online: https://www.researchgate.net/publication/310651741_Exact_Point_Spread_Functions_for_SAR_Imaging (accessed on 6 March 2023).
Ahmed, S.S.; Schiessl, A.; Schmidt, L. A novel active real-time digital-beamforming imager for personnel screening. In Proceedings of the EUSAR 2012; 9th European Conference on Synthetic Aperture Radar, Nuremberg, Germany, 23–26 April 2012; pp. 178–181. [Google Scholar]
Zhou, W.; Bovik, A.C.; Sheikh, H.R.; Simoncelli, E.P. Image quality assessment: From error visibility to structural similarity. IEEE Trans. Image Process. 2004, 13, 600–612. [Google Scholar]

Figure 1. The geometry of imaging system.

Figure 2. The reconstruction results of point scatterer in scale normalized to the maximum. (a) The reflectivity contour projected to xz plane. (b) The 3-D view of the reflectivity projected to the xz plane.

Figure 3. Profiles along x(z) in dB scale normalized to the maximum.

Figure 4. Profiles along x = z in dB scale normalized to the maximum.

Figure 5. Profiles along y in dB scale normalized to the maximum.

Figure 6. The 2-D object under test.

Figure 7. The reconstruction results. (a) Sheen’s method. (b) Meng’s method. (c) The proposed method.

Figure 8. The reconstructed and true normalized magnitude along the line z = 0 m.

Figure 9. The CAD model of three blocks within different ranges.

Figure 10. The front view of the reconstruction results. (a) Sheen’s method. (b) Meng’s method. (c) The proposed method.

Figure 11. The side view of the reconstruction results. (a) Sheen’s method. (b) Meng’s method. (c) The proposed method.

Figure 12. The top view of the reconstruction results. (a) Sheen’s method. (b) Meng’s method. (c) The proposed method.

Figure 13. Profiles along z = 0 m in Figure 10.

Figure 14. Profiles along z = 0 m in Figure 11.

Figure 15. Photograph of the target under test.

Figure 16. Zoomed-in reconstructed result.

Figure 17. Profiles along the red horizontal line in Figure 16.

Figure 18. Profiles along the red vertical line in Figure 16.

Figure 19. The scenario with a volunteer scanned by prototype.

Figure 20. The reconstruction results (a) by Meng’s method and (b) by the proposed method.

Table 1. Simulation parameters adopted in synthetic system.

Parameters	Value
Center frequency	29.9 GHz
Frequency bandwidth	5.8 GHz
The number of sampling frequency points	220
The length of the aperture along x axis	360 mm
The sampling interval along x axis	5 mm
The length of the aperture along z axis	360 mm
The sampling interval along z axis	5 mm

Table 2. Comparison of SSIM and RMSE for synthesis data.

	Sheen’s Method [11]	Meng’s Method [25]	The Proposed Method
SSIM	9.495 × 10⁻¹	9.501 × 10⁻¹	9.544 × 10⁻¹
RMSE	5.957 × 10⁻³	5.921 × 10⁻³	4.316 × 10⁻³

Table 3. Computation time for synthesis data.

Method	Time (s)
Sheen’s method	3.44
Meng’s method	3.54
The proposed method	3.61

Table 4. Comparison of SSIM and RMSE.

	Sheen’s Method [11]	Meng’s Method [25]	The Proposed Method
SSIM	9.2118 × 10⁻¹	9.218 × 10⁻¹	9.270 × 10⁻¹
RMSE	1.0753 × 10⁻³	1.0714 × 10⁻³	8.774 × 10⁻³

Table 5. Computation time for simulation.

Algorithm	Time (s)
Sheen’s method	2.86
Meng’s method	2.93
The proposed method	2.94

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

© 2023 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Nie, X.; Lin, C.; Meng, Y.; Qing, A.; Sykulski, J.K.; Robertson, I.D. An Accurate Millimeter-Wave Imaging Algorithm for Close-Range Monostatic System. Sensors 2023, 23, 4577. https://doi.org/10.3390/s23104577

AMA Style

Nie X, Lin C, Meng Y, Qing A, Sykulski JK, Robertson ID. An Accurate Millimeter-Wave Imaging Algorithm for Close-Range Monostatic System. Sensors. 2023; 23(10):4577. https://doi.org/10.3390/s23104577

Chicago/Turabian Style

Nie, Xinyi, Chuan Lin, Yang Meng, Anyong Qing, Jan K. Sykulski, and Ian D. Robertson. 2023. "An Accurate Millimeter-Wave Imaging Algorithm for Close-Range Monostatic System" Sensors 23, no. 10: 4577. https://doi.org/10.3390/s23104577

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu

An Accurate Millimeter-Wave Imaging Algorithm for Close-Range Monostatic System

Abstract

1. Introduction

2. Theory

3. Results and Discussions

3.1. Point Spread Function

3.2. The Synthetic Data

3.3. The Full-Wave Data

3.4. Experimental Results

4. Conclusions

Author Contributions

Funding

Institutional Review Board Statement

Informed Consent Statement

Data Availability Statement

Acknowledgments

Conflicts of Interest

Appendix A

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI