UAV-Borne_FMCW_InSAR_for_Focusing_Buried_Objects

Download as pdf or txt
Download as pdf or txt
You are on page 1of 5

IEEE GEOSCIENCE AND REMOTE SENSING LETTERS, VOL.

19, 2022 4014505

UAV-Borne FMCW InSAR for Focusing


Buried Objects
Ralf Burr , Student Member, IEEE, Markus Schartel ,
Alexander Grathwohl , Graduate Student Member, IEEE, Winfried Mayer , Senior Member, IEEE,
Thomas Walter , Member, IEEE, and Christian Waldschmidt , Senior Member, IEEE

Abstract— Antipersonnel mines are hidden weapons that are


usually buried close to the surface and are triggered by a foot
step of the victim. A sensor that is able to detect minimum metal
mines is a ground penetrating synthetic aperture radar (GPSAR).
In previous work, an unmanned aerial vehicle (UAV)-based
GPSAR was developed and tested for the detection of buried
landmines. The real topography was neglected and the interface
between air and soil was assumed as a horizontal plane surface.
Since the dielectric properties of the soil reduce the propagation
velocity of the electromagnetic wave, the interface between air
and soil has to be known precisely for subsurface focusing.
Neglecting the surface profile has a negative impact on the
image quality and prevents the detection of buried objects. In
this letter, a novel two-step procedure is used to overcome this
issue. In time-division multiplex mode, the radar transmits two
frequency-modulated ramps in two different frequency bands. Fig. 1. UAV equipped with four TEM-Horn antennas. All antennas are
The data of the upper frequency band are used to generate aligned in vertical polarization.
a digital elevation model (DEM) through interferometry. The
data of the lower frequency band are used for GPSAR focusing
on the basis of the DEM. It is demonstrated that the GPSAR is based on the fact that a part of the electromagnetic wave
image quality is improved on nonideal planar surfaces and the is reflected when an abrupt change of the dielectric properties
probability of detection is increased.
occurs. Thus, it is possible to detect modern minimum metal
Index Terms— Antipersonnel mine, frequency-modulated mines in a reliable way. In case of a down-looking radar,
continuous-wave (FMCW) radar, ground penetrating radar the surface reflection is the first dominant reflection in the
(GPR), radar interferometry, unmanned aerial vehicle (UAV).
radar response. This reflection is greater than that of a buried
minimal metal mine. This strong surface reflection must either
I. I NTRODUCTION be removed by signal processing as demonstrated in [8] or can
be significantly reduced by a side-looking arrangement of the
T HE continuous development of unmanned aerial vehi-
cles (UAVs) as sensor platforms opens up more and more
possibilities in remote sensing. The potential of a completely
antenna [11]. In a side-looking arrangement, measurements
are performed along the synthetic aperture to achieve a high
remotely controlled or even completely autonomous flying resolution [12]. Due to the refractive index, the electromag-
platform is also suited for mine detection [1], [2]. netic wave propagates more slowly in the ground than in the
Ground penetrating radar (GPR) is a technology that has air. This must be considered for subsurface imaging [13].
been used for more than 20 years for the detection of buried In previous publications, the surface was approximated as
objects [3], [4]. In recent years, more and more GPR sensors a horizontal flat surface [11]. Since, in reality, it is rarely
have been integrated on UAVs [5]–[10]. The GPR detection horizontal and flat, a digital elevation model (DEM) is required
to improve the image quality.
Manuscript received October 28, 2020; revised February 3, 2021 and One of the most common standard techniques in Earth
April 13, 2021; accepted May 25, 2021. Date of publication July 15, 2021; observation to generate a DEM is the use of an interferometric
date of current version December 23, 2021. This work was supported by
the Urs Endress Foundation (https://www.uestiftung.org/). (Corresponding synthetic aperture radar (InSAR) [14], [15]. Two 2-D complex
author: Ralf Burr.) SAR images are generated over the same area from data
Ralf Burr and Thomas Walter are with the Laboratory of Microtechnol- recorded by two neighboring receiving antennas. From the
ogy, Ulm University of Applied Sciences, 89081 Ulm, Germany (e-mail:
ralf.burr@thu.de). phase difference of the two SAR images, the height error
Markus Schartel and Winfried Mayer are with the Endress+Hauser SE+Co. related to the reference height of the corresponding pixel can
KG, 79689 Maulburg, Germany. be calculated, and finally, the DEM can be determined [16].
Alexander Grathwohl and Christian Waldschmidt are with the Institute of
Microwave Engineering, Ulm University, 89081 Ulm, Germany. The big advantage of the approach shown here is that the
Digital Object Identifier 10.1109/LGRS.2021.3094165 UAV system with the radar is simultaneously used for both:
1558-0571 © 2021 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission.
See https://www.ieee.org/publications/rights/index.html for more information.

Authorized licensed use limited to: University of Waterloo. Downloaded on January 14,2024 at 08:31:21 UTC from IEEE Xplore. Restrictions apply.
4014505 IEEE GEOSCIENCE AND REMOTE SENSING LETTERS, VOL. 19, 2022

TABLE I
S YSTEM PARAMETER OF THE FMCW R ADAR

InSAR for DEM generation and ground penetrating synthetic


aperture radar (GPSAR) for subsurface sensing. During the
subsequent processing, first the DEM is calculated, which is
then used in the GPSAR processing as an input parameter.
The InSAR-based DEM calculation is designed in a way that
any trajectory can be flown. The DEM is generated with a Fig. 2. Schematic of the geometry.
neighboring frequency band, which in case of geometrically
complex radar signal Sn = F {Sn (t)} multiplied by the
undefined scattering centers leads to a reliable result of the
corresponding phase hypothesis φ [19]
air–ground interface for the GPSAR.
It is evident that the detection of buried mines strongly 
N

depends on an accurate DEM. With our remote sensing system P(x 0 , y0 , z 0 ) = Sn exp(− j φx0 (n),y0 (n),z0 (n) ). (2)
both, an accurate DEM as well as the detection of buried n=1

objects can be performed and evaluated simultaneously. In the case of the flat horizontal focus plane used for
The frequency-modulated continuous-wave (FMCW) radar InSAR, as shown in Fig. 2, the height z = z 0 is the same
used in this publication is published in [17]. The frequency for each computed pixel. The phase hypothesis φ for the pixel
band 1, which spans from 1 to 4 GHz, is used for the P(x 0 , y0 , z 0 ) is calculated with the start frequency f 0 of the
GPSAR images and the frequency band 2, which operates from radar and the distance rh from the phase center of the antenna
6 to 9 GHz, is used for the InSAR images. The broadband to the pixel. In Fig. 2, A1 is the transmitting antenna for the
antennas are shown in [18]. Fig. 1 shows the UAV system receiving antennas A1 and A2
equipped with the radar sensor, extended by one transmitting 2π 
and one receiving antenna, compared to the UAV system φx0 ,y0 ,z0 = 2 f0 (x A − x 0 )2 + (y A − y0 )2 + (z A − z 0 )2 .
c0
shown in [17]. The localization based on Kalman filtering is (3)
presented in [12] with the difference that an RTK GNSS with
centimeter-level accuracy is used instead of a total station. The The SAR image of the second receiving antenna is focused
baseline distance between the two receiving antennas is 35 cm. onto the same focus plane z = z 0 . Thus, for the pixel
The radar parameters are given in Table I. P(x 0 , y0 , z 0 ), the two SAR images are calculated as follows:
This letter is organized as follows. Section II briefly   

describes the SAR signal processing, interferogram generation, P1 (x 0 , y0 , z 0 ) = |Sn1 | exp j f 0 (2rh1 (n) − 2rm1 (n))
c
and DEM reconstruction. In Section III, the algorithms are

n

applied to measurements and a DEM is generated. This DEM 2π
P2 (x 0 , y0 , z 0 ) = |Sn2 | exp j f 0 (rh2 (n) + rh1 (n)
is used for subsurface processing with increased focusability c
of buried objects. Finally, Section IV gives a short conclusion.
n

− rm1 (n) − rm2 (n)) . (4)
II. I N SAR S IGNAL P ROCESSING
With the backprojection algorithm used, a mean-free devi-
The received FMCW radar signal S for a single measure- ation of the antenna position and the resulting error of the
ment is a sum of the received reflections i in the area, which phase hypothesis Δφx0 ,y0 ,z0 are averaged out by integration.
is illuminated by the antennas The interferogram, which contains the phase difference of
N   the two SAR images, is calculated from the argument of their

Sn (t) = Ai exp j f 0 rmi t . (1) conjugate-complex multiplication
c0  
i=1 ϕint = arg P1 · P2∗ . (5)
The amplitude Ai of the target depends on the radar cross
From this relation, the following expression for the phase
section (RCS), the antenna pattern, and the distance to the
difference of the two SAR images and a height error Δh
target. The phase of the reflection is given by the start
(Δh is the difference between the focus height z 0 and the
frequency f0 of the radar and the measured distance rmi . The
real/measured height) can be derived with d being the baseline
distance rmi equals the distance between the scattering center
between the two antennas
at (x mi , ymi , z mi ) and the phase center of the antenna Ai at
2π d
(x Ai , y Ai , z Ai ). ϕint ≈ −h f0 cos(α). (6)
For SAR-image generation, a backprojection algorithm is c0 rm1
used. Simplified for a pixel Pi at the location (x 0 , y0 , z 0 ), Thus, we have to find the height z 0 resulting in a phase
this algorithm integrates the received and range-compressed difference with ϕint = 0.

Authorized licensed use limited to: University of Waterloo. Downloaded on January 14,2024 at 08:31:21 UTC from IEEE Xplore. Restrictions apply.
BURR et al.: UAV-BORNE FMCW InSAR FOR FOCUSING BURIED OBJECTS 4014505

Fig. 4. Picture of the measurement field before the holes for the can lids
are closed. The diameter of the can lid is 5 cm. The perspective points to the
Fig. 3. Phase difference between the two SAR images as a function of height east.
error and distance.

Due to the relation of the wavelength with the periodicity


of φ, the sensitivity of a height error h is higher at the higher
frequency band. The ambiguity free range of the height error is
limited by the interval of ϕint ∈ [−π, π]. For the setup shown
in Fig. 1, the numerically calculated phase difference between
the two receiving antennas for the frequency band 2 with a
flight altitude of h 0 = 4 m is shown in Fig. 3.
The y  -axis represents the range along the focus plane.
The height h = 0 m corresponds to the focus plane z = z 0
of the two SAR images. For a scattering center above or
below the focus plane, a phase difference occurs as a function
of y  with a linear relation regarding (6). For the detection
of buried objects, a circular trajectory is preferred. The 2-D
aperture results in an image that is focused only in a very Fig. 5. Phase difference of two focused SAR images. The trajectory is a
circle with an altitude of 4 m.
small depth range, making clutter and targets more easily
distinguishable and locatable [12]. As a result, a scattering
center is illuminated from different angles and distances on A. InSAR Processing
the circular path. The conversion of the phase difference in
The height h 0 of the two SAR images for each receiving
the interferogram into a height error with respect to the focus
antenna is given as the mean value of the altimeter mea-
plane deviates outside the center of the trajectory from (6).
surement. The two SAR images for the interferogram are
Calculating the interferogram for several focus heights, a linear
calculated using the upper frequency band. The cell size for
relation between focus height and phase difference is obtained.
each pixel is 1 × 1 cm, which corresponds to the ground
Using the least squares approach, a linear approximation
resolution. The interferogram is calculated according to (5)
defines the point for ϕint = 0 rad (z 0 ), defining z 0 as real
and is shown in Fig. 5 after applying a 2-D median filter with
height in the DEM.
a window size of 4 × 4 cm. The y-axis is pointing north.
The positive phase trend toward the top-right edge of the
III. V ERIFICATION AND R ESULTS interferogram corresponds to the negative slope of the terrain.
To verify the overall system including the algorithms, mea- A phase wrap from π to −π can also be seen in the top-right
surements are performed over a meadow with an inclining edge of the interferogram. In addition, the top view of the
profile, as shown in Fig. 4. Can lids with a diameter of 5 cm trajectory is shown in the blue curve.
are buried approximately 5 cm below the grass surface. In In the following, the two SAR images are calculated for a
the following, the can lid in the center is marked with can focus height starting from −40 to +40 cm related to the global
lid 1 and the can lid in the top-left edge of the picture with height H0 with an increment of 1 cm. The resulting phase
can lid 2. After burying the targets, the soil was compacted difference profile at the four circles in Fig. 5 is shown in Fig. 6.
and flattened with a roller. The setup is surrounded by four The colored curves show the linear relation between phase
reference reflectors lying on the surface. The measurements difference and focus height for different positions. The black
are performed by flying in a circle concentrically to the buried curves show the polynomial fits through the phase difference
objects at a height of 4 m above ground level. The diameter of profile as described in Section II. The point of intersection
the circle is 15 m. The control is automated and uses the UAV ϕint = 0 rad corresponds to the measured ground level of the
GPS for (x,y) control. For altitude control, a 24-GHz radar is pixel. It can be seen that the slope at the different positions
fused with a lidar altimeter. After takeoff, the landing gear is is not constant. To verify the calculated height, single points
automatically raised to give the antennas a clear view of the of the surface profile are sampled with an RTK GPS. Based
ground. on these reference points, a surface is interpolated over the

Authorized licensed use limited to: University of Waterloo. Downloaded on January 14,2024 at 08:31:21 UTC from IEEE Xplore. Restrictions apply.
4014505 IEEE GEOSCIENCE AND REMOTE SENSING LETTERS, VOL. 19, 2022

Fig. 8. SAR image for the buried can lids. (a) SAR image for can lid 1 at
z = −7 cm. (b) SAR image for can lid 2 at z = −6 cm.
Fig. 6. Phase difference profile with corresponding line fit for four different
positions.

Fig. 9. Depth profiles of amplitude with (—) and without (- - -) DEM.

image of the region at the depth, in which can lids 1 and 2 are
buried.
The SAR image is normalized to the maximum, which
Fig. 7. DEM from the interferogram. corresponds to the reference reflector on the surface with
an RCS of −16 dBsm. With an amplitude of about −18 dB
area of the test field. The corresponding height for the colored compared to the reference reflector, the buried can lid is very
phase profiles in Fig. 6 is marked with *. It can be seen that clearly visible. The surrounding clutter level is significantly
the calculated altitude with the fitted line agrees very well lower. Fig. 9 shows the depth profile for the pixels in which
with the determined altitude using the RTK GPS. In Fig. 7, the two can lids are buried.
the calculated DEM of the test setup is shown. The focus height z in Fig. 9 is related to the respective
The local ground height z is related to the global height h 0 . height of the DEM. The DEM height is in the area of can lid 1
The declining profile of the terrain toward the x-axis is visible. −2 cm and in the area of can lid 2 −22 cm. The solid line
To verify the calculated DEM, the difference between the represents the amplitude versus focus depth. The amplitude of
DEM and the interpolated RTK GPS surface is calculated. the SAR image has a maximum below the air–ground interface
The standard deviation over the whole DEM is less than 4 cm. for the pixel in which the two cans are buried. For comparison,
With residual systematic localization errors, a relative height a second SAR image was calculated in which the DEM height
error over the DEM of 2.3 cm can be achieved. for can lid 2 corresponds to the local level of can lid 1 and vice
versa. The result is shown with dashed lines in Fig. 9. It can
be seen that the peak is shifted for can lid 1 by 5 cm to the
B. GPSAR Processing position z = −2 cm and the amplitude is reduced by 10 dB.
The lower frequency band is used to focus the objects below The peak of can lid 2 is shifted by −7 cm to position z =
the surface. As an extension of the backprojection algorithm, −13 cm and the amplitude is reduced by 5 dB. The reduction
the law of refraction as well as the DEM is considered [13]. of the amplitude results from the fact that the backscattered
For a focus plane below the DEM, the permittivity is consid- energy cannot be focused. The impact of the ground height
ered. The value of the permittivity is set to the empirical value error on the target peak amplitude is shown in Fig. 10.
r = 8. The focus is swept from −15 to +5 cm, relative to For this purpose, the DEM height z, which refers to the
the air-to-ground interface, with a step size of 1 cm. The pixel global height h 0 , is swept in the range from −32 to 4 cm. From
size for the GPSAR image is 1 × 1 cm. Fig. 8 shows the SAR the resulting depth profile for can lids 1 and 2, the maximum

Authorized licensed use limited to: University of Waterloo. Downloaded on January 14,2024 at 08:31:21 UTC from IEEE Xplore. Restrictions apply.
BURR et al.: UAV-BORNE FMCW InSAR FOR FOCUSING BURIED OBJECTS 4014505

R EFERENCES
[1] C. Castiblanco, J. Gómez, I. Mondragón, C. Parra, and J. Colorado, “Air
drones for explosive landmines detection,” in Proc. 1st Iberian Robot.
Conf. (ROBOT), vol. 253, Jan. 2014, pp. 107–114.
[2] L.-S. Yoo, J.-H. Lee, S.-H. Ko, S.-K. Jung, S.-H. Lee, and Y.-K. Lee,
“A drone fitted with a magnetometer detects landmines,” IEEE Geosci.
Remote Sens. Lett., vol. 17, no. 12, pp. 2035–2039, Dec. 2020.
[3] D. Daniels, Ground Penetrating Radar (Electromagnetics and Radar
Series), no. 1. London, U.K.: Institution of Engineering and Technology,
2004.
[4] W. Wiesbeck and C. Fischer, “New SAR system configuration for
detection of buried objects,” in Proc. IEEE Int. Geosci. Remote Sens.
Symp. (IGARSS), vol. 5, Jun./Jul. 1999, pp. 2471–2473.
[5] M. R. P. Cerquera, J. D. C. Montaño, I. Mondragón, and H. Canbolat,
Fig. 10. Amplitude for the buried can lid by varying the local ground height. “UAV for landmine detection using SDR-based GPR technology,” in
Robots Operating in Hazardous Environments. Rijeka, Croatia: InTech,
as a function of the ground level is displayed. The maximum 2017.
amplitude for can lid 1 is at a local ground level of −2 cm. [6] G. Moussally, K. Breiter, and J. Rolig, “Wide-area landmine survey and
detection system,” in Proc. 10th Int. Conf. Grounds Penetrating Radar
The maximum for can lid 2 is at −19 cm. (GPR), 2004, pp. 693–696.
A phase error due to Δr or Δh would lead to a wrong [7] R. Burr, M. Schartel, P. Schmidt, W. Mayer, T. Walter, and
depth estimation and to a reduction of the object-to-clutter C. Waldschmidt, “Design and implementation of a FMCW GPR for
UAV-based mine detection,” in Proc. IEEE MTT-S Int. Conf. Microw.
ratio (see Figs. 9 and 10 ). Critical values for Δr and Δh can Intell. Mobility (ICMIM), Apr. 2018, pp. 1–4.
be estimated as follows: [8] M. G. Fernández et al., “Synthetic aperture radar imaging system
√ for landmine detection using a ground penetrating radar on board a
c0 c0 r unmanned aerial vehicle,” IEEE Access, vol. 6, pp. 45100–45112, 2018.
Δh  √ ; Δr  . (7)
2 f 0 (1 − r ) 2 f 0 rm soil [9] G. Fasano, A. Renga, A. R. Vetrella, G. Ludeno, I. Catapano, and
F. Soldovieri, “Proof of concept of micro-UAV-based radar imag-
This results for given system parameters in Δh  8 cm ing,” in Proc. Int. Conf. Unmanned Aircr. Syst. (ICUAS), Jun. 2017,
and Δr  4. With respect to (7), the height error is much pp. 1316–1323.
[10] S. Perna, F. Soldovieri, and M. Amin, “Editorial for special, issue radar
more critical than a realistic permittivity error. Thus, for a imaging in challenging scenarios from smart and flexible platforms,”
robust detection of buried targets, a DEM with accuracy in MDPI, Basel, Switzerland, Tech. Rep., 2020, p. 1272, vol. 12, doi:
the centimeter range is required. 10.3390/rs12081272.
[11] M. Schartel, R. Burr, W. Mayer, N. Docci, and C. Waldschmidt, “UAV-
based ground penetrating synthetic aperture radar,” in Proc. IEEE MTT-S
Int. Conf. Microw. Intell. Mobility (ICMIM), Apr. 2018, pp. 1–4.
IV. C ONCLUSION [12] M. Schartel, R. Burr, R. Bähnemann, W. Mayer, and C. Waldschmidt,
In this work, a UAV-based multiband SAR system was “An experimental study on airborne landmine detection using a circular
synthetic aperture radar,” 2020, arXiv:2005.02600. [Online]. Available:
presented, which combines InSAR with GPSAR. Both the data http://arxiv.org/abs/2005.02600
for InSAR and the data for the GPSAR are simultaneously [13] A. Heinzel et al., “Focusing methods for ground penetrating MIMO
recorded with a single FMCW radar supporting two frequency SAR imaging within half-spaces of different permittivity,” in Proc. 11th
Eur. Conf. Synth. Aperture Radar (EUSAR), 2016, pp. 1–5.
bands. It is shown that a reliable DEM can be calculated on the [14] H. A. Zebker and R. M. Goldstein, “Topographic mapping from interfer-
basis of a circular trajectory. The standard deviation between ometric synthetic aperture radar observations,” J. Geophys. Res., vol. 91,
the InSAR DEM compared to the reference points recorded no. 5, pp. 4993–4999, 1986. [Online]. Available: https://agupubs.
onlinelibrary.wiley.com/doi/abs/10.1029/JB091iB05p04993
with the RTK GPS is only a few centimeters. The great advan- [15] R. Bamler and P. Hartl, “Synthetic aperture radar interferometry,” Inverse
tage of the InSAR DEM is that it can be calculated contactless Problems, vol. 14, no. 4, pp. R1–R54, 1998. [Online]. Available:
over the area required for the GPSAR. The calculated DEM http://stacks.iop.org/0266-5611/14/R1
[16] L. C. Graham, “Synthetic interferometer radar for topographic mapping,”
is successfully used as input parameter for the GPSAR as Proc. IEEE, vol. 62, no. 6, pp. 763–768, Jun. 1974.
air–ground interface. The GPSAR subsurface processing is [17] R. Burr, M. Schartel, W. Mayer, T. Walter, and C. Waldschmidt,
improved by an accurate air-to-ground interface obtained from “A broadband UAV-based FMCW GPR and the influence of vege-
tation,” in Proc. 12th German Microw. Conf. (GeMiC), Mar. 2019,
the InSAR DEM. The improvement enabled a clear focusing pp. 210–213.
of can lids with a diameter of 5 cm buried in grass covered [18] R. Burr, M. Schartel, W. Mayer, T. Walter, and C. Waldschmidt,
soil. The target response over depth of the can lids versus depth “Lightweight broadband antennas for UAV based GPR sensors,” in Proc.
15th Eur. Radar Conf. (EuRAD), Sep. 2018, pp. 245–248.
increased by 10 dB compared to processing with unknown air- [19] M. Soumekh, Synthetic Aperture Radar Signal Processing With
to-ground interface and thus enables detection. MATLAB Algorithms. Hoboken, NJ, USA: Wiley, 1999.

Authorized licensed use limited to: University of Waterloo. Downloaded on January 14,2024 at 08:31:21 UTC from IEEE Xplore. Restrictions apply.

You might also like