Millimeter-wave compressive holography
Christy Fernandez Cull,
1,
* David A. Wikner,
2
Joseph N. Mait,
2
Michael Mattheiss,
3
and David J. Brady
1
1
Fitzpatrick Institute for Photonics and Department of Electrical and Computer Engineering,
Duke University, 129 Hudson Hall, Durham, North Carolina 27708, USA
2
United States Army Research Laboratory, 2800 Powder Mill Road,
Adelphi, Maryland 20783, USA
3
Department of Electrical and Computer Engineering, University of Maryland,
2405 A.V. Williams Building, College Park, Maryland 20742, USA
*Corresponding author: caf11@ee.duke.edu
Received 16 December 2009; revised 26 March 2010; accepted 16 April 2010;
posted 19 April 2010 (Doc. ID 121568); published 12 May 2010
We describe an active millimeter-wave holographic imaging system that uses compressive measurements
for three-dimensional (3D) tomographic object estimation. Our system records a two-dimensional (2D)
digitized Gabor hologram by translating a single pixel incoherent receiver. Two approaches for compres-
sive measurement are undertaken: nonlinear inversion of a 2D Gabor hologram for 3D object estimation
and nonlinear inversion of a randomly subsampled Gabor hologram for 3D object estimation. The object
estimation algorithm minimizes a convex quadratic problem using total variation (TV) regularization for
3D object estimation. We compare object reconstructions using linear backpropagation and TV minimi-
zation, and we present simulated and experimental reconstructions from both compressive measurement
strategies. In contrast with backpropagation, which estimates the 3D electromagnetic field, TV minimi-
zation estimates the 3D object that produces the field. Despite undersampling, range resolution is
consistent with the extent of the 3D object band volume. © 2010 Optical Society of America
OCIS codes: 090.1995, 100.3200, 100.6950, 110.1758, 110.3010, 110.3200.
1. Introduction
Various methods exist for concealed weapons detec-
tion [1]. These methods aim to penetrate common
obstructions such as clothing or plastics. X-ray [2]
and millimeter-wave (MMW) [3] imaging systems
are technologies capable of penetrating these bar-
riers for imaging suicide bomb vests or weapons com-
posed of metals, nonmetals, or plastics. While x-ray
imaging capabilities are highly effective, questions
about health risks impair the feasibility of such sys-
tems for real-time imaging. MMW for low-power (on
the order of mW) imaging systems do not present a
health hazard and therefore enable real-time
imaging of targets with high contrast and high
resolution.
Several studies have explored both active and pas-
sive MMW imagers for concealed weapons detection
[4,5] where the system limitation is the detector ar-
ray cost. Some systems include portal, or handheld
devices [6] operating in close range to the target.
Other systems are holographic [7,8]. These MMW fo-
cal and interferometric systems map object informa-
tion onto a two-dimensional (2D) or a linear array
that is typically scanned for image formation. These
scanning systems [9–11] are plagu ed by their asso-
ciated data acquisition times. For these systems,
there is a trade-off between scan time and measure-
ment signal-to-noise ratio (SNR). Therefore, rapid
scanning of conceal ed weapons is challenging for cur-
rent MMW systems.
For stand-off explosive detection, rapid scanning of
the target is a necessity. To overcome the bottleneck
associated with current MMW scanning systems, we
consider compressive sensing (CS). Recent studies in
0003-6935/10/190E67-16$15.00/0
© 2010 Optical Society of America
1 July 2010 / Vol. 49, No. 19 / APPLIED OPTICS E67
CS reveal that an N-point image can be restored
from M measurements, where M ≪ N [5,12–14].
Chan et al. [13] used a focal system to randomly sam-
ple spatial frequencies in the Fourier plane for 2D
object estimation at 100 GHz. We are further moti-
vated to investigate CS for MMW imaging based on
similar work in 633 nm compressive holography [15].
Holography, which measures a limited set of spatial
frequencies in the Fourier domain, is a compressive
encoder, because it compresses three-dimensional
(3D) spatial information into a single interferometric
planar field. Since the entire extent of an object’s3D
spatial frequency band volume can not be captured in
a single exposure, mul tiangle illumination or object
rotation is typically used to improve 3D object esti-
mation. The results in [15] suggest that 3D tomo-
graphic estimation can be achieved from a 2D
hologram recording.
This paper extends compressive holography to
millimeter wavelengths. The subject matter in this
paper differs from compressive holography at visible
wavelengths, because sparse holographic sampling is
implemented to minimize the data acquisition scan
cost associated with imaging at millimeter wave-
lengths. Also, the MMW holography system operates
in a completely different wavelength region with cor-
responding differences in optics and detectors, so it
required a completely new system design compared
to work at visible wavelengths. Further, in this pa-
per, we holographically sample a subset of spatial fre-
quencies for 3D object estimation. Al so, we randomly
subsample a 2D hologram to further analyze the im-
pact of fewer measurements on 3D object estimation.
Our holographic technique is similar to [8], as we are
not band limited by a lens aperture and phase infor-
mation is preserved. We differ in our nonlinear inver-
sion approach for 3D object estimation. Our method
optimizes a convex quadratic problem using total
variation (TV) regularization. Unlike classical recon-
struction with backpropagation, we demonstrate
that undiffra cted fields, overlaid in the frequency
domain of a Gabor hologram, can be separated by
imposing a TV sparsity constraint.
Although other contributions in the literature em-
body a mathematical framework similar to compres-
sive holography [16–18], there exists a fundamental
difference in philosophy. Compressive holography
exploits encoding and undersam pling for 3D object
estimation, whereas techniques in diffraction tomo-
graphy are designed to overcome sampling limita-
tions imposed by the data collection process. Also,
recent work by Denis et al. [19] presents a similar
twin-image suppression method; however, a spar-
sity-enforcing prior in a Bayesian framework com-
bined with l1 regulariz ation is used for object
estimation. Our work represents the confluence of
MMW digitized holographic measurement and TV
minimization for 3D object estimation with minimal
error. We adapted the algorithm framework of [15]
for sparse holographic sampling and data inversion.
This paper is organized as follows. In Section 2,we
describe the theoretical background for diffraction to-
mography and holographic measurement. Hologram
recording geometry and resolution metrics are also
discussed in this section. Section 3 summarizes our
TV minimization algorithm used for 3D object estima-
tion from a 2D digitized composite hologram. Also, we
present simulated data of randomly subsampled 2D
holograms to analyze the impact of fewer measure-
ments on 3D object estimation. Section 4 describes
the experimental platform. Section 5 presents TV
minimization and backpropagation reconstructions.
Finally, in Section 6, we provide a summary of the
results and concluding remarks.
2. Theory
Our ultimate goal is to make the smallest number of
measurements about a 3D (x
0
; y
0
; z
0
) object f
o
ðr
0
Þ,
where r
0
is a 3D spatial vector, such that it is possible
to reconstruct f
o
ðr
0
Þ with minimal error. Rather than
attempting to form an image of f
o
ðr
0
Þ point by point,
our approach is based on making measurements in
the far field where spatial frequencies (u
x
and u
y
)
are measured. To do this, we record a hologram. A ho-
logram gðr
h
Þ is a record of the interference between
two wave fields, a reference field E
r
ðrÞ and an object
scattered field E
o
ðrÞ. To record a hologram, a square-
law detector in the hologram plane r
h
¼ðx; y; z
h
Þ mea-
sures a time-averaged intensity of the interference:
gðr
h
Þ¼∣E
r
ðr
h
ÞþE
o
ðr
h
Þ∣
2
¼ ∣E
r
ðr
h
Þ∣
2
þ ∣E
o
ðr
h
Þ∣
2
þ 2∣E
r
ðr
h
ÞE
o
ðr
h
Þ∣ cos½θ
r
ðr
h
Þ − θ
o
ðr
h
Þ; ð1Þ
where θ
r
ðr
h
Þ represents the phase associated with the
propagated reference wave field and θ
o
ðr
h
Þ represents
the phase associated with the propagated object wave
field. We assume our object field is generated by illu-
minating a 3D object f
o
ðr
0
Þ by an on-axis plane wave
expð−2πiu
o
· r), where u
o
¼ðu
x
o
; u
y
o
; u
z
o
Þ. The 3D ob-
ject f
o
ðr
0
Þ represents an object scattering amplitude
where after reference plane wave illumination, a frac-
tion of the energy is either transmitted or reflected at
a point in 3D space. We do not assume the object in-
duces any phase change in an incident wave field
through polarization or birefringence.
If the object is transmissive and located at z
h
dis-
tance away from the hologram plane, under the Born
approximation the scattered field is
E
o
ðr
h
Þ¼
−π
λ
2
Z
E
r
ðr
0
Þf
o
ðr
0
Þhðr
h
− r
0
Þdr
0
; ð2Þ
where hðr
h
− r
0
Þ is the shift-invariant impulse
response and E
r
ðr
0
Þ is the reference plane wave.
For scalar waves in homogeneous space, the impulse
response is
hðr
h
− r
0
Þ¼
expð2πi∣r
h
− r
0
∣=λÞ
∣r
h
− r
0
∣
: ð3Þ
E68 APPLIED OPTICS / Vol. 49, No. 19 / 1 July 2010
We can reformulate the convolution integral in Eq.
(2) using the Fourier convolution theorem. The Four-
ier transform of the scattered field along the trans-
verse axes in the recording plane is
^
E
0
ðu
x
; u
y
; z
h
Þ¼
1
iπλ
^
f
o
u
x
− u
x
o
; u
y
− u
y
o
;
×
ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
1
λ
2
− u
2
x
− u
2
x
r
− u
z
0
G
2D
ðu
x
; u
y
; zÞ; ð4Þ
where
G
2D
ðu
x
; u
y
; zÞ¼
exp
2πiz
ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
1
λ
2
− u
2
x
− u
2
y
q
ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
1
λ
2
− u
2
x
− u
2
y
q
; ð5Þ
^
f
o
is the 3D Fourier-trans form of the object density,
and the exponential term represents a propagation
transfer function. Under the small angle approxi-
mation, u
z
¼ 1=λ and u
x
, u
y
≤ 1=λ. The frequency-
domain scattered field is then approximated by
^
E
o
ðu
x
; u
y
; zÞ¼
1
iπλ
^
f
o
u
x
− u
x
o
; u
y
− u
y
o
; −
λ
2
ðu
2
x
þ u
2
y
Þ
× exp
2πiz
ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
1
λ
2
− u
2
x
− u
2
y
r
: ð6Þ
As discussed in Section 3, digital processing of the
Gabor hologram aims to isolate the scattered field,
E
o
ðr
h
Þ, signal term from background and conjugate
terms. If we assume the recorded hologram measures
E
o
ðr
h
Þ directly and that E
r
ðr
h
Þ¼1, Eq. (6) demon-
strates that a 2D hologram captures a 3D parabolic
slice of the object’s band volume. Figure 1 describes
tomographic sampling of a 3D band volume in a
Gabor geometry. Typically, the illumination (or the
object) must be rotated to fully sample the 3D band
volume. To increase longitudinal resolution, the sys-
tem may alternatively be scanned in frequency. In-
stead of scanning in the frequency domain on the
surface of a sphere, this approach allows one to scan
a spherical shell with radii corresponding to a wave-
length range. In [15], we showed that one may apply
CS theory to estimate the 3D distribution of f ðr
0
Þ
from a single holographic image without scanning
in frequency or rotating the object. Instead, we
exploit sparsity. This work extends 3D tomographic
estimation from 2D holographic measurements to
millimeter wavelengths.
A. Hologram Recording Geometry
The geometry used to record the hologram impacts
the postdetection signal processing and the ability
to reconstruct the image. In an off-axis geometry
[20], the signal and its conjugate are separated from
each other in frequency space and from the on-axis
undiffracted energy [see Fig. 2(a)].
Note that the maximum spatial frequency u
max
that the detector can record is limited by the sam-
pling pitch (dx ) of the detector:
u
max
¼
1
2dx
: ð7Þ
With our detector, the pixel pitch is set by the WR-08
waveguide size (2:32 mm × 1:08 mm). The maximum
spatial frequency recorded in the vertical direction is
0:463 mm
−1
, and in the horizontal direction it
is 0:216 mm
−1
.
Figure 2(a) shows that the information content of
the object and the pixel pitch of the detector impose
minimum and maximum limits on the angle θ
cz
of the
off-axis reference. For simplicity, we assume the re-
ference beam has no y component. To separate the
object from its squared magnitude without ambigu-
ity due to detector aliasing or from confusion with un-
diffracted terms, the angle of the off-axis beam must
satisfy
θ
c min
≤ θ
cz
≤ θ
c max
; ð8Þ
where
θ
c max
¼ sin
−1
λ
2
1
dx
− u
B
; ð9Þ
Fig. 1. Fourier-transform domain sampling of the object band volume in a transmission geometry. (a) 2D slice of a 3D sphere where the
dotted curve represents the measurement from single plane wave illumination. (b) Rectilinear pattern represents wave vectors sampled by
the hologram due to a finite detector plane sampling. (c) Wave normal sphere cross section for spatial and axial resolution analysis.
1 July 2010 / Vol. 49, No. 19 / APPLIED OPTICS E69
θ
c min
¼ sin
−1
3
λu
B
2
; ð10Þ
and u
B
is the spatial frequency bandwidth of the
object.
Use of an off-axis reference beam simplifies our re-
construction because digital signal processing allows
one to yield an estimate of the object field
~
g
offaxis
¼ E
r
ðr
h
Þ
E
o
ðr
h
Þ: ð11Þ
In an on-axis geometry, it is more difficult to separate
the object field from the undiffracted, zero-order
fields and from its conjugate. The overlap of these
three fields degrades resolution and contrast in
the object reconstruction. One can apply DC suppres-
sion techniques to enhance object reconstructions
[21–23], and measurements of the energy in the re-
ference beam alone can be made and subtracted from
the hologram:
~
g
onaxis
¼ ∣E
o
ðr
h
Þ∣
2
þ E
r
ðr
h
ÞE
o
ðr
h
ÞþE
r
ðr
h
ÞE
o
ðr
h
Þ:
ð12Þ
In this paper, we record a hologram in an on-axis geo-
metry because the need for an increased bandwidth
in the off-axis case outweighs the complexity for
on-axis object isolation.
B. Hologram Plane Sampling and Resolution Metrics
In our implementation, the field gðr
h
Þ is sampled and
digitized into a 2D matrix by translating a point de-
tector in x and y at the hologram plane z
h
. An ana-
lytical discussion of the discrete model is detailed
in [15] and addressed in Section 3.
Holographic measurements captured digitally, by
a scanning detector, are related to measurements
made in the spatial frequency domai n. The total
number of detector measurements N is
N ¼ n
x
n
y
; ð13Þ
where n
x
and n
y
represent the number of measure-
ments along each spatial dimension in x and y. Given
the detector sampling pitch (dx), the number of
measurements (n
x
) in the hologram plane along
the horizontal dimension is given by
n
x
¼
W
x
dx
; ð14Þ
where
Fig. 2. (Color online) (a) Spectrum for an off-axis hologram recording depicting an inherent increase in bandwidth for adequate object
separation from undiffracted terms. (b) Spectrum for a Gabor hologram recording, depicting the overlay of undiffracted, object, and con-
jugate terms. (c) Transverse slices from linear inverse propagation results at various z planes.
E70 APPLIED OPTICS / Vol. 49, No. 19 / 1 July 2010
W
x
¼
λz
Δx
o
; ð15Þ
z is the distance between the object and the detector,
Δx
o
is the one-dimensional spatial resolution with
which we wish to image the object, and W
x
refers
to the spatial extent of a diffracted object (Δx
o
).
The inverse scaling relationship arises from the
conjugate relationship between the object and the
hologram [24].
Detector sampling over a finite field size affects
sampling resolution in the frequency domain. The
sampling resolution, Δ
u
, in the frequency domain
along both the horizontal and vertical dimensions
(u
x
and u
y
), assuming n
x
¼ n
y
, is determined by the
sampling field size at the detector plane, Δ
u
¼
1=ð2n
x
dxÞ. The maximum spatial frequency sampled
by the detector is equal to u
max
in Eq. (7).
Resolution metrics, lateral (Δ
x
) and axial (Δ
z
), for
the Gabor geometry are determined by the illumina-
tion wavelength and the system numerical aperture
(NA). Based on the Gabor recording geometry, the
NA is defin ed by
n sin θ
u
¼
W
x
2z
¼
λ
2Δx
o
; ð16Þ
where n is the refractive index of air and θ
u
is the
half-angle subtended by the object to half the spatial
extent of the hologram plane (W
x
=2). Recall that the
spatial resolution is related to the inverse scaling re-
lationship in the frequency domain. The half-angle,
θ
u
, is also defined as the angular bandwidth sam-
pling on the wave normal sphere due to system
NA. Thus, the NA can also be described in the spatial
frequency domain. The wave vector geometry for
hologram recording is shown in Figs. 1(b) and 1(c).
Considering the geometry in Fig. 1(c), we write
sinðθ
u
Þ¼
Δu
x
∣u∣
; ð17Þ
where under the small angle approximation
∣u∣θ
u
¼ Δu
x
: ð18Þ
If we assume that NA ≈ θ
u
and ∣u∣≈1=λ, the spatial
resolution is equal to
Δ
x
¼
λ
NA
: ð19Þ
Similarly, from the wave vector geometry, the spatial
frequency resolution along zðΔu
z
Þ is determined by
Δu
z
¼ Δu
z;max
− Δu
z;min
¼ ∣u∣ð1 − cosðθ
u
ÞÞ ¼ ∣u∣θ
2
u
:
ð20Þ
Under the small angle approximation, the axial
resolution is
Δ
z
¼
λ
NA
2
: ð21Þ
After substituting the expression for NA from Eq.
(16) into Eqs. (19) and (21), we see that lateral reso-
lution is also defined as Δ
x
≈ 2Δx
o
and range resolu-
tion is defined as Δ
z
≈ 4Δx
2
o
=λ. Defining the lateral
and axial resolution using NA describes resolution
in terms of system geometry (a function of object dis-
tance), whereas the second metric is modeled as a
function of feature size, Δx
0
. The maximum of the
two measures for lateral and axial resolution pro-
vides a baseline metric for resolution. We use these
metrics for evaluating resolution from TV minimiza-
tion object reconstructions in Section 5.
This section provided motivation for implementing
a Gabor geometry instead of an off-axis approach.
The impact of detector sampling at the hologram
plane was addressed and a relation between object
sampling and frequency domain sampling was dis-
cussed. Finally, theoretical resolution metrics were
derived.
3. Reconstruction Methods and Simulations
In this section, we discuss two reconstruction meth-
ods: 3D object estimation from a Gabor hologram and
3D object estimation from randomly subsampled
Gabor holographic measurements. Subsampling is
implemented to further analyze the impact of com-
pressive measurement on 3D object estimation.
The continuous model for Gabor holography, mod-
eled under the first Born approximation, is shown in
Eq. (2). The detector plane is located at the z
h
¼ 0
plane in the rðx; y; z
h
Þ coordinate system. The object
data, f , are located in the r
0
ðx
0
; y
0
; z
0
Þ coordinate sys-
tem. The recorded hologram in Eq. (1) can be reformu-
lated if we assume that E
r
ðr
h
Þ¼1 and if operations on
f
o
ðr
0
Þ in the convolution integral in Eq. (2) are ex-
pressed using an operator, H. After squared-reference
field subtraction, we can represent the recorded holo-
gram in algebraic notation using
g ¼ ∣Hf ∣
2
þ Hf þ H
f þ n; ð22Þ
where g is an N × 1 vectorized detector measurement,
H is a 2D discrete system matrix, f is an M × 1 vector-
ized object representation [f
o
ðr
0
Þ], and n is the noise
associated with the measurement. If we ignore the ob-
ject-squared-field contribution in Eq. (22), we can es-
tablish a linear relationship between the detector
measurement and object field distribution, g ¼ Hf .
Once we record a digital hologram, our goal is to
estimate the object distribution, f . From [15], we
know that the digitized holographic measurement is
g
n
1
;n
2
¼
X
l
−1
2D
×
^
f
m
1
;m
2
;l
exp
ılΔ
z
ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi
1
λ
2
−m
2
1
Δ
2
u
−m
2
2
Δ
2
u
r
n
1
;n
2
;ð23Þ
1 July 2010 / Vol. 49, No. 19 / APPLIED OPTICS E71