Tracking Infiltration Front Depth Using Time-lapse Multi-offset Gathers Collected with Array Antenna Ground Penetrating Radar

Summary

Here we present a Ground Penetrating Radar (GPR) system based on a ground-coupled, densely populated antenna array for monitoring the dynamic process of subsurface water infiltration. A time-lapse radar image of the infiltration process allowed estimating the depth of the wetting front during the course of the infiltration process.

Abstract

A Ground Penetrating Radar (GPR) system based on a ground-coupled, densely populated antenna array was used to collect data during an infiltration experiment conducted at a test site near the Tottori Sand Dune, Japan. The antenna array used in this study consists of 10 transmitting antennas (Tx) and 11 receiving antennas (Rx). For this experiment, the system was configured to use all possible Tx-Rx pairings, resulting in a Multi-Offset Gather (MOG) consisting of 110 Tx-Rx combinations. The array was left stationary at a position directly above the infiltration area and data were collected every 1.5 seconds using a time-based trigger. Common-Offset Gather (COG) and Common Mid-Point (CMP) data cubes were reconstructed from the MOG data during post-processing. There have been few studies that used time-lapse CMP data to estimate changes in velocity of propagation. In this study, electromagnetic (EM) wave velocity was estimated heuristically at 1-minute intervals from the reconstructed CMP data through curve fitting, using the hyperbola equation. We then proceeded to calculate the depth of the wetting front. The evolution of the wetting front over time obtain through this method is consistent with the observations from a soil moisture sensor which was placed at a depth below 20 cm. The results obtained in this study demonstrate the ability of such array GPR system to monitor a subsurface dynamic process like water infiltration accurately and quantitatively.

Introduction

Understanding mass and energy transport processes in the vadose zone are important for many applications in agricultural and environmental disciplines. Among those processes, variably saturated water flow is the essential process as many of the other processes, such as physical, geochemical, biological, and even mechanical processes, are usually coupled with water flow. Recent developments in geophysical techniques have allowed one to monitor hydrological processes in the vadose zone non-invasively. Among many geophysical techniques, ground penetrating radar (GPR) is one of the most widely used techniques to monitor and to characterize soil water dynamics because the propagation of electromagnetic (EM) waves emitted and received by GPR antennas is characterized by soil moisture content^1,2,3,4. Among available systems, on-ground surface GPR (referred to as surface GPR in the remainder of the manuscript) is the most common to use in a field. Traditional surface GPR systems with one transmitter and one receiver (bistatic radar systems) are commonly used to scan the subsurface with a constant transmitter-receiver separation (offset). Data sets collected in this configuration are also known as common offset gathers (COG). Radar data are displayed as time series based on total travel time between the transmitter, eventual reflectors and back to the receiver. In order to convert the travel time to depth information, EM wave velocity in the subsurface needs to be estimated. For example, this can be done through the analysis of multi-offset gather (MOG) data sets⁵.

Although there have been a number of studies using GPR to monitor subsurface infiltration processes^6,7,8,9, none of them determined directly the location of the wetting front or the EM wave velocity structure that is changing with time during infiltration. The common approach is to use objects buried at known depths as reference reflectors to determine average EM wave velocity and wetting front depth. Since the wetting front changes dynamically during infiltration, time-lapse MOG must be collected at short time intervals to exploit changes in the EM wave velocity structure without using reference objects. With common bistatic surface GPR antennas, the collection of time-lapse MOG at short intervals between each other is difficult or impossible as it requires to manually move antennas to set up the various offset configurations. Recently, a family of antenna array GPR (referred to as array GPR hereafter) has been used extensively to image the subsurface fast and accurately¹⁰. The basic concept of the array GPR is to provide dense swaths with a minimal effort by electronically switching multiple antennas mounted in a single frame. Array GPR systems have been mainly used to generate 3D subsurface images of wide areas quickly. Some examples of typical applications for these systems are road and bridge inspection¹¹, archeological prospection¹² and UXO and landmine detection^13,14. For such purposes, the array GPR is mainly used to scan the subsurface with a constant antenna separation configuration to collect COG. Although it has been demonstrated that MOG collected with an array GPR could be used for velocity estimation¹⁵, the practical application of this methodology has been limited to just a few cases. By placing the antenna array at a fixed location, time-lapsed MOG can be easily collected. As demonstrated in our recent publication¹⁶, time-lapse radargrams collected with the array GPR system portray rather clearly the reflection signals from the wetting front as it gradually moves downwards during a vertical infiltration experiment carried out on a sand dune. The main objective of the paper was to demonstrate how to use the array GPR to collect time-lapse MOG during the infiltration test and how to analyze such data for tracking the depth of the wetting front.

In this study, we used an antenna array consisting of 10 transmitting (Tx0 - Tx9) and 11 receiving (Rx0 - Rx10) bow-tie monopole antennas. The displacement of the antenna elements within the array is shown in Figure 1 (see Table of Materials). The antenna array is controlled by a Step-Frequency Continuous Waveform (SFCW) radar unit operating over the frequency range from 100 MHz to 3,000 MHz. The array GPR switches through a user-defined sequence of Rx-Tx pairs using radio frequency (RF) multiplexers in the antenna array¹⁰. The maximum number of Tx-Rx combinations for this particular system is 110. For this experiment, we configured the array GPR to use all 110 combinations, programming the scan sequence so that each transmitter, from Tx0 to Tx9, was paired sequentially with all 11 receivers from Rx0 to Rx10. The time necessary to perform a scan through all 110 combinations is less than 1.5 seconds. The offset between transmitter and receiver was calculated based on the distance between the feed points of the antenna elements, where the vertical offset is 85 mm as shown in Figure 1.

Protocol

The main objective of this study was to use an array GPR system to visualize subsurface dynamic processes, such as water infiltration, it is also important to describe how the infiltration test was conducted.

1. Field Infiltration Test

1. Select a site with a bare surface for the infiltration test.
   Note: For this study, the infiltration test was conducted in a large greenhouse facility of the Arid Land Research Center of Tottori University, Japan. See the recent publication for more details¹⁶.
2. Position on the surface six 2.5-m long porous tubes parallel at a distance of 15-cm between each other, so as to ensure uniform infiltration (Figure 2).
3. Connect one end of the tubes to a valve regulated inlet, which is connected to a water reservoir or to a faucet, and the other end to an outlet in order to remove excess water.
4. Place a thin 910 x 1802 mm² wood panel (slightly larger than the antenna array), to cover the porous tubes.
5. Check the effect of the wood panel on the radar signal prior to the infiltration test.
6. Install soil moisture sensors right next to the wood panel to monitor changes in the soil moisture content at various depths.
   Note: The purpose of these sensors is to provide reference data used to cross-check the results obtained from the analysis of the GPR data. In this study, a rod-type soil moisture sensor (see Table of Materials) was installed next to the antenna to monitor the soil dielectric constant at 10, 20, 30, 40, 60, and 100 cm depths.

2. GPR Measurement

1. Place the array GPR antenna on the wood panel.
   Note: In this study, the antenna was left stationary to ensure consistency in the recorded signals.
2. Connect the antenna to the controller using coaxial cables.
3. Connect the controller to a laptop PC with an Ethernet cable.
4. Decide the antenna combinations and their sequence to ensure both common offset gather (COG) and multi-offset gather (MOG) are collected seamlessly.
   Note: In this study, all the possible combinations for the array GPR at our disposal was used (110). It took less than 1.5 s to perform a full scan through all combinations.
5. Start the recording of the GPR data using the data collection software.
   Note: It is important to start the data acquisition before water is injected into the porous tubes.
6. Open the valve to apply water (i.e., begin the infiltration process).
   Note: The water injection rate in the subsoil can be determined from the hydraulic properties (e.g., hydraulic conductivity) of the target soil. The flow rate used in this study was 7,000 cm³/min, which was less than the saturated hydraulic conductivity of the dune sand.
7. After a predetermined amount of water is injected, close the valve.
   Note: Water was injected continuously for 4 h, and a total of 1,680 L of water were used for the test.
8. Keep scanning with the array GPR for an additional period of time, in order to monitor the redistribution process for further analysis.
   Note: The array GPR has been kept collecting data for 4 h after the valve had been closed. Nonetheless, these data have not been analyzed for the study presented here.

3. Data Analysis for Velocity Estimation

1. Reconstruct COG and CMP data by extracting the relative Tx-Rx combinations from the overall data cube. This was achieved through some custom code developed by the authors which rearranged the radar profiles so to group together Tx-Rx pairs at identical offset (COG) and the ones with the same mid-point (CMP).
   Note: Because a SFCW GPR unit was used for this study, the raw data is in the frequency domain; the operation of frequency to time conversion was performed using the software provided by the manufacturer of the GPR system.
2. Estimate the optimum EM wave velocity in the wetting zone at a given elapsed time by fitting the computed two-way travel time, t_c, of the EM wave reflected at the infiltration front to the observed signals. The velocity with the best curve fit was determined heuristically.
       (1)
   where d₀ is the depth to the reflection point, x is the distance between Tx and Rx, v_r is the root-mean-square velocity to the reflection plane and t₀ is the two-way travel time at zero offsets determined from the CMP radargram.
   Note: The most common velocity estimation method using MOG is semblance analysis^5,15,17. In semblance analysis, the root mean square velocity is picked using a velocity spectrum. Semblance analysis was not used for velocity estimation in this study because the low signal-to-noise ratio observed in the array GPR data did not allow us to develop a reliable velocity spectrum. Finding an appropriate filtering method to generate a reliable velocity spectrum from these data is being investigated in other ongoing research. For the scope of this study, the curve fitting method was used to fit Eq (1) to the CMP radargrams.
3. Using Eq. (1), compute the wetting front depth d₀ at a given elapsed time using the estimated EM velocity v_r.
   Note: For this part of the data analysis, the COG data set was used at an offset x equal to 113 mm.

Results

Figure 3 shows a time-lapse panel diagram of 110 Tx-Rx combinations acquired every 1.5 s during the first 60 min of the infiltration experiment, after converting the data from frequency to time domain. In order to enhance reflections at depth, a bandpass filter was applied followed by gain compensation. The panel diagram can be divided into 10 sections, each section corresponding to a particular Tx. The location of the Tx is indicated by a white triangle, and each slice corresponds to a time-lapse signal recorded for an Rx. The vertical axis shows the two-way travel time or TWT, i.e. the time necessary for the EM wave to travel from the transmitter, to the eventual reflector and back to the receiver. The profiles display the signal amplitude in a grayscale color map. A large contrast in color indicates high amplitude in the recorded radar signal. EM wave reflections are produced at the interface between layers at different dielectric constant or by objects with different electrical characteristics than the surrounding medium; a high dielectric contrast will determine a high amplitude reflection. During the water infiltration process, there is a zone referred to as the transition zone where the water content gradually increases from the edge of the wetting zone, that is where water penetrates into the initially dry soil. The EM wave is likely to be reflected not at the very edge but within the transition zone, as observed in studies of water table detection¹⁸. In the remainder of the manuscript, this reflection area is referred to as the wetting front. In Figure 3, a high amplitude signal appears and moves steadily downwards as time advances during the experiment. This reflection is indeed produced by the wetting front as the water gradually penetrates downwards in the subsurface. From this diagram, COG and CMP can be reconstructed as shown in Figure 2 of Iwasaki et al.¹⁶

Velocity analysis was performed on the CMP data obtained every 1 min. For every set of CMP data, the two-way travel time given by Eq. (1) was fitted to the reflection from the wetting front by adjusting t₀ and v_r, assuming a uniform layer in the wet zone. Time zero was corrected by fitting the air wave with the velocity of 0.3 m/ns. Figure 4 shows CMP data at 5 min intervals from elapsed time t_e = 5 min to t_e = 50 min along with the best-fit curves displayed as white lines (solid lines for the reflected wave and dashed line for the air wave). Since the air wave is not a reflected wave but the direct signal between transmitter and receiver, the travel time increases linearly with the offset. All the curves were fitted to the positive peaks (in white color) of the reflected waves. All the curves fitted well to the observed reflected curves shown in the CMP radargrams at all offsets, which means that the estimated values for t₀ and v_r are good. For the infiltration experiment, a dry wood panel was placed between the antenna and the porous tubes. Because the panel has a much lower dielectric constant than that of wet soil, its effects on the EM wave propagation may be not negligible, even though it is thin. A two-layer model was then considered in addition to the aforementioned uniform layer model, assuming a value of 3 for the dielectric constant of the top 5 cm. Also for this second model, the EM wave velocity v_r was estimated by curve fitting at the reflection produced by the wetting front.

In Figure 5, estimated wetting front depths are plotted as a function of t_e for both the single and two-layer models. It can be appreciated that the wetting front moves downwards almost linearly with time for both models, with the exception of a slow-down between t_e = 10 min and t_e = 20 min. Differences between two models are initially not significant, but as time elapses the estimate for the uniform model moves slightly faster compared to the two-layer model. In Figure 5, diamond symbols are used to mark the times when readings from the moisture sensors started to increase and when they later became steady; these are connected with a solid line for every sensor depth. As mentioned above, the reflection of the EM wave does not necessarily occur at the very edge of the wet zone; in other words, considering a certain depth, one cannot expect this reflection to be matching the point in time when the readings from a sensor start to increase. In this sense, the reflection is attributable to a depth level at the infiltration front where a certain water saturation has been reached, compared to the area immediately below. Considering the sensors at 30, 40, and 60 cm in depth, the estimate of the wetting front depth obtained from the GPR data falls well in the range shown by the solid lines on the timeline. The time when the GPR estimated wetting front arrived at 20 cm depth corresponds to the time when the sudden increase in the sensor reading was observed, while the GPR estimate reached the depth of 10 cm much faster than what was produced by the moisture sensor, although the reflection signal from the wetting front is clearly observed after t_e = 5 min (Figure 4). Also, it has to be mentioned that the extrapolation of the GPR estimated curve does not pass through the origin. Although it is not clear what caused this discrepancy at shallower depths, there may be a few possible explanations. It can be attributed to heterogeneity in soil properties, or it may be due to non-uniformity in water application. If that is indeed the case, this would have a greater effect earlier during the infiltration process than at a later stage. Another explanation could be that surface roughness affects the determination of time zero. In addition to the effect of the wood panel and the porous tubes, the effect of surface roughness should be taken into account.

Figure 1: Schema of array ground penetrating radar antenna configuration used in this study. V-shaped structures are bowtie monopole antennas. There are 10 transmitting antennas (Tx) and 11 receiving antennas (Rx) aligned horizontally. Please click here to view a larger version of this figure.

Figure 2: Schemas of the infiltration experiment. (A) a top view and (B) a side view in which the array antenna was placed on top of six 250-cm porous tubes aligned 15 cm apart. Please click here to view a larger version of this figure.

Figure 3: Time-lapse radargram obtained during the first 60 min of the infiltration experiment. Data consist of recorded signals for 110 antenna combinations. One slice corresponds to time-lapse data collected with a single Tx-Rx combination. A different color is used for signal amplitude. Vertical black lines separate data for each Tx. White triangles indicate the locations of Tx. Please click here to view a larger version of this figure.

Figure 4: Radargrams of the CMP at t_e = 5 min to t_e = 55 min at 5 min intervals. The white solid lines represent manually fitted two-way travel time of the reflection from the wetting front, while the white dashed lines represent the travel time of the air wave. Please click here to view a larger version of this figure.

Figure 5: Wetting front depths. Wetting front depths estimated from array GPR in time-lapse multi-offset gather as a function of the elapsed time both for the uniform (triangles) and the two-layer (square) models. Black lines with diamonds at both ends show the time between the initial increase in readings and when these reached a steady level for each sensor depth (i.e., the duration of transition zone). Please click here to view a larger version of this figure.

Discussion

In this study, array ground penetrating radar (GPR) was used to track the depth of the wetting front during an infiltration experiment conducted at an experimental field near the Tottori Sand Dune, Japan. The array GPR system used in this study consists of 10 transmitting antennas (Tx) and 11 receiving antennas (Rx). The system can be configured to use up to 110 different Tx-Rx combinations. During the infiltration experiment, all 110 combinations were scanned continuously at 1.5 seconds intervals, leaving the array stationary at the location where water was applied through some porous tubes placed on the surface. Common-offset gather (COG) and common mid-point data (CMP) were reconstructed from the time-lapse data cube. It would be practically impossible to collect CMP data at the same rate with conventional bistatic GPR systems. It is very important to leave the antenna in stationary recording during the experiment in order to obtain reproducible and meaningful time-lapse data.

Although array GPR data have been used to estimate EM wave velocities¹⁵, there are only a few studies that analyzed time-lapse array GPR data to estimate EM wave velocities for transient processes like water infiltration. In this study, the electromagnetic (EM) wave velocity structure was estimated from the time-lapse CMP data. Instead of performing semblance analysis, the hyperbola curve for the two-way travel time was fitted heuristically to the reflected signals in the CMP radargrams in order to estimate the average EM wave velocity in the wetting zone because of low signal-to-noise (S/N) ratio in the data. When the S/N ratio is low, semblance analysis cannot be used to generate a reliable velocity spectrum. An appropriate filtering method would need to be developed for using the semblance analysis method. A rod-type soil moisture sensor was installed next to the GPR antenna to measure the changes in soil moisture content during the infiltration experiment; the sensors were displaced at the depths of 10, 20, 30, 40, and 60 cm and worked independently.

Using the estimated EM wave velocity, the depth of the wetting front was computed at 1-minute intervals of the infiltration process. The evolution in time of the estimated wetting front agrees well with the observations from the soil moisture sensors at depths below 20 cm. At shallower depths, the GPR estimate of the wetting front depth shows a discrepancy with the readings from the soil moisture sensors.

Overall, this study demonstrates that the array GPR system is capable of tracking the evolution of the wetting front depth during water infiltration into the soil, by collecting time-lapse common mid-point (CMP) data. As this type of data was not easily collected from conventional surface GPR before, the data obtained in this study are the very first ones that actually show how the wetting front evolved over time in the subsoil. Future work will explore the possibility of using data inversion to estimate the hydraulic parameters of the soil from the data obtained during this experiment.

Materials

List of materials used in this article

Name	Company	Catalog Number	Comments
GeoScope Radar Unit	3D Radar AS
DXG1820 antenna	3D Radar AS
PR2/6 Profile Probe	Delta-T