# Three-Dimensional Microwave Imaging for Concealed Weapon Detection Using Range Stacking Technique

Authors: Weixian Tan; Pingping Huang; Zengshu Huang; Yaolong Qi; Wenqin Wang
Journal: International Journal of Antennas and Propagation
 (2017)

---

## Abstract

Three-dimensional (3D) microwave imaging has been proven to be well suited for concealed weapon detection application. For the 3D image reconstruction under two-dimensional (2D) planar aperture condition, most of current imaging algorithms focus on decomposing the 3D free space Green function by exploiting the stationary phase and, consequently, the accuracy of the final imagery is obtained at a sacrifice of computational complexity due to the need of interpolation. In this paper, from an alternative viewpoint, we propose a novel interpolation-free imaging algorithm based on wavefront reconstruction theory. The algorithm is an extension of the 2D range stacking algorithm (RSA) with the advantages of low computational cost and high precision. The algorithm uses different reference signal spectrums at different range bins and then forms the target functions at desired range bin by a concise coherent summation. Several practical issues such as the propagation loss compensation, wavefront reconstruction, and aliasing mitigating are also considered. The sampling criterion and the achievable resolutions for the proposed algorithm are also derived. Finally, the proposed method is validated through extensive computer simulations and real-field experiments. The results show that accurate 3D image can be generated at a very high speed by utilizing the proposed algorithm.


---

## Body

## 1. Introduction

Recently, increasing terrors have been motivating the development of various concealed weapon detection techniques for better security monitoring at airports, courthouses, and other buildings [1–15]. The ability to detect dangerous hidden items is seen as a significant tool in defeating terrorism in today’s environment. However, it is monumentally difficult to timely detect concealed weapons carried on an individual with conventional metal detectors. Since microwave possesses a unique property of passing transparently through materials such as heavy clothing, it is possible to identify the presence of hidden threatening objects under clothing with microwave imaging technique. In addition, microwave imaging systems eliminate the use of ionizing radiation such as those present in X-ray systems and, therefore, pose no known health hazard at a low power. In the past few years, an increased interest has been addressed to three-dimensional (3D) microwave imaging [16–22]. The system has the capability of mapping the microwave reflectivity of an individual in 3D space and providing the basis for detection, identification, and further automatic recognition. Lately, the Smith Detection Inc. in England developed a flat-panel millimeter-wave technique for concealed weapon detection for the increasing people screening market [4]. Two-dimensional (2D) planar scanning geometries are adopted in all the above systems and, consequently, a 3D reflectivity image of the human body can be formed.As for the personnel surveillance applications, most of the previous image reconstruction algorithms based on 2D planar synthetic aperture mainly focus on the 3D range migration algorithm (RMA) [23], which requires a 3D STOLT interpolation before the 3D inverse discrete Fourier transform (DFT). The formulation of the 3D RMA is derived by using the stationary phase [24] or the Fourier decomposition of the 3D free space Green function. Here, we discuss yet another interpolation-free reconstruction algorithm which accounts precisely for not only the free space propagation losses but also the wavefront curvature. The algorithm, known as range stacking algorithm (RSA), is usually derived based on wave equations and has been widely and successfully used in seismic signal processing for oil exploration for decades [25]. To date, it has only been used in its 2D version in the field of synthetic aperture radar (SAR) imaging. The discussed algorithm is an extension of the work of Soumekh, who primarily developed the 2D RSA for the conventional SAR imaging [26]. Here we introduce an alternative viewpoint to reconstruct the image and derive the interpolation-free reconstruction algorithm from wavenumber domain instead of wave equations. This algorithm uses a different matched filtering reference for each range bin and then numerically integrates over the fast-time frequencies before inverse Fourier transform to get the final result. By explicitly performing these steps, the interpolation step required in existing algorithms is not required and, moreover, the final image does not suffer from the truncation and DFT wrap-around errors in the wavenumber domain. The algorithm has been successfully extended to 3D microwave imaging for concealed weapon detection and the reflectivity images at several range bins within the illuminated range swath can be easily formed.The remaining of the paper is organized as follows. An introduction to 3D microwave imaging geometry and an alternative viewpoint for image reconstruction under 2D planar aperture condition are described in Section2. The extended 3D RSA for implementing the image reconstruction is presented, and the corresponding procedures are addressed in Section 3. Section 4 is dedicated to numerical simulations and experimental results with realistic measurement data. Finally, the conclusions are provided in Section 5.


## 2. Formulation of the Problem

### 2.1. Planar Scanning Geometry Signal Model and Image Reconstruction

We begin by considering a 3D imaging geometry for the detection of concealed weapons or other hidden items carried on personnel, as shown in Figure1. A wide-band stepped-frequency (or frequency-modulated) continuous-wave (CW) signal is radiated with a linear array of antennas mounted on a fast mechanical scanner along the y-axis or z-axis. When the scanner is moving along the corresponding axis, the sequentially switched transmitter receivers scan quickly over a large linear aperture to actively illuminate the target, so that a full 2D planar rectangular aperture of data can be gathered. The equivalent sample point of the quasi-monostatic transmit and receive antenna pair is assumed to be at the position (0,u,v), and the target is assumed to be at the position (x,y,z) with a spatial distribution of reflectivity function f(x,y,z). Taking the free space propagation losses into consideration, the signal reflected from a given target Pn at position (x,y,z) can be written as follows when measured at position (0,u,v) [2, 14, 17]:(1)s0Kw,u,v∣x,y,z=fx,y,zexp⁡-j2KwRR2,where Kw∈[Kwmin,Kwmax] is the fast-time frequency wavenumber and related to the temporal angular frequency w and speed of light C via Kw=w/c; Kwmin and Kwmax are the wavenumbers at the minimum and maximum angular frequencies, respectively; R is the distance between the target Pn and the measurement point (0,u,v); that is,(2)R=Xc+x2+y-u2+z-v2,where x∈[-X0/2,X0/2], Xc is the reference range, and X0 is the range swath along the x-axis.Figure 1
Planar scanning geometries for concealed weapon detection.LetV denote the support band of the 3D illuminated area, and then the response measured at the planar aperture transceiver positions at frequency wavenumber Kw can be formulated as follows:(3)sKw,u,v=∭Vs0Kw,u,v∣x,y,zdr=∭Vfx,y,zexp⁡-j2KwRR2dr,where dr=dxdydz and |y-u|≤LY/2 and |z-v|≤LZ/2 with LY and LZ representing the aperture length along the y-axis and z-axis, respectively. Note that either of them subtends the corresponding half-power beamwidth in that direction.The signal model (3) indicates that the measured response comprises echoes of different targets, any one of which has an amplitude inversely proportional to the square of the distance to the target. This amplitude, which decays with the range, is referred to as propagation loss. If left unhandled, it will degrade the image focusing quality. However, according to (3), the propagation losses are difficult to be directly removed in the frequency wavenumber domain. In this paper, we compensate by first transforming s(Kw,u,v) to the slant range spatial domain and then multiplying with an amplitude factor r2; that is,(4)s′Kw,u,v=∫r∫KwsKw,u,vexp⁡jKwrdKw×r2exp⁡-jKwrdr,and, except for some additional computations, (4) brings about the benefit of being able to reduce the coupling between transmitter and receiver and undesired noise beyond the region of interest by a simple range gating in the spatial domain. The computational load and the data volume can be further reduced by subsampling s(Kw,u,v) without aliasing, and then the useful signal reflected from the human body could be picked up for subsequent processing at the subsampled set.After the propagation loss compensation, the 3D microwave reflectivity image of the human body can be obtained by an integral resembling the inverse transform of (4) with respect to the variables Kw, u, and v; namely,(5)f^x,y,z=12π3∫∫S∫Kws′Kw,u,v×exp⁡j2KwRdKwdudv,where S denotes the surface of the planar synthetic aperture and f^(x,y,z) represents the reconstructed values of the ideal target function f(x,y,z). The image reconstruction can be derived by using the stationary phase. For instance, both Sheen and Soumehk decomposed the spherical-wave exponential term in (5) into a superposition of plane-wave components with Fourier transform. To facilitate the understanding of the image reconstruction, we give an alternative derivation instead of that of (5) in the following subsection.


### 2.2. An Alternative Viewpoint for Image Reconstruction

Firstly, let us consider the target function reconstruction in the 3D Fourier domain:(6)f^x,y,z=∭KFKx,Ky,Kz×exp⁡jKxx+jKyy+jKzzdKxdKydKz,where Kx, Ky, and Kz are the frequency wavenumber variables along the x-, y-, and z-axis, respectively; K denotes the support band of this 3D wavenumber domain. In order to reconstruct the 3D reflectivity image of the illuminated target area, the relationship between F(Kx,Ky,Kz) and the measured response s(Kw,u,v) can be found as follows:(7)FKxKw,Ku,Kv,KyKw,Ku,Kv,KzKw,Ku,Kv=SKw,Ku,KvHKw,Ku,Kv,where S(Kw,Ku,Kv) is the Fourier transform of the propagation loss compensated signal s′(Kw,u,v) with respect to (u,v) and H(Kw,Ku,Kv) is the matched filter for correcting the wavefront curvature of all target points with the same reference range Xc, which has the following form:(8)HKw,Ku,Kv=exp⁡j4Kw2-Ku2-Kv2Xc.When the 3D RMA is used, the wavefront curvature errors can be compensated by the known STOLT transformation; namely,(9)Kx=4Kw2-Ku2-Kv2,Ky=Ku,Kz=Kv,where computationally intensive 3D interpolation is required to recover the rectilinear samples of S(Kx,Ky,Kz) from available S(Kw,Ku,Kv) data samples. To avoid this, substituting (7)~(9) into (6) and making a change of variables from dKxdKydKz to dKwdKudKv yield(10)f^x,y,z=∫Kv∫Ku∫KwJKw,Ku,KvSKw,Ku,Kv×exp⁡j4Kw2-Ku2-Kv2Xc×exp⁡j4Kw2-Ku2-Kv2x+jKuy+jKvzdKwdKudKv,where JKw,Ku,Kv is the Jacobian determinant to accommodate for the change of coordinates which can be formulated as(11)JKw,Ku,Kv=4Kw4Kw2-Ku2-Kv2.The Jacobian determinant is a slowly fluctuating function and as a result its impact on the 3D image reconstruction is negligible. In the following discussions, we assume that the Jacobian function, or any other amplitude factor, is absorbed in the measured response spectrumS(Kw,Ku,Kv). Suppose also that we are interested in forming the 3D image at range bins xns, where xn∈[-X0/2,X0/2]. Consider the above reconstruction for f(x,y,z) at a fixed range x=xn:(12)f^x=xn,y,z=∫Kv∫Ku∫KwSKw,Ku,Kv×exp⁡jKuy+jKvz×exp⁡j4Kw2-Ku2-Kv2Xc+xndKwdKudKv.It can be derived that the second exponential term in (12) is related to a reference signal sxnKw,u,v (also propagation loss compensated) at the reference range bin x=xn via the following expression [25]:(13)exp⁡j4Kw2-Ku2-Kv2Xc+xn=∫v∫uexp⁡-j2KwXc+xn2+u2+v2exp⁡jKuu+jKvvdudv∗=∫v∫usxnKw,u,vexp⁡jKuu+jKvvdudv∗=Sxn∗Kw,Ku,Kv,where the reference signal sxnKw,u,v is the response from a reflector located at (x,y,z)=(xn,0,0) with Sxn(Kw,Ku,Kv) representing its 2D Fourier transform with respect to (u,v) and ∗ denotes the complex conjugate operator. Then the target function reconstruction at (x=xn,y,z) can be rewritten as(14)f^x=xn,y,z=∫Kv∫Ku∫KwSKw,Ku,Kv×Sxn∗Kw,Ku,Kvexp⁡jKuy+jKvzdKwdKudKv=∫Kv∫Ku∫KwSKw,Ku,KvSxn∗Kw,Ku,KvdKw×expjKuy+jKvzdKudKv.(14) indicates that the spatial distribution of target reflectivity at a specific range bin x=xn could be obtained by multiplying the measured response with the reference signal at the same reference range bin xn, followed by integrating the result over the available fast-time frequency wavenumbers Kw, as well as a 2D inverse Fourier transform with respect to (Ku,Kv), rather than mapping the (Kw,Ku,Kv) domain into the (Kx,Ky,Kz) domain through a computationally intensive interpolation as that of (9). It is worth noting that, before transforming to s(Kw,Ku,Kv), the propagation loss compensated signal s′(Kw,u,v) is often upsampled along the Y/U- and Z/V-axes through zero-padding to avoid any possible aliasing. Thereafter, a parallel implementation of the above procedures can simultaneously form the target function at the individual range bins xns and then stack the outcomes together to form a 3D image of the target function. Next, we develop the detailed procedures of the 3D range stacking algorithm (3D RSA) for implementing the above reconstructions.


## 2.1. Planar Scanning Geometry Signal Model and Image Reconstruction

We begin by considering a 3D imaging geometry for the detection of concealed weapons or other hidden items carried on personnel, as shown in Figure1. A wide-band stepped-frequency (or frequency-modulated) continuous-wave (CW) signal is radiated with a linear array of antennas mounted on a fast mechanical scanner along the y-axis or z-axis. When the scanner is moving along the corresponding axis, the sequentially switched transmitter receivers scan quickly over a large linear aperture to actively illuminate the target, so that a full 2D planar rectangular aperture of data can be gathered. The equivalent sample point of the quasi-monostatic transmit and receive antenna pair is assumed to be at the position (0,u,v), and the target is assumed to be at the position (x,y,z) with a spatial distribution of reflectivity function f(x,y,z). Taking the free space propagation losses into consideration, the signal reflected from a given target Pn at position (x,y,z) can be written as follows when measured at position (0,u,v) [2, 14, 17]:(1)s0Kw,u,v∣x,y,z=fx,y,zexp⁡-j2KwRR2,where Kw∈[Kwmin,Kwmax] is the fast-time frequency wavenumber and related to the temporal angular frequency w and speed of light C via Kw=w/c; Kwmin and Kwmax are the wavenumbers at the minimum and maximum angular frequencies, respectively; R is the distance between the target Pn and the measurement point (0,u,v); that is,(2)R=Xc+x2+y-u2+z-v2,where x∈[-X0/2,X0/2], Xc is the reference range, and X0 is the range swath along the x-axis.Figure 1
Planar scanning geometries for concealed weapon detection.LetV denote the support band of the 3D illuminated area, and then the response measured at the planar aperture transceiver positions at frequency wavenumber Kw can be formulated as follows:(3)sKw,u,v=∭Vs0Kw,u,v∣x,y,zdr=∭Vfx,y,zexp⁡-j2KwRR2dr,where dr=dxdydz and |y-u|≤LY/2 and |z-v|≤LZ/2 with LY and LZ representing the aperture length along the y-axis and z-axis, respectively. Note that either of them subtends the corresponding half-power beamwidth in that direction.The signal model (3) indicates that the measured response comprises echoes of different targets, any one of which has an amplitude inversely proportional to the square of the distance to the target. This amplitude, which decays with the range, is referred to as propagation loss. If left unhandled, it will degrade the image focusing quality. However, according to (3), the propagation losses are difficult to be directly removed in the frequency wavenumber domain. In this paper, we compensate by first transforming s(Kw,u,v) to the slant range spatial domain and then multiplying with an amplitude factor r2; that is,(4)s′Kw,u,v=∫r∫KwsKw,u,vexp⁡jKwrdKw×r2exp⁡-jKwrdr,and, except for some additional computations, (4) brings about the benefit of being able to reduce the coupling between transmitter and receiver and undesired noise beyond the region of interest by a simple range gating in the spatial domain. The computational load and the data volume can be further reduced by subsampling s(Kw,u,v) without aliasing, and then the useful signal reflected from the human body could be picked up for subsequent processing at the subsampled set.After the propagation loss compensation, the 3D microwave reflectivity image of the human body can be obtained by an integral resembling the inverse transform of (4) with respect to the variables Kw, u, and v; namely,(5)f^x,y,z=12π3∫∫S∫Kws′Kw,u,v×exp⁡j2KwRdKwdudv,where S denotes the surface of the planar synthetic aperture and f^(x,y,z) represents the reconstructed values of the ideal target function f(x,y,z). The image reconstruction can be derived by using the stationary phase. For instance, both Sheen and Soumehk decomposed the spherical-wave exponential term in (5) into a superposition of plane-wave components with Fourier transform. To facilitate the understanding of the image reconstruction, we give an alternative derivation instead of that of (5) in the following subsection.


## 2.2. An Alternative Viewpoint for Image Reconstruction

Firstly, let us consider the target function reconstruction in the 3D Fourier domain:(6)f^x,y,z=∭KFKx,Ky,Kz×exp⁡jKxx+jKyy+jKzzdKxdKydKz,where Kx, Ky, and Kz are the frequency wavenumber variables along the x-, y-, and z-axis, respectively; K denotes the support band of this 3D wavenumber domain. In order to reconstruct the 3D reflectivity image of the illuminated target area, the relationship between F(Kx,Ky,Kz) and the measured response s(Kw,u,v) can be found as follows:(7)FKxKw,Ku,Kv,KyKw,Ku,Kv,KzKw,Ku,Kv=SKw,Ku,KvHKw,Ku,Kv,where S(Kw,Ku,Kv) is the Fourier transform of the propagation loss compensated signal s′(Kw,u,v) with respect to (u,v) and H(Kw,Ku,Kv) is the matched filter for correcting the wavefront curvature of all target points with the same reference range Xc, which has the following form:(8)HKw,Ku,Kv=exp⁡j4Kw2-Ku2-Kv2Xc.When the 3D RMA is used, the wavefront curvature errors can be compensated by the known STOLT transformation; namely,(9)Kx=4Kw2-Ku2-Kv2,Ky=Ku,Kz=Kv,where computationally intensive 3D interpolation is required to recover the rectilinear samples of S(Kx,Ky,Kz) from available S(Kw,Ku,Kv) data samples. To avoid this, substituting (7)~(9) into (6) and making a change of variables from dKxdKydKz to dKwdKudKv yield(10)f^x,y,z=∫Kv∫Ku∫KwJKw,Ku,KvSKw,Ku,Kv×exp⁡j4Kw2-Ku2-Kv2Xc×exp⁡j4Kw2-Ku2-Kv2x+jKuy+jKvzdKwdKudKv,where JKw,Ku,Kv is the Jacobian determinant to accommodate for the change of coordinates which can be formulated as(11)JKw,Ku,Kv=4Kw4Kw2-Ku2-Kv2.The Jacobian determinant is a slowly fluctuating function and as a result its impact on the 3D image reconstruction is negligible. In the following discussions, we assume that the Jacobian function, or any other amplitude factor, is absorbed in the measured response spectrumS(Kw,Ku,Kv). Suppose also that we are interested in forming the 3D image at range bins xns, where xn∈[-X0/2,X0/2]. Consider the above reconstruction for f(x,y,z) at a fixed range x=xn:(12)f^x=xn,y,z=∫Kv∫Ku∫KwSKw,Ku,Kv×exp⁡jKuy+jKvz×exp⁡j4Kw2-Ku2-Kv2Xc+xndKwdKudKv.It can be derived that the second exponential term in (12) is related to a reference signal sxnKw,u,v (also propagation loss compensated) at the reference range bin x=xn via the following expression [25]:(13)exp⁡j4Kw2-Ku2-Kv2Xc+xn=∫v∫uexp⁡-j2KwXc+xn2+u2+v2exp⁡jKuu+jKvvdudv∗=∫v∫usxnKw,u,vexp⁡jKuu+jKvvdudv∗=Sxn∗Kw,Ku,Kv,where the reference signal sxnKw,u,v is the response from a reflector located at (x,y,z)=(xn,0,0) with Sxn(Kw,Ku,Kv) representing its 2D Fourier transform with respect to (u,v) and ∗ denotes the complex conjugate operator. Then the target function reconstruction at (x=xn,y,z) can be rewritten as(14)f^x=xn,y,z=∫Kv∫Ku∫KwSKw,Ku,Kv×Sxn∗Kw,Ku,Kvexp⁡jKuy+jKvzdKwdKudKv=∫Kv∫Ku∫KwSKw,Ku,KvSxn∗Kw,Ku,KvdKw×expjKuy+jKvzdKudKv.(14) indicates that the spatial distribution of target reflectivity at a specific range bin x=xn could be obtained by multiplying the measured response with the reference signal at the same reference range bin xn, followed by integrating the result over the available fast-time frequency wavenumbers Kw, as well as a 2D inverse Fourier transform with respect to (Ku,Kv), rather than mapping the (Kw,Ku,Kv) domain into the (Kx,Ky,Kz) domain through a computationally intensive interpolation as that of (9). It is worth noting that, before transforming to s(Kw,Ku,Kv), the propagation loss compensated signal s′(Kw,u,v) is often upsampled along the Y/U- and Z/V-axes through zero-padding to avoid any possible aliasing. Thereafter, a parallel implementation of the above procedures can simultaneously form the target function at the individual range bins xns and then stack the outcomes together to form a 3D image of the target function. Next, we develop the detailed procedures of the 3D range stacking algorithm (3D RSA) for implementing the above reconstructions.


## 3. 3D Range Stacking Algorithm

### 3.1. Reconstruction Procedure

LetFy,z(x,Ky,Kz) denote the marginal Fourier transform of the reconstructed target function f^(x,y,z) with respect to the variables (y,z), and then the target function at range x=xn can be recovered by a 2D inverse Fourier transform of Fy,z(x=xn,Ky,Kz) with respect to (Ky,Kz), which is(15)f^x=xn,y,z=∫Kz∫KyFy,zx=xn,Ky,Kz×exp⁡jKyy+jKzzdKydKz,where the constant amplitude factor in the Fourier integral is neglected. Comparing (15) with (14), we can get that the reconstruction for Fy,z(x=xn,Ky,Kz) is(16)Fy,zx=xn,Ky,Kz=∫KwSKw,Ku,Kv×Sxn∗Kw,Ku,KvdKw,where(17)Ky=Ku,Kz=Kv.Practically, the above one-dimensional integral can be implemented by summing over the available discrete values of Kw and the digital implementation of the Fourier transform or its inverse; that is, the forward/inverse DFT can be realized through efficient fast Fourier transform (FFT) and inverse FFT (IFFT) algorithms. Thus, according to (14)~(17), the discrete form of the target function reconstruction at (x=xn,y,z) is expressed as(18)f^x=xn,y,z=IFFTKy=KuKz=Kv∑KwFFTu,vSKw,u,vSxn∗Kw,Ku,KvΔKw.In summary, the reconstruction procedure of the 3DRSA is illustrated in Figure2, which contains the following steps.Figure 2
Flowchart of 3D range stacking algorithm for concealed weapon detection.Step 1.
Perform an IFFT on the measured responseS(Kw,u,v) with respect to the fast-time frequency wavenumbers Kw. Compensate the propagation losses, and reduce the coupling and noises beyond the region of interest by range gating, and then the compensated signal is upsampled along the y- and z-axes through zero-padding to avoid any aliasing, followed by transforming to the 3D Fourier domain via a 3D FFT.Step 2.
The resultant 3D Fourier domain signal from Step1 is first matched filtered with the reference signal spectrum at the range bin xn, Sxn∗(Kw,Ku,Kv). The result is summed over the available fast-time frequency wavenumbers Kw to yield the marginal Fourier transform of the target function Fy,z(x=xn,Ky,Kz) at the range xn with Ky=Ku and Kz=Kv. Then a 2D IFFT with respect to (Ky,Kz) is performed on Fy,z(x=xn,Ky,Kz) to obtain the target reconstruction f^(x=xn,y,z) at the range bin xn. For all available range bins xns, the above calculations are implemented.Step 3.
Stack all outcomes of Step2 together to form a 3D image of the target function, f^(x,y,z); then apply possible postprocessing procedures such as smoothing filtering, geometry regridding, and object recognition and classification [27, 28] and finally display the image.


### 3.2. Sampling Criteria and Resolutions

According to the Nyquist sampling criterion [29], the sampling interval along some specific direction is inversely proportional to the spatial frequency domain support band of the signal in that direction.As for the wave propagation direction, the frequency domain support band of the transmitted signal is determined by(19)ΩKw=2Kwmin,2Kwmax.Thus the following Nyquist sampling constraint should be met in order to avoid aliasing along the slant range direction:(20)RmaxKwmax-Kwmin/ΔKw≤2π2Kwmax-Kwmin⟹ΔKw≤πRmax⟹Δw≤CπRmax,where Δw is the frequency sampling interval and Rmax denotes the maximum observed distance within the antenna beam pattern.As for the sampling criterion along theY/U direction, the samples are realized via a real or synthetic antenna array along this direction, and the frequency wavenumber support band ΩKy depends on the antenna beamwidth of the individual antenna ΦY through(21)ΩKy=-2Kwsin⁡ΦY2,2Kwsin⁡ΦY2.Thus the Nyquist sampling spacing Δu along the Y/U axis for gathering the echoed signal must satisfy the following inequality:(22)Δu≤2π4Kwsin⁡ΦY/2⟹Δu≤2π4Kwmaxsin⁡ΦY/2=λmin4sin⁡ΦY/2,where λmin is the smallest wavelength at the highest frequency Kwmax and the antenna beamwidth ΦY is determined by the antenna size along the corresponding direction. The worst-case situation corresponds to the case of ΦY = 180° for which the Nyquist sampling spacing in the Y/U dimension is(23)Δu≤λmin4.As for the case along theZ/V direction, a similar conclusion can be drawn that the Nyquist sampling criterion in the Z/V dimension takes the form (24)Δv≤λmin4sin⁡Φz/2,where Φz denotes the antenna beamwidth along the Z/V direction.The frequency wavenumber support band and the derived sampling constraints depend on the dimensions of the volume occupied in the 3D wavenumber domain. Therefore, the spatial resolutions of the resulting 3D reflectivity imagef^(x,y,z) are determined by the signal bandwidth, the center frequency, and the effective synthetic aperture. The achievable −3 dB resolutions along the x-, y-, and z-axes, are approximated, respectively, by(25)δx≈0.886πKwmax-Kwmin,(26)δy≈0.886λc4sin⁡ΦY/2,(27)δz≈0.886λc4sin⁡ΦZ/2,where λc is the wavelength at the center frequency fc and the units of δx, δy, and δz are all meters. In practice, due to the sidelobe structure of the radiation patterns, the final resolutions are usually slightly poorer than those theoretical values given by (25)~(27). In spite of this, both simulation and experimental results provided in the next section show that the 3D RSA is well suited to reconstruct 3D microwave imaging in concealed weapon detection at a high-precision manner.


## 3.1. Reconstruction Procedure

LetFy,z(x,Ky,Kz) denote the marginal Fourier transform of the reconstructed target function f^(x,y,z) with respect to the variables (y,z), and then the target function at range x=xn can be recovered by a 2D inverse Fourier transform of Fy,z(x=xn,Ky,Kz) with respect to (Ky,Kz), which is(15)f^x=xn,y,z=∫Kz∫KyFy,zx=xn,Ky,Kz×exp⁡jKyy+jKzzdKydKz,where the constant amplitude factor in the Fourier integral is neglected. Comparing (15) with (14), we can get that the reconstruction for Fy,z(x=xn,Ky,Kz) is(16)Fy,zx=xn,Ky,Kz=∫KwSKw,Ku,Kv×Sxn∗Kw,Ku,KvdKw,where(17)Ky=Ku,Kz=Kv.Practically, the above one-dimensional integral can be implemented by summing over the available discrete values of Kw and the digital implementation of the Fourier transform or its inverse; that is, the forward/inverse DFT can be realized through efficient fast Fourier transform (FFT) and inverse FFT (IFFT) algorithms. Thus, according to (14)~(17), the discrete form of the target function reconstruction at (x=xn,y,z) is expressed as(18)f^x=xn,y,z=IFFTKy=KuKz=Kv∑KwFFTu,vSKw,u,vSxn∗Kw,Ku,KvΔKw.In summary, the reconstruction procedure of the 3DRSA is illustrated in Figure2, which contains the following steps.Figure 2
Flowchart of 3D range stacking algorithm for concealed weapon detection.Step 1.
Perform an IFFT on the measured responseS(Kw,u,v) with respect to the fast-time frequency wavenumbers Kw. Compensate the propagation losses, and reduce the coupling and noises beyond the region of interest by range gating, and then the compensated signal is upsampled along the y- and z-axes through zero-padding to avoid any aliasing, followed by transforming to the 3D Fourier domain via a 3D FFT.Step 2.
The resultant 3D Fourier domain signal from Step1 is first matched filtered with the reference signal spectrum at the range bin xn, Sxn∗(Kw,Ku,Kv). The result is summed over the available fast-time frequency wavenumbers Kw to yield the marginal Fourier transform of the target function Fy,z(x=xn,Ky,Kz) at the range xn with Ky=Ku and Kz=Kv. Then a 2D IFFT with respect to (Ky,Kz) is performed on Fy,z(x=xn,Ky,Kz) to obtain the target reconstruction f^(x=xn,y,z) at the range bin xn. For all available range bins xns, the above calculations are implemented.Step 3.
Stack all outcomes of Step2 together to form a 3D image of the target function, f^(x,y,z); then apply possible postprocessing procedures such as smoothing filtering, geometry regridding, and object recognition and classification [27, 28] and finally display the image.


## 3.2. Sampling Criteria and Resolutions

According to the Nyquist sampling criterion [29], the sampling interval along some specific direction is inversely proportional to the spatial frequency domain support band of the signal in that direction.As for the wave propagation direction, the frequency domain support band of the transmitted signal is determined by(19)ΩKw=2Kwmin,2Kwmax.Thus the following Nyquist sampling constraint should be met in order to avoid aliasing along the slant range direction:(20)RmaxKwmax-Kwmin/ΔKw≤2π2Kwmax-Kwmin⟹ΔKw≤πRmax⟹Δw≤CπRmax,where Δw is the frequency sampling interval and Rmax denotes the maximum observed distance within the antenna beam pattern.As for the sampling criterion along theY/U direction, the samples are realized via a real or synthetic antenna array along this direction, and the frequency wavenumber support band ΩKy depends on the antenna beamwidth of the individual antenna ΦY through(21)ΩKy=-2Kwsin⁡ΦY2,2Kwsin⁡ΦY2.Thus the Nyquist sampling spacing Δu along the Y/U axis for gathering the echoed signal must satisfy the following inequality:(22)Δu≤2π4Kwsin⁡ΦY/2⟹Δu≤2π4Kwmaxsin⁡ΦY/2=λmin4sin⁡ΦY/2,where λmin is the smallest wavelength at the highest frequency Kwmax and the antenna beamwidth ΦY is determined by the antenna size along the corresponding direction. The worst-case situation corresponds to the case of ΦY = 180° for which the Nyquist sampling spacing in the Y/U dimension is(23)Δu≤λmin4.As for the case along theZ/V direction, a similar conclusion can be drawn that the Nyquist sampling criterion in the Z/V dimension takes the form (24)Δv≤λmin4sin⁡Φz/2,where Φz denotes the antenna beamwidth along the Z/V direction.The frequency wavenumber support band and the derived sampling constraints depend on the dimensions of the volume occupied in the 3D wavenumber domain. Therefore, the spatial resolutions of the resulting 3D reflectivity imagef^(x,y,z) are determined by the signal bandwidth, the center frequency, and the effective synthetic aperture. The achievable −3 dB resolutions along the x-, y-, and z-axes, are approximated, respectively, by(25)δx≈0.886πKwmax-Kwmin,(26)δy≈0.886λc4sin⁡ΦY/2,(27)δz≈0.886λc4sin⁡ΦZ/2,where λc is the wavelength at the center frequency fc and the units of δx, δy, and δz are all meters. In practice, due to the sidelobe structure of the radiation patterns, the final resolutions are usually slightly poorer than those theoretical values given by (25)~(27). In spite of this, both simulation and experimental results provided in the next section show that the 3D RSA is well suited to reconstruct 3D microwave imaging in concealed weapon detection at a high-precision manner.


## 4. Experimental Results

### 4.1. Numerical Simulation

In the computer simulations, a human body model consisting of 816 unit-amplitude point scatterers is numerically generated; as in Figure3, the model has pixel spacings of 2.5 cm and 5.0 cm on the horizontal plane and in the height direction, respectively; therefore it is practically confined within a box of dimensions 25.0 cm (x-axis) × 50.0 cm (y-axis) × 170.0 cm (z-axis) (the definition of the coordinate is in accordance with that of Figure 1, with x-, y-, and z-axes representing the range, azimuth, and height directions, resp.). Assuming that the human body is facing the planner scanning antenna aperture and considering that the microwave signals are not capable of penetrating the human body, therefore only 486 points from the model representing the surface of the front part of the body are picked out to generate the simulated radar echoes. The measurement parameters used in the simulations are listed in Table 1. Though a rectangular radiation pattern is assumed for each antenna element for simplicity, we will see in the latter real data experiments that the effects of the sidelobes from the practical antenna radiation pattern has a trivial influence on the performance of the proposed imaging algorithm.Table 1
Measurement parameters used in the simulation experiment.
Parameters
Value
Center frequency
35 GHz
Sweeping frequency points
101
y-axis samples number
151
y-axis sampling interval
6 mm
y-axis aperture length
0.9 m
y-axis antenna beamwidth
30°
Theoreticaly-axis resolution
6.8 mm
Frequency bandwidth
10 GHz
Theoreticalx-axis resolution
1.3 cm
z-axis samples number
301
z-axis sampling interval
6 mm
z-axis aperture length
1.8 m
z-axis antenna beamwidth
30°
Theoreticalz-axis resolution
6.8 mmFigure 3
Numerically generated human body model with an ensemble of 816 point scatterers, where 486 of the scatterers from the front of the body (red points) are utilized to generate the simulated radar echoes, whereas the left ones (blue crosses) represent the body part that cannot be penetrated by microwave: (a) 3D view; (b) side view; (c) front view.
(a)
(b)
(c)After performing the processing chain of Figure2, a 3D image of the human body is reconstructed, as in Figure 4(a). The proposed 3D imaging algorithm shows good imaging performance; for example, the reflectivity values as well as the positions of the point targets are all in perfect agreement with those of the simulated model. To quantitatively evaluate the performance, a point target with the coordinate (x,y,z) = (0.5 m, 0, 0) is simulated independently, the reconstructed point target is shown in Figures 5 and 6, and the figures of merits like resolution, peak sidelobe ratio (PSLR), and integrated sidelobe ratio (ISLR) are evaluated and compared with their theoretic values, which are indicated in Table 2. All the above measured parameters indicate that the developed algorithm provides a proper and accurate implementation of the 3D reconstruction.Table 2
The measured values of PSLR, ISLR and resolutions.
Measured parameter
x-axis/range
y-axis/azimuth
z-axis/height
PSLR (dB)
-
13.015
-
13.102
-
13.092
ISLR (dB)
-
10.355
-
10.347
-
10.336
Resolution (cm)
1.342
0.684
0.685Figure 4
(a) Isosurfaces (−3 dB) of the 3D reconstructed image processed by the proposed algorithm; maximum values projections of (b) side view and (c) front view.
(a)
(b)
(c)Figure 5
The impulse response on the (a) azimuth-range plane and (b) azimuth-height plane. The displayed dynamic range is 40 dB.
(a)
(b)Figure 6
(a) Range, (b) azimuth, and (c) height profiles of the reconstructed target impulse response.
(a)
(b)
(c)Actually, as for the people screening system, the resolution along the azimuth and height direction can be easily achieved, while being compared with that along the wave propagation (range) direction, so that the main lobe of the reconstructed impulse response along the range direction is much wider than the main lobes along other directions. Accordingly, if the 3D image is projected onto a 2D plane directly, the energy of the main lobe along the range direction would possibly be spread out across other directions due to the shift-varying property of the impulse response. To overcome this problem, the 3D image is collapsed using a maximum-value projection into a fully focused 2D image of the target for display, rather than being directly projected. This is illustrated in Figures4(b) and 4(c), where the maximum values of the 3D reconstructed image are orthogonally projected onto the side and front planes for display. The scatterers on both views are basically distinguishable; the front view of Figure 4(c) is usually the image that is really being of interest, where concealed weapons and any threat can be detected effectively in this image.


### 4.2. Real-Field Experiment Result and Discussion

A wide-band imaging system has been developed that utilizes a scanned linear array of antennas to quickly gather the image data for subsequent reconstruction and display. The system consists of a vertically scanned linear antenna array which consists of 281 virtual elements which are formed by the multiple-input multiple-output (MIMO) technology [30, 31]. The antenna array is driven by a wide-band millimeter-wave transceiver which transmits and receives step-frequency continuous-wave (SF-CW) signals and is demodulated to baseband through analogue dechirping. The array antenna is mounted on a large fast horizontal mechanical scanner, which allows a full-image data set to be gathered in less than 1 s. During a scan, the array is sequentially switched at a very high speed to gather each column of data as the scanner is moving from left or right, thereby gathering a full 2D aperture of data.The system operates in Ka band and each antenna element operates in the HH polarization state and has a highest sidelobe level of −25 dB. A typical set of measurement parameters are listed in Table3, where the current target scene consists of a human body model carrying a concealed knife near the left waist and a small metallic sphere around the neck. The wide-band data are sampled vertically across the array and horizontally over the synthetic aperture and are digitized by an A/D converter for subsequent storage in the computer. The proposed 3D range stacking algorithm is programmed with CUDA C language to run on an NVIDIA TITAN X GPU platform parallel, which allows it to be executed at very high speed with an average period of 0.7 seconds to generate a frame of 3D image of the human body, compared to the average period of 5 seconds or more for a typical STOLT-interpolated counterpart. So the proposed algorithm is more suitable for the real-time application. It is worse noting that the computation complexity (or the processing time) of the imaging algorithm is a monotone increasing function of the frequency bandwidth, array element number, and the image samples number, and the practical processing speed also depends on the specific hardware platform to be used.Table 3
Measurement parameters of the real data experiment.
Parameters
Value
Center frequency
35 GHz
Array element number
281
Array element spacing
5 mm
Array length
1.4 m
Frequency bandwidth
10 GHz
A
z
i
m
u
t
h
∗ samples number
181
Azimuth element spacing
5 mm
Azimuth aperture length
0.9 m
T
h
e
∗ azimuth denotes the direction along the y-axis in Figure 1, whereas the array direction denotes the height direction along the z-axis.Figure7(a) shows the optical image of the human body model to be imaged. A 2 cm diameter metallic sphere fixed by the adhesive tape is exposed near the neck in the front of the body; however, a knife composed of a stain-less steel blade and a plastic handle is also concealed under the model’s jacket near its left waist. Figure 7(b) is the front view of the 3D microwave image result. Due to the high sensitivity of the transceiver design, as well as the high-precision focus of the imaging algorithm, a significantly high dynamic range is apparent in the microwave image. It can be seen that microwaves readily pass through the clothing materials and are reflected by the body and the concealed knife. In this experiment, the 3D range stacking imaging techniques incorporating large apertures and wide angle illumination achieve the high-resolution reconstruction of the targets, and the theoretical resolution for the image is approximately 5 mm, where the differences in the shapes and reflectivity of the targets, such as the metallic sphere, the knife handle, blade, and the model’s body, are readily apparent in the reconstructed imagery. Given the remarkable quality of the reconstructed image, it can also be concluded that the effects of the sidelobes of the antenna radiation pattern on the final image are negligible.Figure 7
Photograph of the human body model carrying a concealed knife and a metallic sphere (a); orthogonal projection of the 3D microwave image result onto the azimuth-elevation (Y-Z) plane for display (b).
(a)
(b)It is worth noting that, given a sufficient small point-like target (as compared to the system resolution), the resolution, PSLR, and ISLR can be measured from the reconstructed image, and this kind of point target usually works as a calibrator to calibrate any errors (such as the multichannel amplitude-phase errors of the array) of the imaging system. Other metrics of the image may include the image contrast relating the level of the signal to noise-clutter ratio, which is important to enhance the succeeding concealed weapon detection performance (such as the detection probability and the precision).This experiment demonstrates microwave imaging can effectively view objects through barriers that are opaque in the optical band. It also demonstrates that the proposed 3D imaging algorithm can achieve the high resolution and image fidelity.


## 4.1. Numerical Simulation

In the computer simulations, a human body model consisting of 816 unit-amplitude point scatterers is numerically generated; as in Figure3, the model has pixel spacings of 2.5 cm and 5.0 cm on the horizontal plane and in the height direction, respectively; therefore it is practically confined within a box of dimensions 25.0 cm (x-axis) × 50.0 cm (y-axis) × 170.0 cm (z-axis) (the definition of the coordinate is in accordance with that of Figure 1, with x-, y-, and z-axes representing the range, azimuth, and height directions, resp.). Assuming that the human body is facing the planner scanning antenna aperture and considering that the microwave signals are not capable of penetrating the human body, therefore only 486 points from the model representing the surface of the front part of the body are picked out to generate the simulated radar echoes. The measurement parameters used in the simulations are listed in Table 1. Though a rectangular radiation pattern is assumed for each antenna element for simplicity, we will see in the latter real data experiments that the effects of the sidelobes from the practical antenna radiation pattern has a trivial influence on the performance of the proposed imaging algorithm.Table 1
Measurement parameters used in the simulation experiment.
Parameters
Value
Center frequency
35 GHz
Sweeping frequency points
101
y-axis samples number
151
y-axis sampling interval
6 mm
y-axis aperture length
0.9 m
y-axis antenna beamwidth
30°
Theoreticaly-axis resolution
6.8 mm
Frequency bandwidth
10 GHz
Theoreticalx-axis resolution
1.3 cm
z-axis samples number
301
z-axis sampling interval
6 mm
z-axis aperture length
1.8 m
z-axis antenna beamwidth
30°
Theoreticalz-axis resolution
6.8 mmFigure 3
Numerically generated human body model with an ensemble of 816 point scatterers, where 486 of the scatterers from the front of the body (red points) are utilized to generate the simulated radar echoes, whereas the left ones (blue crosses) represent the body part that cannot be penetrated by microwave: (a) 3D view; (b) side view; (c) front view.
(a)
(b)
(c)After performing the processing chain of Figure2, a 3D image of the human body is reconstructed, as in Figure 4(a). The proposed 3D imaging algorithm shows good imaging performance; for example, the reflectivity values as well as the positions of the point targets are all in perfect agreement with those of the simulated model. To quantitatively evaluate the performance, a point target with the coordinate (x,y,z) = (0.5 m, 0, 0) is simulated independently, the reconstructed point target is shown in Figures 5 and 6, and the figures of merits like resolution, peak sidelobe ratio (PSLR), and integrated sidelobe ratio (ISLR) are evaluated and compared with their theoretic values, which are indicated in Table 2. All the above measured parameters indicate that the developed algorithm provides a proper and accurate implementation of the 3D reconstruction.Table 2
The measured values of PSLR, ISLR and resolutions.
Measured parameter
x-axis/range
y-axis/azimuth
z-axis/height
PSLR (dB)
-
13.015
-
13.102
-
13.092
ISLR (dB)
-
10.355
-
10.347
-
10.336
Resolution (cm)
1.342
0.684
0.685Figure 4
(a) Isosurfaces (−3 dB) of the 3D reconstructed image processed by the proposed algorithm; maximum values projections of (b) side view and (c) front view.
(a)
(b)
(c)Figure 5
The impulse response on the (a) azimuth-range plane and (b) azimuth-height plane. The displayed dynamic range is 40 dB.
(a)
(b)Figure 6
(a) Range, (b) azimuth, and (c) height profiles of the reconstructed target impulse response.
(a)
(b)
(c)Actually, as for the people screening system, the resolution along the azimuth and height direction can be easily achieved, while being compared with that along the wave propagation (range) direction, so that the main lobe of the reconstructed impulse response along the range direction is much wider than the main lobes along other directions. Accordingly, if the 3D image is projected onto a 2D plane directly, the energy of the main lobe along the range direction would possibly be spread out across other directions due to the shift-varying property of the impulse response. To overcome this problem, the 3D image is collapsed using a maximum-value projection into a fully focused 2D image of the target for display, rather than being directly projected. This is illustrated in Figures4(b) and 4(c), where the maximum values of the 3D reconstructed image are orthogonally projected onto the side and front planes for display. The scatterers on both views are basically distinguishable; the front view of Figure 4(c) is usually the image that is really being of interest, where concealed weapons and any threat can be detected effectively in this image.


## 4.2. Real-Field Experiment Result and Discussion

A wide-band imaging system has been developed that utilizes a scanned linear array of antennas to quickly gather the image data for subsequent reconstruction and display. The system consists of a vertically scanned linear antenna array which consists of 281 virtual elements which are formed by the multiple-input multiple-output (MIMO) technology [30, 31]. The antenna array is driven by a wide-band millimeter-wave transceiver which transmits and receives step-frequency continuous-wave (SF-CW) signals and is demodulated to baseband through analogue dechirping. The array antenna is mounted on a large fast horizontal mechanical scanner, which allows a full-image data set to be gathered in less than 1 s. During a scan, the array is sequentially switched at a very high speed to gather each column of data as the scanner is moving from left or right, thereby gathering a full 2D aperture of data.The system operates in Ka band and each antenna element operates in the HH polarization state and has a highest sidelobe level of −25 dB. A typical set of measurement parameters are listed in Table3, where the current target scene consists of a human body model carrying a concealed knife near the left waist and a small metallic sphere around the neck. The wide-band data are sampled vertically across the array and horizontally over the synthetic aperture and are digitized by an A/D converter for subsequent storage in the computer. The proposed 3D range stacking algorithm is programmed with CUDA C language to run on an NVIDIA TITAN X GPU platform parallel, which allows it to be executed at very high speed with an average period of 0.7 seconds to generate a frame of 3D image of the human body, compared to the average period of 5 seconds or more for a typical STOLT-interpolated counterpart. So the proposed algorithm is more suitable for the real-time application. It is worse noting that the computation complexity (or the processing time) of the imaging algorithm is a monotone increasing function of the frequency bandwidth, array element number, and the image samples number, and the practical processing speed also depends on the specific hardware platform to be used.Table 3
Measurement parameters of the real data experiment.
Parameters
Value
Center frequency
35 GHz
Array element number
281
Array element spacing
5 mm
Array length
1.4 m
Frequency bandwidth
10 GHz
A
z
i
m
u
t
h
∗ samples number
181
Azimuth element spacing
5 mm
Azimuth aperture length
0.9 m
T
h
e
∗ azimuth denotes the direction along the y-axis in Figure 1, whereas the array direction denotes the height direction along the z-axis.Figure7(a) shows the optical image of the human body model to be imaged. A 2 cm diameter metallic sphere fixed by the adhesive tape is exposed near the neck in the front of the body; however, a knife composed of a stain-less steel blade and a plastic handle is also concealed under the model’s jacket near its left waist. Figure 7(b) is the front view of the 3D microwave image result. Due to the high sensitivity of the transceiver design, as well as the high-precision focus of the imaging algorithm, a significantly high dynamic range is apparent in the microwave image. It can be seen that microwaves readily pass through the clothing materials and are reflected by the body and the concealed knife. In this experiment, the 3D range stacking imaging techniques incorporating large apertures and wide angle illumination achieve the high-resolution reconstruction of the targets, and the theoretical resolution for the image is approximately 5 mm, where the differences in the shapes and reflectivity of the targets, such as the metallic sphere, the knife handle, blade, and the model’s body, are readily apparent in the reconstructed imagery. Given the remarkable quality of the reconstructed image, it can also be concluded that the effects of the sidelobes of the antenna radiation pattern on the final image are negligible.Figure 7
Photograph of the human body model carrying a concealed knife and a metallic sphere (a); orthogonal projection of the 3D microwave image result onto the azimuth-elevation (Y-Z) plane for display (b).
(a)
(b)It is worth noting that, given a sufficient small point-like target (as compared to the system resolution), the resolution, PSLR, and ISLR can be measured from the reconstructed image, and this kind of point target usually works as a calibrator to calibrate any errors (such as the multichannel amplitude-phase errors of the array) of the imaging system. Other metrics of the image may include the image contrast relating the level of the signal to noise-clutter ratio, which is important to enhance the succeeding concealed weapon detection performance (such as the detection probability and the precision).This experiment demonstrates microwave imaging can effectively view objects through barriers that are opaque in the optical band. It also demonstrates that the proposed 3D imaging algorithm can achieve the high resolution and image fidelity.


## 5. Conclusions

This paper presented a 3D imaging algorithm especially tailored for accurately reconstructing the 3D image of the human body for concealed weapon detection under planner scanning geometry. The algorithm, which is referred to as 3D range stacking algorithm, allows for the propagation losses compensation, wavefront reconstruction, and aliasing mitigating and without interpolations. Moreover, it suffers no approximations and truncation errors. The essence of the algorithm is the multiplication of the reference signal spectrum at each specific range bin, followed by the coherent summation over the whole frequency band, instead of the 3D frequency interpolation. In addition, the sampling criteria and the resolutions for planner scanning geometry are also derived. The imaging results using simulated data and experimental data demonstrated the high resolution and good image fidelity of the proposed algorithm.



---

*Source: 1480623-2017-08-07.xml*