Array-Based Ground Penetrating Synthetic Aperture Radar on Board an Unmanned Aerial Vehicle for Enhanced Buried Threats Detection

The usage of unmanned aerial vehicles (UAVs)-based ground penetrating radar (GPR) systems has gained interest over the last years thanks to advantages over ground-based systems, such as contactless inspection and capability, to reach difficult-to-access areas. The former is of paramount importance concerning the detection of buried threats, such as improvised explosive devices (IEDs) and landmines. Current state-of-the-art UAV-based GPR systems are able to provide centimeter-level resolution thanks to the use of GPR-synthetic aperture radar (SAR) processing techniques. One of the challenges to keep improving these systems is the scanning throughput, that is, the area that can be scanned in a given time. This contribution presents an array-based GPR-SAR system for subsurface imaging, aiming at maximizing the scanning throughput without jeopardizing the imaging capabilities of the system. First, the antenna array is mounted on a portable setup to evaluate its performance and imaging capabilities. Next, the antenna array is integrated into the UAV platform, and the UAV-based GPR-SAR system with the array is tested in realistic scenarios with different kinds of buried targets. Results show that the scanning throughput is significantly improved and, furthermore, the coherent combination of all transmitting–receiving channels of the array provides enhanced detection capabilities.

In the last years, GPR systems and unmanned aerial vehicles (UAVs) have been synergized enabling the contactless inspection of the subsoil [9], [10]. This is of special interest in the case of landmine and IED detection, where groundbased GPR systems can accidentally detonate buried threats. In this context, the capability of UAV-based GPR systems to detect either metallic and nonmetallic targets enables to overcome the limitation of airborne prototypes making use of metal detectors [11] or magnetometers [12], which are constrained to the detection of metallic targets. These systems also overcome the limitations of hyperspectral systems for landmine detection [13], which are mainly devoted to the detection of partially buried or shallow targets.
Reviews of the state-of-the-art UAV-based GPR systems have been presented in [14], [15], [16], describing the main features and challenges faced by these systems. Although most of the UAV-based GPR system prototypes developed in the last years have been devoted to detecting landmines and IEDs, they have been also introduced in other application areas like snow moisture and ice thickness measurement [17], [18], soil moisture characterization [19], or search and rescue missions [20], [21], among others.
Some of the current state-of-the-art UAV-based GPR prototypes provide 3-D images of the underground with centimeterlevel resolution [22], [23], [24], [25], [26]. To achieve this goal, these prototypes are equipped with cm-level accuracy positioning and geo-referring subsystems integrating laser rangefinder or LIDAR sensors and global navigation satellite systems (GNSS)-real-time kinematics (RTK) receivers. They also rely on advanced GPR-synthetic aperture radar (SAR) processing algorithms that overcome the issues intrinsic to UAV-based systems, such as nonuniform acquisition grids or strong clutter contributions from the air-soil interface, which pose significant challenges for the detection of shallow targets. A description of some of these processing techniques, which are devoted to enhancing the quality of the obtained GPR images, reducing the presence of clutter, and improving the detection capabilities of the system, is presented in [27]. This   One of the main challenges to improving the performance of UAV-based GPR systems is the increase of the scanning throughput, that is, the area (or volume) that can be inspected in a given time. Current prototypes, like the ones presented in [23] and [24], are able to scan up to 35 m 2 in a single flight. Wide-area GPR surveys using airborne-based systems have been conducted with manned aircrafts (e.g., planes and helicopters) [28], [29], but at the expense of increasing the cost, and complexity with respect to UAV-based systems. Increasing the scanning throughput would require either: 1) a radar module capable of performing faster acquisitions, so that the survey speed of the UAV can be increased, or 2) the use of an antenna configuration that enables to inspect a greater area per each of the sweeps performed by the UAV. The former strategy relies on the use of more complex (and thus, more expensive) hardware. Furthermore, it could also imply bulkier and/or heavier radar modules. The latter strategy has been adopted in several ground-based GPR systems by integrating an antenna array to increase the width of the scanned strip (see Table I). Nevertheless, the implementation of these solutions is not straightforward in the context of UAVs. Ground-based GPR systems usually rely on large vehicles and trucks to carry large antenna arrangements (e.g., [30], [31], [32], [33], [34]). However, the integration of a GPR within a UAV is heavily constraint by the available space and weight in the platform, which directly impact the maximum flight time. Moreover, the payload integrated within the UAV, especially the GPR antennas, severely affects its flight stability and its maneuverability. Therefore, the design of an antenna array for UAV-based GPR systems must take into account these issues to ensure that the scanning throughput is effectively improved without significantly affecting the flight time, stability, and maneuverability.
In addition, another challenge faced by array-based GPR systems that cannot be overlooked is the time required to perform a complete sweep of all the available radiofrequency channels. More specifically, if these systems are mounted on moving platforms (e.g., ground vehicles or trolleys), the time required to perform an acquisition would limit the moving speed of the system (e.g., in [35] the system stops, performs a full acquisition and then moves again). Otherwise, several zones of the surveyed area could be missed or, in the case of combining measurements acquired at different positions, subsampling would take place, leading to a degradation of the system performance.
Recently, there have been some attempts to overcome these issues. In this regard, the prototype presented in [36] implemented a linear array with one Tx and three Rx (separated 30, 90, and 120 cm from the transmitter) though at a very early stage (single roundtrip sweep of about 12 m). Also, the system proposed in [22] takes advantage of the two receiving channels of the employed radar module, enabling a one Tx and two Rx configuration. In addition, the use of polarimetric SAR measurements (one Tx and two Rx) was considered in [37], aiming at improving the detection capabilities of the system. Their main features are summarized in Table I, along with those of other ground-based GPR systems using antenna arrays. Most of these systems are conceived to be mounted on a ground vehicle or at an early stage of development. Thus, the weight and size of the different components are not a major concern, and laboratory equipment is frequently used. Most of them make use of (one or more) vector network analyzers (VNAs) working in a frequency band that, in general, ranges from a few hundreds of MHz up to a few GHz. In addition, practically all of these systems have, at least, 1500-MHz bandwidth (BW), resulting in a free-space range resolution ( r ) equal or better than r = c/(2BW) = 10 cm. To achieve this large fractional bandwidth, Vivaldi antennas are usually selected thanks to their good trade-off in terms of compact size and wideband frequency response.
The contributions listed in Table I employ different architectures. For example, [30], [38], and [39] use several Tx and Rx antenna pairs that are located at different distances with respect to a common midpoint. This is called the common midpoint (CMP) technique, and the scanned area is illuminated under different incident angles. Moreover, the combination of the measurements from the different Tx-Rx pairs improves the signal-to-clutter ratio. In the case of [40] and [41], a single Tx is used, so the midpoint position changes when switching between Rx antennas. In contrast, a side-looking architecture is considered in [31], [32], and [33]. This system benefits from a significant separation between Txs and Rxs, which are also at different heights, to obtain 3-D images. This is enabled by a large structure elevated with a lifter. Such kind of antenna arrangement cannot be integrated within a UAV, making it necessary to resort to other antenna configurations and several sweeps at different heights to achieve 3-D images.
This contribution focuses on a down-looking GPR (DLGPR) architecture [22], [23], in which adjacent samples should not be separated greater than half a wavelength (λ min /2) at the highest working frequency to avoid aliasing in the resulting SAR images. Therefore, the scanning throughput is severely limited by the large number of along-track sweeps required to scan a particular area, and minimizing this number is essential to facilitate the usage of this kind of system in real scenarios (drastically reducing the time needed to scan a given area).
In particular, this contribution presents a novel UAV-based GPR-SAR system, which makes use of an antenna array to significantly improve the scanning throughput without jeopardizing its performance in terms of detection capabilities. For this purpose, a particular switching scheme of the transmitting and receiving channels, together with a specific flight path configuration has been conceived. To the best of our knowledge, this is the first time that an array with multiple transmitting channels and multiple receiving channels is implemented on board a UAV for a GPR-SAR system. Furthermore, extensive experimental validation has been carried out, showing that the proposed approach results not only in a significant increase in the scanning throughput but also in enhanced detection capabilities. The latter is due to the fact that the proposed multistatic system allows illuminating the scanned area from different angles and the combination of radar images from all the Tx-Rx channels results in a higher signal-to-clutter ratio.
This article is structured as follows. Section II describes the architecture of the UAV-based GPR system, as well as the GPR-SAR processing technique to retrieve the 3-D images of the subsoil. In addition, this section also includes a theoretical analysis of the system and some practical considerations. On-ground validation with a portable setup that emulates the behavior of the payload prior to its integration in the UAV is presented in Section III. Section IV shows the results corresponding to in-flight tests of the array. Finally, conclusions will be drawn in Section V.

II. SYSTEM DESCRIPTION
As aforementioned, this contribution focuses on the DLGPR architecture presented in [22] and [23]. Under such configuration, the UAV follows a predefined zigzag path to cover the area of interest, as shown in https://youtu.be/HDUwgka8Dns. However, as the GPR system relies on SAR techniques to enable high-resolution images of the subsoil, consecutive sweeps of the UAV should not be separated greater than half a wavelength (λ min /2) at the highest working frequency, which, in this case, is 5 cm as the highest frequency considered for processing is 3 GHz. This entails a significant constraint in the required across-track spacing to avoid aliasing, which in turn results in a large number of sweeps to inspect a given area, limiting the scanning throughput of the system. Furthermore, it is worth noting that at the end of each along-track sweep, the UAV has to slow down, perform a lateral displacement of a few cm (e.g., not greater than λ min /2 = 5 cm), and then increase the speed again, resulting in additional flight time. Consequently, minimizing the number of along-track sweeps to scan a particular area would result in significant flight time savings.
Although the elements in an array are usually separated λ min /2 to avoid the presence of grating lobes, the information from multiple acquisition positions can be combined so that the interelement spacing requirement can be relaxed. In this contribution it is proposed to exploit the UAV movement to achieve a proper sampling, enabling the use of SAR processing techniques, whilst the separation between Tx and Rx antennas can be increased. This facilitates their integration within the UAV platform and enlarges the area surveyed per flight.

A. Architecture
The antenna array, including the switches employed to select the different Tx-Rx combinations, will be integrated into the UAV platform employed in [22] and [23], which provides a maximum take-off weight (MTOW) of 11 kg. The UAV-based GPR system consists of the following subsystems: 1) Flight control subsystem, which comprises a microcomputer (also acting as the UAV flight controller). 2) Positioning and geo-referring subsystem, formed by the standard devices and sensors used by most UAVs (inertial measurement unit (IMU), barometer, compass, and a conventional GNSS receiver), and a laser rangefinder [46] and a triple-band multiconstellation GNSS-RTK receiver [47] to enable the high accuracy navigation and localization required to apply GPR-SAR processing. 3) Communications subsystem, which includes the radio controller (R/C) transmitter and receiver modules, and an ad hoc wireless network to establish connectivity between the UAV and the ground station (processing laptop) to transfer measurement data and position corrections. 4) Radar subsystem, composed by a UWB radar module [48] with one transmitting channel and two receiving channels, and the radar antennas. Concerning the positioning and geo-referring subsystem, it is worth highlighting the importance of geo-referring the radar measurements with high accuracy to retrieve wellfocused radar images. In this regard, the RTK integrated into the UAV is configured to track the following constellations and signals: Global Positioning System (GPS) L1, L2, L2C, and L5; GLONASS G1, G2, and G3; and Galileo E1 and E5. The accuracy of this receiver is 5 mm+0.5 ppm × baseline in the horizontal plane and 10 mm + 0.8 ppm × baseline in the vertical direction (where the baseline is the distance between the RTK receiver and the RTK base station) [47]. Both a network of RTK base stations and a dedicated RTK base station have been used to provide corrections, although the latter approach is preferred in order to minimize the baseline and avoid the dependence on network coverage. In this case, our dedicated RTK base station is deployed next to the ground control station (which, in the tests shown in this contribution, is located less than 50 m away from the RTK receiver on board the UAV).
Regarding the radar subsystem, it should be noted that the UWB radar mounted on board the UAV transmits a pseudorandom sequence. Among other advantages, this technology offers a lower transmit power, which helps to reduce possible interferences. The fact that the radar antennas are pointing downward also helps to mitigate them. Nevertheless, due to the importance of high-accuracy GNSS-based positioning, a dual-notch filter at the GNSS frequency bands is connected at the output port of the radar to further reduce the risk of interference into the GNSS receivers. No interferences have been detected between the radar and any other systems mounted on board the prototype (GNSS receivers, wireless network, and R/C transceiver).
In this contribution, an array of antennas is employed to increase the overall scanning throughput of the GPR-SAR system. The design and initial assembly of the antenna array, which is not within the aim and scope of this contribution, were conducted by the research group Teoría do Sinal e Comunicacions, Universidade de Vigo, Vigo, Spain, which had previously designed another array for GPR applications [49]. The following requirements and restrictions are of special relevance: 1) Maximum weight of 1200 g. This value results from the subtraction of the MTOW and the sum of the weight of all the UAV components except for the GPR antennas and switches.
2) The maximum length of the array should be 80 cm approximately. A larger array could have an impact on the in-flight performance of the UAV (e.g., less stability).
3) The working frequency band of previous prototypes [22] should be maintained. That means that antennas and switches should be able to operate, approximately from 600 MHz to 6 GHz. Given these design specifications, a seven-element antenna array was provided, consisting of two subarrays, with three Tx and four Rx antennas, respectively. This layout contributes to balancing the weight of the payload, thus resulting in a more stable flight. Besides, the choice of three Tx and four Rx results in 3 × 4 = 12 virtual pairs of Tx-Rx antennas, where the placement of each virtual Tx-Rx antenna is the midpoint between the Tx and the Rx. The UWB Vivaldi antennas described in [50] were selected as a trade-off between the size and weight limitations and radiation performance. It must be pointed out that, in the study conducted in [23], these antennas were found to provide acceptable detection capabilities when compared with larger UWB Vivaldi antennas.
A picture of the antenna array is shown in Fig. 1(a), whereas the position of the antennas and the spacing between them is depicted in Fig. 1(b). The overall weight of the antenna array, including radiofrequency switches and cables, was 1150 g. As observed in Fig. 1, a gap was left between the Tx and Rx subarrays so that the landing gear of the UAV could be deployed and retracted. One radiofrequency switch is used to commute between the three-element Tx subarray, and a second switch controls the four-element Rx subarray. The latter has two output ports to take advantage of the two receiving channels of the employed radar module (i.e., two receiving antennas are activated at a time). The CMP of each Tx-Rx pair is also depicted in Fig. 1(b) with gray dots. The time required by the GPR subsystem to acquire the measurements corresponding to the three Tx × four Rx combinations is 280 ms.
Taking these features into account, the switching scheme between the antenna array elements has been designed to ensure that the sampling of the surveyed area is as uniform as possible while fulfilling the Nyquist sampling rate. This switching scheme is conditioned by the across-track spacing between two consecutive along-track sweeps of the array and also by the flight speed of the UAV (as previously explained, the maximum UAV flight speed is limited by the maximum spacing between consecutive samples in the along-track direction, i.e., λ min /2). In particular, the switching scheme has six stages: first, the receiving switch activates Rx1 and Rx3 and the transmitter switch cycles activating the different transmitters (Tx1, Tx2, and Tx3); then, the receiving switch activates Rx2 and Rx4 and the transmitter switch cycles again through all transmitters. It is worth noting that although some Tx-Rx pairs have the same CMP, in practice, there is not redundant information as the acquisitions are performed at different positions (e.g., channel Tx1-Rx2 and channel Tx2-Rx1 are not measured simultaneously). Furthermore, even if all the channels were acquired at the same UAV position (simultaneously), Tx-Rx channels with the same CMP would provide complementary information (due to the different illumination angles). After an along-track sweep is completed, the UAV moves 20 cm across-track to continue the survey of the area under inspection. For the sake of clarity, an animation illustrating how the scanning is conducted with the antenna array can be watched in the supplemental multimedia file "AnimArraySweep." The animation also illustrates how the proposed scheme results in a significantly uniform sampling of the investigation domain. A picture of the implemented UAV prototype with the antenna array is shown in Fig. 2. As can be seen, a 3-D-printed piece was designed to integrate the array structure in the UAV platform and also to give stability to the overall structure.

B. GPR-SAR Imaging
Some of the existing UAV-based GPR systems make use of SAR processing techniques to coherently combine the measurements collected over the area of interest. Like physical apertures, spacing between radar measurements must not be greater than λ min /2 (e.g., 5 cm for a working frequency of 3 GHz). Besides, geo-referring uncertainty should be small enough so that it does not have a significant impact on GPR-SAR images. In [22], a masking procedure was introduced to define the set of measurements that can be coherently combined given a certain value of geo-referring uncertainty.
In the UAV-based GPR system presented in [22] and [23], radar measurements are processed using the delay-and-sum (DAS) or backprojection method. It consists of adding all measurements coherently at one focal point of the imaging or investigation domain, then repeating this procedure for all the points of the investigation domain [51]. Under free-space considerations (that is, the relative permittivity of the soil, ε r , is assumed to be 1), the reflectivity on a point located underground, ρ(r ′ ) can be calculated as where E scatt (r n , f m ) is the scattered field measured at the nth measurement position, r n , and the mth discrete frequency, f m . The terms φ Tx and φ Rx account for the phase shift due to the wave propagation from the transmitter to the focusing point r ′ and from the focusing point to the receiver, respectively, and are given by and φ Rx,n,m = k 0,m ||r ′ − r Rx,n || 2 k 0,m is the free-space wavenumber at the mth discrete frequency, and r Tx and r Rx are the positions of the Tx and the Rx, respectively. This formulation can be modified to take into account the relative permittivity of the soil, as explained in [52]. The reason why the DAS technique is chosen for the processing of UAV-based GPR measurements is due to its capability of handling measurements that do not follow a regular acquisition pattern. Microwave tomography can also be used to process irregularly acquired measurements, but it can be computationally demanding, especially for inspecting large scenarios [53]. There are other methods more computationally efficient than DAS or tomography, like the phase shift migration (PSM) [54], [55], but they require the acquisition domain to be evenly sampled. Furthermore, the DAS algorithm can be directly adapted to the use of an antenna array taking into account the position of each Tx-Rx pair at each acquisition point. For this purpose, the position of the UAV is available whenever an acquisition with a Tx-Rx combination is performed.
Although the working frequency band of the Vivaldi antennas selected for the antenna array ranges from 1000 to 6000 MHz [50], only the 1000-3000-MHz frequency band is considered for processing. The reason behind this is that frequencies above 3000 MHz do not provide relevant information for subsurface imaging due to the fast attenuation in the subsoil, as discussed in [23] and [27]. Therefore, the free-space range or depth resolution of the presented system is r = 7.5 cm.
As mentioned earlier, the DAS algorithm can be modified by adopting a masking procedure, so that only those measurements in the vicinity of the point at which the reflectivity is computed (r ′ ) are considered for processing. In this method, called masked SAR, the size of the mask is given by the integration length (in both the across-track and along-track directions), which is the length along which the measurements can be coherently combined. In this regard, it must be noted that in real conditions the integration length is limited due to several factors, such as the size of the area illuminated by the antennas and the accuracy of the positioning system. For the proposed system, the integration length in the alongtrack axis has been set to L along-track = 2 m, in agreement with the projection of the antenna 3-dB main beamwidth at the central frequency considered (2 GHz). However, in the across-track axis, the integration length is further limited by the positioning accuracy as the geo-referring uncertainties are less correlated among different sweeps [22]. Based on the experience gathered with all our prototypes, the integration length across-track has been set to L across-track = 1 m to guarantee that the positioning uncertainties do not distort the radar images.
It should be noted that the GPR-SAR algorithm (either conventional DAS or masked SAR) is the most time-consuming step within the processing chain. The processing time depends mainly on the number of measurements, the number of voxels in which the investigation domain is discretized (which, in turn, is given by the size of the area under inspection and the voxels size), and the characteristics of the machine used to process the data. In this system, the measurements are processed on a conventional laptop with an i7-8750H CPU, Although the core algorithm in the processing of the radar measurements is the masked SAR method, complementary techniques must be also implemented in order to mitigate the clutter and enhance the image focusing. In particular, before the SAR processing, the following clutter mitigation techniques are applied: average subtraction, height correction, and distance-based singular value decomposition (SVD) filtering. Furthermore, after the SAR processing, the resulting radar images are equalized and co-registered to improve their quality. Finally, to further improve the signal-to-clutter ratio, the radar images of each Tx-Rx channel are coherently combined to obtain the final 3-D multichannel radar image.

C. Theoretical Analysis
Before conducting the experimental validation, the performance of the proposed array-based radar system was evaluated by computing the point spread function (PSF). This analysis has been carried out taking into account the conditions faced by the system when the radar is mounted on board the UAV. A 4 m × 4 m observation domain has been considered, in which measurements are gathered by performing a set of parallel along-track (y-axis) sweeps (with a separation between consecutive sweeps of 20 cm). The switching scheme described in Section II-A has also been adopted. In addition, as the radar module has two receiving ports, it has been assumed that at each acquisition position two different channels are measured. The acquisition positions are shown in Fig. 3(a), and a zoomed-in view depicting the position of the active transmitter at the different acquisition positions in a small area of the observation domain is depicted in Fig. 3(b). As can be seen, the scene is fully illuminated.
The height has been set at h array = 1.3 m, in agreement with the average height of the array during the UAV flights. It should be noted that the UAV is configured to fly at a constant distance from the soil surface (i.e., at a constant height). In particular, the UAV height is set at h = 1.5 m, as it is measured by the on board laser rangefinder (which is placed slightly closer to the UAV frame than the array). This height has been selected as a trade-off between the following factors. On the one hand, the lower the height is, the better the cross-range resolution is and the lower the attenuation losses are. On the other hand, UAV navigation is affected by an aerodynamic phenomenon called the ground effect when flying close to the soil surface [56], [57]. Besides this issue, a stand-off distance from the soil surface must also be kept to guarantee the safety of the system, so that the UAV can successfully avoid obstacles. Therefore, according to these navigation and safety requirements, the stand-off distance has been set to h = 1.5 m.
A point-scatterer target placed at the origin of the system of coordinates has been considered (i.e., at x = y = z = 0 m) and the reflectivity has been computed for each Tx-Rx combination using the DAS algorithm in an investigation domain of 2 m × 2 m × 2 m. The total number of measurements considered is 6762, and the investigation domain is discretized in voxels of 0.025 m × 0.025 m × 0.025 m size, resulting in 531 441 voxels. The processing time for the DAS algorithm is approximately 320 s. Nevertheless, it should be noted that this time can be significantly reduced if a coarser discretization of the investigation domain is considered and/or if the size of the investigation domain is reduced. For instance, if only a 2-D cut (e.g., z = 0 m) is considered in the reconstruction, the time is reduced to only slightly over 13 s.
The main cuts of the resulting PSFs for the Tx placed in the middle of the array (Tx2) and each Rx are shown in the left column of Fig. 4. The PSFs are almost identical for all these combinations, and the main differences are observed in the x-axis cut (y = z = 0 m, Fig. 4(a) are slightly different due to the different position of each Rx across the x-axis. In the right column of Fig. 4, the main cuts of the PSFs for each Tx (considering all the Rxs) have been depicted, together with the PSF obtained when combining all Tx-Rx channels. It can also be observed that the PSFs are almost identical for each Tx, with only some discrepancies in the x-axis cut (Fig. 4(d)), due to the different locations of each Tx across this axis. Furthermore, the sidelobes level is greatly reduced when all the Tx-Rx channels are combined. This is because when only one Tx is used the sampling is too coarse, whereas when all the Txs are considered the observation domain is fully sampled.
To further analyze the behavior of the array, two of the main 2-D cuts (the horizontal cut z = 0 m, and the acrosstrack cut y = 0 m) are shown in Fig. 5. When only one Tx and all Rx are considered (Fig. 5(a) and (d)), the sidelobes level is significantly higher than when all Tx-Rx combinations are considered (Fig. 5(c) and (f)), as expected from the conclusions extracted from Fig. 4(d).
Another configuration of interest is based on selecting the following six Tx-Rx channels: Tx1-Rx1, Tx1-Rx2, Tx2-Rx2, Tx2-Rx3, Tx3-Rx3, and Tx3-Rx4. As observed in Fig. 2(b), this corresponds to a multi-quasi-monostatic architecture and, besides, these channels have different CMPs. When only these six Tx-Rx channels are considered, the sidelobes level increases slightly (Fig. 5(b) and (e)), compared with using all Tx-Rx combinations (Fig. 5(c) and (f)). Nevertheless, this increase is below 3 dB. It should also be highlighted again that, although several channels produce the same CMPs (e.g., Tx1-Rx3, Tx2-Rx2, and Tx3-Rx1), they are not acquired simultaneously and the look angles are different, so they can provide complementary information. Furthermore, as already mentioned in previous studies [58], in antenna array architectures, such as the one proposed in this contribution, combining the radar images of the different channels benefits from the spatial diversity of the antennas, which in turn contributes to producing a decorrelation effect on the clutter.

1) On the Coherent Combination of All Tx-Rx Channels:
To improve the detection capabilities of the system, it is necessary to increase the signal-to-clutter ratio of the GPR images. For this purpose, the information from all the Tx-Rx channels is coherently combined. The coherence of the measurements can be degraded due to positioning and calibration errors. The former mainly depend on the uncertainty of the positioning sensors, whilst the latter can be mitigated with an appropriate calibration procedure.
As aforementioned, the radar module has two receiving channels, which are connected through a dual switch to the Rx antennas of the array. In addition, the Tx of the radar module is connected to the Tx antennas of the array using another switch. In order to calibrate the system, a first step consisting of retrieving the impulse response function of the two different radar channels was performed. Then, a calibration term was computed so that the two radar-receiving channels have the same amplitude levels and can be added coherently. To avoid the need to individually calibrate each Tx-Rx channel of the array, the same equal-length cables were used to connect each antenna of the array to the switches. In order to verify this approach, a test target was placed over the soil surface and the radar images obtained with the different channels were computed (considering the same time zero for all channels). Comparing these images, it was concluded that the target was located at its actual position regardless of the particular Tx-Rx channel considered and, thus, there was no need for an extra calibration term.
The impact of the calibration and of the coherent/incoherent sum of the channels can be observed by comparing the GPR-SAR images of a VS-1.6 plastic antitank landmine obtained from in-flight measurements (Fig. 6). As can be observed when no calibration is applied and the Tx-Rx channels are summed coherently ( Fig. 6(a)), the shape of the target is not well reconstructed and there is a significant amount of clutter. When an incoherent combination is performed (Fig. 6(b)), although the circular shape of the target is retrieved, the clutter level is significantly high. However, when the calibration is considered and all the Tx-Rx channels are coherently combined (Fig. 6(c)), not only the target is accurately reconstructed but also the clutter level is greatly reduced compared with the previous cases. Therefore, it can be concluded that the targetto-clutter ratio improves drastically.
2) Height Estimation From Radar Measurements: As aforementioned, the coherence of the measurements can be severely affected by positioning errors, especially in the range direction (i.e., z-axis in this case). In this system, the horizontal position, the height, and the orientation of the UAV are obtained from the GNSS-RTK receiver, the laser rangefinder, and the IMU, respectively. Then, the position of each antenna is computed taking into account the actual position of the UAV and the location of each antenna with respect to the UAV origin of coordinates. In the previous monostatic system [27], a technique to estimate the UAV height from the GPR measurements was proposed, showing that it allows retrieving better focused radar images. This technique has been extended to the arraybased architecture presented in this contribution. In particular, every time an acquisition is performed a height estimation is computed. As previously explained, two Tx-Rx channels are simultaneously measured at each acquisition position. The average height of each Tx-Rx combination can be estimated from the arrival time of each radar signal [27]. Assuming that the pth transmitter and the qth and lth receivers are active, the estimated heights, denoted as h Txp-Rxq and h Txp-Rxl , respectively, are given by Then, assuming the soil surface is locally planar and constant in the along-track direction, namely the y-axis, it is possible to define a local plane to model the soil surface at each measurement position. This plane contains the vectorsŷ and Fig. 7. Picture of the portable setup [59] and the targets selected for testing purposes: two 10 cm × 10 cm metallic plates placed 6 cm above the ground, a 16.5 cm × 16.5 cm metallic plate placed on the ground, and a plastic VS-1.6 antitank landmine buried 6 cm in the ground. r CMP,Txp-Rxq − r CMP,Txp-Rxl (where r CMP,Txp-Rxq is the position of the CMP between the pth Tx and the qth Rx, and r CMP,Txp-Rxl is the position of the CMP between the pth Tx and the lth Rx). Finally, evaluating this plane in the across-track coordinates (namely, the x-axis) of the pth Tx, and the qth and the lth Rx, it is possible to find their respective heights.

III. VALIDATION WITH A PORTABLE SETUP
The first step toward the integration of the antenna array in the UAV frame was the assessment of its imaging capabilities. For this goal, the array was mounted on a portable setup used to validate the UAV-based GPR-SAR subsystems before integrating them into the UAV platform. The setup consists of a plastic frame that has a sliding arm which includes a cart made of a plastic box attached to it. The sliding arm can be moved along-track using two ropes, and manually displaced in the across-track direction. As a result, it allows scanning an area of approximately 1 m × 1 m, as described in [59].
A picture of the array attached to the bottom side of the cart is shown in Fig. 7. The rest of the payload placed inside the cart replicates that of the UAV platform, and it consists of the UWB radar module, the GNSS-RTK receiver, the microcomputer, and a battery to power up all the devices.
As previously explained, the across-track spacing between consecutive along-track sweeps was set to 20 cm. This value was selected as a trade-off between minimizing along-track sweeps while ensuring a proper sampling of the surveyed area. At this point, it should be remarked that, even though the flight path performed by the UAV system follows a set of predefined waypoints, several factors inherent to the UAV operation (such as wind and the UAV dynamics) result in a nonperfect acquisition grid. This should be taken into account as deviations of the actual flight with respect to the original flight plan are to be expected. To illustrate the sampling distribution resulting from this scheme, six along-track sweeps were performed to inspect the 1 m × 1 m area under the portable setup, which is enclosed by a blue dash-dotted line in Fig. 8(a). Also in Fig. 8(a) the positions where acquisitions were performed during the scan are depicted in black dots. As can be observed, the different along-track sweeps deviate slightly from straight lines (as it would be the case during a UAV flight). The resulting sampling distribution, considering all Tx-Rx combinations, is shown in Fig. 8(b). In particular, the acquisition domain was discretized in 5 cm × 5 cm cells, i.e., λ min /2 × λ min /2. Thus, Fig. 8(b) shows a heatmap of the number of measurements performed within each cell.
In order to provide a clear sight of the sampled areas within the acquisition domain, a binary representation of the sampling distribution is depicted in Fig. 8(c). In particular, sampled cells are colored in yellow. As can be observed, the area that can be covered by the array with the portable setup (enclosed by a black and white line) is fully sampled. The measurements that lie outside the area under inspection (enclosed by a black and white line), were acquired due to the length of the array (crossrange dimension of the scanning scheme), which partly falls outside the area delimited by the frame of the portable setup when the array is at its edges. However, they are not relevant for this analysis as they are outside the area of interest. The total number of measurements considered is 474.
Three targets, depicted in Fig. 7, were placed in the scanned area: 1) two 10 cm × 10 cm metallic plates on top of a piece of foam to place them 6 cm above the air-soil interface; 2) a 16.5 cm × 16.5 cm metallic plate placed on the air-soil interface, and 3) a VS-1.6 antitank plastic landmine, buried 6 cm deep with dimensions 22 cm of diameter × 9.2 cm thick. The loamy soil of this scenario had an estimated relative permittivity of approximately ε r = 7.
The investigation domain has been defined as a volume of 1.5 m × 1.5 m × 0.3 m, and the voxels size is set to 0.03 m × 0.03 m × 0.01 m. This size has been found as a good tradeoff to facilitate the qualitative inspection of the resulting images (a finer discretization would require more computational resources, and a coarser discretization would make it more difficult to distinguish small details). As a result, the total number of voxels is 80 631, and the processing time of the GPR-SAR imaging algorithm is slightly less than 5 s.
The GPR-SAR imaging results obtained using the antenna array and considering the measurement samples depicted in Fig. 8 are shown in Fig. 9(a)-(c). As can be observed, the three targets are properly imaged. In particular, the z-cut shown in Fig. 9(a) corresponds to the two metallic plates located over foam, placed 6 cm above ground (profile enclosed with a black and white line). It should be noted that there is a second area with a high reflectivity level (around x ≈ 0.5 m and y ≈ 0.8 m), which comes from the second target, a bigger metallic plate placed right at the soil surface. The image of the z-cut corresponding to this second target is depicted in   Fig. 9(b), where the squared shape of the plate, which is the most reflective target of this scenario, can be observed. As can be seen, the dimensions of the high reflectivity area match the real size of the target (black and white line). Finally, the z-cut corresponding to the top interface of the last target, the VS-1.6 antitank plastic landmine, is shown in Fig. 9(c). In this case, as the target was buried in the ground, the SVD filtering technique described in [27] was applied to mitigate the clutter coming from the air-soil interface. Analogously to the previous two z-cuts corresponding to the other targets, the size of the plastic landmine is well-reconstructed (see the black and white line), and two high-reflectivity areas caused by the other two metallic targets can be observed around (x ≈ 0.5 m and y ≈ 0.8 m), and (x ≈ 0.6 m and y ≈ 0.1 m).
It should be noted that all targets are properly imaged, with more than 10-dB dynamic range over the clutter level, thus validating the imaging capabilities of the array-based GPR-SAR imaging system. In addition, it should be remarked that no echoes due to grating lobes are observed within the imaging domain, which is in agreement with the results depicted in Fig. 8(b) and (c) about the sampling of the investigation domain.
In terms of scanning time savings, only six along-track sweeps spaced 20 cm were needed to scan the area of interest while fulfilling the Nyquist sampling rate (5 cm spacing between consecutive samples). This means four times fewer along-track sweeps with respect to the GPR-SAR architecture that does not make use of the array, i.e., a 75% time reduction.

IV. IN-FLIGHT VALIDATION
Once the operation of the GPR-SAR imaging system with the antenna array was tested on the ground, then it was mounted in the UAV frame. The same UAV platform as in previous prototypes (see [22] or [23]), a DJI S1000+, was considered. The antenna array was mounted on a customized carbon fiber-made frame that was reinforced to minimize vibrations during in-flight operations. An ad hoc interface was designed and manufactured using 3-D printing to adapt the antenna array frame to the UAV.

A. Initial Tests
The first flight tests with the array mounted on the UAV were conducted in a scenario located at the airfield of the University of Oviedo, Gijón, Spain. The size of this scenario was 5.25-m along-track × 1.5-m across-track. Two metallic targets, a 1-m-long × 9-cm-wide metallic bar, and a 16.5 cm × 16.5 cm size metallic plate, were placed on the ground for testing the imaging capabilities of the array-based GPR-SAR UAV system. A picture of the prototype while performing the flight as well as the two targets is shown in Fig. 10.  The same scanning strategy adopted in Section III was followed in the flight conducted with the UAV: the predefined flight path consisted of a zigzag pattern where along-track sweeps were spaced 20 cm. The position of the performed along-track sweeps is plotted in Fig. 11(a). As can be observed, the actual flight path deviates slightly from the predefined one due to the presence of wind gusts during the flight and the dynamics of the UAV. As a consequence, the spacing between consecutive sweeps is not exactly 20 cm. Nonetheless, as can be seen in Fig. 11(b) and (c), almost all the 5 cm × 5 cm cells have, at least, one measurement.
Flight tests were conducted at a speed of 50 cm/s, at a height of 1.5 m above the ground, resulting in 2522 measurements to be processed. As observed in Figs. 11(b) and (c), there are almost no gaps in the along-track direction, thus fulfilling the Nyquist sampling rate. Concerning the performance and operation of the UAV-based GPR-SAR system, no issues were found due to the use of the antenna array on board the UAV with respect to the previous prototypes based on larger UWB Vivaldi antennas.
Measurements were processed with the DAS backpropagation technique described in Section II-B. The two metallic targets were placed on the ground, so no SVD filtering was applied in the processing. In [27], different sources of information were considered to obtain the height above the ground, concluding that the use of the height information extracted from the radar measurements provided more accurate reflectivity images compared with considering the laser rangefinder height information. This technique has been further developed to estimate the height of each antenna in the array, as explained in Section II-D. A comparison between the reflectivity images obtained using the height information given by the laser rangefinder and extracted from radar measurements is depicted in Fig. 12(a) and (b), respectively. As can be observed, the latter results exhibit a higher signal-to-clutter ratio, and the shape of the metallic bar is fully reconstructed. This is in agreement with the analysis presented in [27], confirming that the method proposed in Section II-D can be used to estimate the height of each antenna in the array from the radar measurements. The metallic plate shows a reflectivity level comparable to that of some areas of the surrounding terrain which, as in the example of Section III, was a bare wet loamy soil (see Fig. 10). This, together with the fact that the metallic plate is placed practically at the same z level as the soil, contributes to the presence of clutter in some parts of the image, which is especially noticeable around y ≈ 5.5 m.

B. In-Flight Validation With Buried Targets
The capability of the array-based GPR-SAR system on board the UAV for subsurface imaging is presented in this section. This experimental validation campaign was conducted in October 2021 at the Spanish military training and shooting range "El Palancar," located north of Madrid [60]. A picture of the area scanned with the array-based GPR-SAR system on board the UAV is shown in Fig. 13. The scenario consisted of a 12-m along-track × 4.5-m across-track section of a dirt road, where ten targets of different shapes and compositions were buried. These targets, which mainly comprise inert replicas of IEDs and landmines, are depicted in Fig. 14, identified with roman numerals, and described in Table II.
The targets were buried by experts in counter-IED techniques [1] to ensure that they were placed as similar as possible to realistic scenarios. These experts also annotated the placement of each target within the 4.5 m × 12 m validation scenario prior to the realization of the flights with the UAV, but the position of the targets was not disclosed to the research team conducting the flights. In fact, the validation of the detection results followed a blind procedure in which the research team who implemented the prototype had to estimate the position of the targets from the analysis of the GPR-SAR images obtained with the prototype. Then, these results were forwarded to the counter-IED team who buried the targets, comparing the estimated positions with the true ones. Finally, they provided the research team who implemented the prototype with a report indicating the targets correctly detected, the targets missed, and the number of false alarms.
A video of the UAV conducting the scanning of the validation scenario can be watched in the supplemental multimedia file "VideoArrayFlightMadrid." In a similar fashion to the previous in-flight validation example, the flight speed of the UAV was set to 50 cm/s, and the flight height was 1.5 m, resulting in 23 048 measurements to be processed.
The flight path followed by the UAV is plotted in Fig. 15(a). As in the previous examples, along-track sweeps were spaced 20 cm in the across-track direction, requiring 24 along-track sweeps to cover the 4.5 m × 12 m investigation domain. The different sweeps can be seen in Fig. 15(a), where it can be observed that there is a significant resemblance between the actual spacing between the along-track-sweeps and that of the predefined flight path. The overall flight time to complete the scanning was 10 min. The heatmap showing the  Table II. distribution of the measurements on each 5 cm × 5 cm cell in which the scanned area was divided is depicted in Fig. 15(b). When converting this heatmap into a binary plot, Fig. 15(c), it can be observed that the overall sampling of the area under scan is sufficient. It should be noted that there is a 15-cm offset in the y-axis between the array center and the origin of the system of coordinates, which contributes to having a small nonsampled area around y = 0 m.
Once the measurements with the antenna array were collected, the GPR-SAR processing algorithm was applied to obtain the 3-D GPR-SAR images of the subsoil and the targets buried in it. In this scenario, the relative permittivity of the soil was estimated to be around ε r ≈ 4.
The recovered GPR-SAR images are depicted in Fig. 16. In particular, those z cuts where the buried targets were identified were selected. Concerning the qualitative identification of the targets, the antitank mine TM-62 (ii) exhibits the strongest reflectivity due to its size and metallic composition (the echo associated with this target is observed in the four depicted z cuts). The two VS-1.6 antitank plastic landmines (i) and (x)  are also clearly noticed in Fig. 16(c), z = −12 cm. The elongated shape of the 81-mm mortar shell (ix), buried parallel to the air-soil interface, can be identified in Fig. 16(a) and (b). The plastic jug filled with ammonium nitrate (viii) also creates a strong reflection (more than 15 dB above the surrounding clutter) due to its large size and the dielectric contrast between the soil (ε r = 4) and the ammonium nitrate (ε r around 7.1 [7]).
The plastic polyvinyl chloride (PVC) pipe (iii), the antipersonnel (AP) plastic landmine (iv), the wooden trunk-like IED (vi), and the wooden pressure plate (PP) (vii), exhibit lower reflectivity contrast with respect to their surrounding area (equal or less than 10 dB) than the previous targets, making their detection more challenging from a visual inspection of the GPR-SAR images. The reasons why these targets are more difficult to detect are mainly their size and composition (no metal content) and their proximity to the air-soil interface, which partially masks their responses.
Concerning the 120-mm mortar shell (v), this target was buried tilted with respect to the air-soil interface. Thus, despite its size (74 cm long) and metallic composition, the echo of this target is around 10 dB weaker than the echo of the 81-mm mortar shell (ix) buried horizontally.
The range or depth resolution achieved by the radar subsystem ( r = 3.75 cm taking into account the permittivity of this soil) is sufficient to identify some details of the buried targets. Thus, the representation of the vertical cuts of the 3-D GPR-SAR images helps to resolve unclear detections observed on the horizontal cuts. To illustrate this, a vertical cut of the 3-D GPR-SAR image is depicted in Fig. 17. In this vertical cut, the echo associated with the wooden PP (vii) is better observed than in Fig. 16, and also the one corresponding to the plastic jug (viii). In the case of the latter target, it had a motorbike battery next to it, whose upper side is also observed in the vertical cut of Fig. 17.
The cross-range resolution achieved by the GPR-SAR system also helps to distinguish some features of large targets in the horizontal cuts. For example, in the case of the 81-mm mortar shell (ix), the two reflectivity peaks observed in Fig. 16(a) correspond to the projectile body and the tail. Another large target is the antitank mine TM-62 (ii) whose upper side is not flat: the fuze and the center part of the case are approximately 4 cm more elevated than the outer case (a picture is shown in Fig. 14). A detail of the GPR-SAR image centered on this target is depicted in Fig. 18: the fuze is imaged in the z = −8 cm cut, whereas the circular shape imaged in the z = −12 cm cut corresponds to the outer case of the TM-62.
Apart from the qualitative analysis of the GPR-SAR images, they were also processed with a cell-averaging CFAR detector (CA-CFAR) [61]. In brief, this CA-CFAR detector estimates the noise variance for the cell under test by selecting a range of neighboring cells within each cut of the GPR-SAR image. The detector is based on the a priori assumption that the noise in the neighboring cells follows the same distribution as in the cell under test, and that there are no targets in these neighboring cells. The parameters of the CA-CFAR detector considered in this example are the same as for the GPR-SAR images of [27,. When applying the CA-CFAR detector to the images depicted in Fig. 15, all the targets were detected. The number of false alarms of the CA-CFAR detector for this scenario was only one.

C. Comparison of Fully Multistatic and Multiquasi-Monostatic Configurations
As mentioned in Section II-A, the proposed multistatic GPR system benefits from the spatial diversity of the antennas (which contributes to reducing the clutter level) and from the different illumination angles (which helps to obtain complementary information). In order to illustrate how the proposed approach helps to improve the detection capabilities, the GPR-SAR images of the smallest targets in the previous scenario are compared when all Tx-Rx channels are used (fully multistatic configuration, left column of Fig. 19) and when only those Tx-Rx corresponding to a quasi-monostatic arrangement are considered (multi-quasi-monostatic configuration, right column of Fig. 19). As explained in Section II-C, the latter includes the following channels: Tx1-Rx1, Tx1-Rx2, Tx2-Rx2, Tx2-Rx3, Tx3-Rx3, and Tx3-Rx4.
In particular, the targets analyzed are the AP landmine P5 (iv), the wooden PP (vii), and the plastic PVC pipe (iii). The plastic AP landmine is detected with both configurations (Fig. 19(a) and (d)), but with the fully multistatic configuration the reflectivity is higher and the clutter is slightly lower. In the case of the wooden PP, its top face can only be detected with the fully multistatic configuration, as observed when comparing Fig. 19(b) and (e). Finally, in the case of the PVC pipe (Fig. 19(c) and (f)), the fully multistatic configuration provides an image in which the target exhibits a significantly better contrast with the surrounding clutter. Therefore, it can be concluded that when all the channels Fig. 19. GPR-SAR images when all channels are considered (fully multistatic configuration, on the left) and when only the quasi-monostatic channels are considered (multi-quasi-monostatic configuration, on the right) for the targets: (a) and (d) plastic AP landmine P-5, (b) and (e) wooden PP, and (c) and (f) piece of PVC pipe. are considered (fully multistatic configuration), these small shallow targets are better detected as their reflectivity is higher and the surrounding clutter is lower.

V. CONCLUSION
In this contribution, the use of an antenna array to improve the scanning capabilities of a UAV-based GPR-SAR system has been presented, focusing on increasing the scanning throughput without jeopardizing the detection of threats. Extensive tests have been conducted in different scenarios to show the performance of the proposed system for subsurface imaging while keeping most features of the previous prototypes in terms of detection capabilities of buried targets. The switching scheme, flight speed, and across-track spacing have been optimized to achieve the maximum scanning throughput while ensuring that the sampling rate is fulfilled in almost the entire investigation domain.
The scanning throughput achieved with the antenna array is four times greater than the one of previous DLGPR prototypes without the array [23], [27], where the spacing between along-track sweeps was 5 cm to fulfill the sampling rate at 3000 MHz. This means that with the use of the array, a 75% reduction of the time required to scan the investigation domain is achieved, which for scenarios like the one of Section IV-B, allows completing the scan of the full area with a single flight. For this scenario, 24 along-track sweeps are now required to scan the 4.5 m × 12 m area (instead of 91 sweeps as with the previous prototype). Given a flight speed of 50 cm/s, the use of the array results in a scanning throughput of, approximately, 54 m 2 in 10 min. This is a throughput significantly higher than the one achieved by other UAV-based GPR architectures with similar resolution, such as the one based on a side-looking architecture presented in [24] (with a scanning throughput of around 33 m 2 in slightly more than 12 min).
Finally, it should be highlighted that the usage of an antenna array brings further benefits in addition to the increase in the scanning throughput. Combining the radar images of all Tx-Rx channels coherently allows us to reduce the clutter level (thanks to the spatial diversity of the antennas). Furthermore, the proposed multistatic configuration allows illuminating the scanned scene from different angles, which helps to detect challenging threats (in particular, small, nonmetallic, and shallowly buried targets). Therefore, the proposed system paves the way for the usage of UAV-mounted GPR technology for the detection of buried threats in real scenarios.