Published online 25 February 2008
Published in Vadose Zone J 7:215-226 (2008)
DOI: 10.2136/vzj2006.0137
© 2008 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
SPECIAL SECTION: GROUND PENETRATING RADAR IN HYDROGEOPHYSICS
Calibration of a Vadose Zone Model Using Water Injection Monitored by GPR and Electrical Resistance Tomography
Rita Deianaa,
Giorgio Cassianib,*,
Alberto Villaa,
Andrea Bagliania and
Vittorio Brunoa
a Dipartimento di Scienze Geologiche e Geotecnologie, Univ. of Milano-Bicocca, Milan, Italy
b Dipartimento di Geoscienze, Univ. of Padua, Padova, Italy
* Corresponding author (giorgio.cassiani{at}unipd.it).
All rights reserved. No part of this periodical may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, recording, or any information storage and retrieval system, without permission in writing from the publisher.
Received 17 September 2006.
 |
ABSTRACT
|
|---|
The noninvasive characterization of the vadose zone has been the objective of intense research aiming at collecting data for the calibration of flow and transport models. This interest is motivated by the difficulty of accessing the vadose zone with direct methods without causing disturbance to the natural conditions. On the other hand, characterizing the vadose zone requires detailed knowledge of space and time variations of water distribution. Geophysical methods capable of imaging the moisture content, such as ground penetrating radar (GPR) and electrical resistance tomography (ERT), provide data with high resolution in space and time: natural infiltration and artificial injection (tracer) tests can be imaged, by repeating measurements over time (time-lapse mode). We conducted a forced-injection experiment at a test site in Gorgonzola, east of Milan (Italy). The site is characterized by Quaternary sand and gravel sediments that house an extensive unconfined aquifer. We used ERT and GPR in two-dimensional cross-hole configuration and time-lapse mode during a period of several days preceding and following a tracer experiment with the injection of about 20 m3 of fresh water in a purposely excavated trench. The calibration of a three-dimensional infiltration model on the basis of the geophysical data provides useful estimates of vertical saturated hydraulic conductivity; however, we observed a discrepancy in mass balance between calibrated simulations and measurements. This fact cannot be explained using a simple isotropic model for hydraulic properties; therefore, the method requires further investigation before being fully established.
Abbreviations: ERT, electrical resistance tomography GPR, ground penetrating radar MOG, multioffset gather ZOP, zero-offset profile
 |
INTRODUCTION
|
|---|
The hydrologic characterization of the vadose zone is technically challenging, particularly when the investigation extends deeper than a couple of meters below ground. All direct measurements are invasive, based on drilling and sample collection, and consequently difficult to perform without causing major disturbance to the natural in situ conditions, particularly volumetric water content. In addition, invasive measurements cannot be repeated over time on the same soil sample, thus giving only a single snapshot of volumetric water content conditions. Also, the support volume of invasive measurements is small (at most on the cubic decimeter scale) and localized, while soil water content can have large spatial variations and predictive models usually can represent reality only at a much larger scale. All these facts have contributed to the increasing use of noninvasive—geophysical—methods, particularly GPR (Annan, 2005) and ERT (Binley and Kemna, 2005), often in cross-borehole configuration, to investigate the vadose zone (Hubbard et al., 1997; Slater et al., 1997; Binley et al., 2001, 2002b; Alumbaugh et al., 2002; French et al., 2002; Cassiani et al., 2004; Chang et al., 2004; Schmalholz et al., 2004, among others; for reviews see Huisman et al., 2003; Cassiani et al., 2006a).
A direct link between measurable geophysical quantities (e.g., dielectric constant and electrical resistivity) and hydrologic model parameters such as hydraulic conductivity and unsaturated flow parameters is generally not available and in the best case it is highly site specific. More reliable links can be established between geophysical quantities and hydrological variables such as water content and solute concentration, generally in the form of empirical or semiempirical relationships—consider the classical relationships proposed by Archie (1942) for electrical conductivity and by Topp et al. (1980) and Roth et al. (1990) for the dielectric constant. Using such relationships, it is possible, albeit not always straightforward, to obtain quantitative estimates of hydrologic data to be used for the calibration of hydrologic models, thus identifying the parameters of interest.
This successful use of geophysical data for hydrologic investigations requires (Fig. 1
)
- that data are collected in time-lapse mode, i.e., repeated with time on the same geometry, to capture the interaction between static geological signals and dynamic hydrological processes;
- that the collected geophysical data have a clear, identifiable, and quantitative petrophysical relationship to environmental variables of interest; and
- that the resolution and sensitivity of geophysical methods in space and time is fully understood and appropriately taken into account in the calibration of the hydrologic model.

View larger version (23K):
[in this window]
[in a new window]
|
FIG. 1. A conceptual scheme for the use of geophysical data to calibrate unsaturated water flow models. The key issues related to the process are highlighted on the right.
|
|
For this use of geophysical data, experience shows that cross-hole acquisition provides the most informative results, particularly thanks to the superior resolution of this acquisition geometry. This is, of course, particularly true in the deep vadose zone, where surface-to-surface measurements fail to provide suitable information. A number of studies have been performed in recent years to try to exploit the information content of cross-hole time-lapse ERT and GPR data during natural and forced infiltration to characterize the vadose zone (e.g., Binley et al., 2002a; Binley and Beven, 2003; Cassiani et al., 2004; Cassiani and Binley, 2005; Winship et al., 2006).
Intuitively, the calibration of unsaturated flow models requires that the system is adequately stressed so that different states of the variables are explored—i.e., volumetric water content varies significantly with time. This condition is not guaranteed to be met under natural conditions in the deep vadose zone (e.g., Binley et al., 2002b; Cassiani and Binley, 2005). Consequently, it is faster and more informative to conduct water injection experiments that can be monitored in a time span of hours, days, or weeks depending on the dynamics of the system.
There were two objectives of this study. The first was to conduct a water injection experiment in the unsaturated zone of a Quaternary sediment formation in the Po River plain, monitored via cross-hole GPR and ERT. While the methodology adopted is similar to that presented by Binley et al. (2002a), these Quaternary sediments were expected to be orders of magnitude more permeable than the sandstone described by those researchers, with strong consequences on the timing and accuracy of time-lapse monitoring. The second objective was to calibrate a three-dimensional, unsaturated flow model on the data derived from the cross-hole geophysical monitoring, with the aim of estimating unsaturated hydraulic parameters at the field scale, assessing the data value and limitations.
 |
Field Site Description
|
|---|
The Gorgonzola experimental site was specifically designed for noninvasive cross-hole monitoring of vadose zone dynamics and consequently for the assessment of the vulnerability of the underlying aquifer. The site is located few kilometers east of Milan, in the Po River valley, northern Italy. This is a heavily industrialized region, where important underlying aquifers exist and are potentially threatened by contamination. The aquifer and the unsaturated zone are composed of Quaternary sediments with a fairly coarse sand–gravel grain size distribution. The water table oscillates yearly between 13.5 and 20 m below ground as a consequence of seasonal irrigation. The site was initially equipped (Fig. 2
) in October 2004 with three boreholes (internal diameter 0.0762 m [3 inches]) that extend to a depth of 20 m, i.e., cover the entire vadose zone under all water table conditions. These boreholes (A, B, and C in Fig. 2) are permanently equipped with a set of 24 stainless steel borehole electrodes for ERT imaging, spaced 0.8 m vertically, between 2- and 20-m depth. Contact with the surrounding formation, made of loose to poorly cemented sand and gravel, was ensured by backfilling the boreholes with a mixture of drilling cuttings and a small percentage of concrete powder between 2 and 15 m below ground. A gravel pack around the slotted section was placed between 15 and 20 m to ensure that the changing elevation of the water table could be measured at all times. Analysis of a soil core extracted in July 2005 from a fourth borehole (D in Fig. 2) has provided direct knowledge of the site stratigraphy (Fig. 3
). Several features can be noted that have a direct impact on the site's hydrology: (i) the presence of a 3-m-thick surface layer composed of rubble and agricultural soil, partly overconsolidated by construction operations in a nearby parking lot—this layer is present only locally and is of no interest in terms of characterization of the vadose zone in the area; (ii) a prevalence of Quaternary sand and gravel sediments of fluvial and glacial origin, highly permeable; and (iii) the existence of a 2-m-thick layer of cemented gravel and sand between 12 and 14 m below ground.

View larger version (30K):
[in this window]
[in a new window]
|
FIG. 2. Scheme of the Gorgonzola experimental site. All four boreholes were drilled to 20-m depth. Three boreholes (A, B, and C) were equipped with electrodes for electrical resistance tomography (ERT) measurements, while Borehole D was used for core retrieval only. The trench was dug to 2-m depth for the purpose of water injection during the tracer experiments. Cross-hole ERT and ground penetrating radar time-lapse data were collected between Boreholes A and B.
|
|

View larger version (12K):
[in this window]
[in a new window]
|
FIG. 3. Description of the core retrieved from Borehole D and results of the corresponding laboratory analyses for saturated hydraulic conductivity and grain size distribution (as 10, 50, and 90% mass fractions); b.g.l. = below ground level.
|
|
During drilling of the cored Borehole D, a Lefranc permeability test (Hvorslev, 1951) was performed at a depth of 6 m. The resulting field-saturated hydraulic conductivity value was 1.5 m/h. Given the test characteristics, this value corresponds to an estimate of horizontal hydraulic conductivity.
The grain size distribution of the core material was measured in the laboratory using dry and wet sieving (ASTM International, 2004, D421 and D422); the fraction passing the 75-µm sieve (ASTM International, 2004, no. 200) was characterized using densimeters (ASTM International, 2004, 151H and 152H) (using Stokes' law). Permeability was measured in the laboratory using a constant-head water permeameter (ASTM International, 2004, D2434). Laboratory analyses of material from the retrieved core are summarized in Fig. 3: both grain size distribution and saturated hydraulic conductivity above 6-m depth are characterized by a large fraction of fine to very fine sediments that is not observed in the deeper sediments. As discussed below, we believe that this fine material comes from washing out of the upper 2 m, composed mostly of allochtonous soil, during mud-drilling operations: neither grain size distribution nor hydraulic conductivity are representative of actual in situ conditions between 2- and 6-m depth. In July 2005, the site was equipped also with a 2.6-m-wide, 2-m-deep trench (Fig. 2) located halfway between Boreholes A and B. The trench was designed to cut the shallowest rubble and agricultural soil layer and allow injection of water directly into the natural sand and gravel formation below.
 |
Cross-Hole Geophysical Imaging
|
|---|
The main data set collected at the Gorgonzola site is composed of repeated (time-lapse) cross-hole measurements performed using both borehole GPR and ERT. Both measurements were collected periodically starting in January 2005 over the three two-dimensional planes defined by Boreholes A, B, and C (Fig. 2). Only Boreholes A and B were used during the water injection experiment.
Cross-Hole Electrical Resistance Tomography
Electrical resistance tomography data were collected using an IRIS Syscal Pro system (IRIS Instruments, Orléans, France) with 10 physical channels and 48 multiplexed channels. The acquisition was based on a dipole–dipole measurement scheme built as a mixture of AB–MN and AM–BN schemes (A and B refer to current electrodes, M and N refer to potential electrodes, see Binley and Kemna, 2005), to benefit both from the good resolution of the AB–MN scheme (where current and potential dipoles are small, with both electrodes of the same type in the same borehole) and the good signal/noise ratio of the AM–BN scheme. The latter scheme (advocated, e.g., by Bing and Greenhalgh, 2000) proved practically redundant in this case, however, as the resistivity images produced by using only the AB–MN data are identical to the images using the full data set. During the water injection experiment, in addition to the borehole electrodes, eight electrodes were placed at the ground surface along the line connecting Boreholes A and B, with 0.9-m spacing: these electrodes were necessary to cover the upper part of the A–B section. To accommodate these surface electrodes in the 48-channel capability of the instrument, eight borehole electrodes (the six deepest and two electrodes with poor contact) were removed from the measurement scheme. A total of 3168 measurements were taken for each ERT image, including direct and reciprocal configurations. The latter are made by interchanging the electrode pairs used for voltage measurement and the ones used for current injection. In theory, the resistance for direct and reciprocal configuration shall be identical. Therefore, a comparison between direct and reciprocal resistance provides an estimate of data error (Daily et al., 2005). The total time required for each acquisition, optimizing the use of the 10 physical channels of the IRIS Syscal Pro, was about 20 min.
Inversion of ERT time-lapse data requires special care, as changes in resistivity cannot, in general, be successfully imaged by carrying out independent data inversions and subtracting pixel-by-pixel values. In most cases, changes with time are small relative to the natural spatial variability within the region of interest. In addition, when using two-dimensional acquisition and inversion while monitoring three-dimensional processes, such as a tracer injection, off-plane changes may induce a considerable level of data error, which changes with time. Successful approaches to overcome these problems have been proposed by Daily et al. (1992) (ratio inversion—e.g., Binley et al., 2002a; Cassiani et al., 2006b) and by LaBrecque and Yang (2000) (difference inversion—e.g., Kemna et al., 2002).
Inversion of the Gorgonzola ERT data was performed using the CRTomo 2D code (Kemna, 2000), which adopts difference inversion for time-lapse data. The vertical plane defined by Boreholes A and B was discretized into a grid of 56 cells in the y direction and 35 cells in the x direction, accommodating also outer padding regions (allowing for the flow of current away from the interborehole region) of relatively limited extent, thanks to third-type boundary conditions defined at the lower and lateral grid boundaries. For the absolute images, such as the preinjection background, the reciprocal analysis was set to exclude data that exceed a 10% error level, and a corresponding error level was set as the target for the Occam inversion procedure (LaBrecque et al., 1996) adopted in CRTomo. For difference inversion, the error analysis cannot be conducted on the basis of the reciprocals, and tighter error levels are generally adopted (e.g., Kemna et al., 2002). In our case, error levels in the range of 2 to 5% were investigated, and a 2% error level was finally chosen for all difference inversions, based on a necessarily rather subjective judgment of the quality of the resulting images of resistivity change. Conversion from resistivity images into volumetric water content was performed by using Archie's law parameters derived in the laboratory from core data (Fig. 4
).
A number of important limitations must be kept in mind when interpreting cross-hole ERT images of flow and transport processes in the vadose zone. First, regularization is necessary to constrain the number of possible solutions to the nonunique ERT inverse problem. Choice of the regularization method is not defined by the data, but is instead based on a priori information and, generally, on experience. The obtained images depend strongly on the selected regularization, leading to potential misinterpretation of the results (Kemna et al., 2002). Second, for cross-hole electrode arrangements, data sampling and spatial resolution are markedly enhanced in the vertical compared with the horizontal direction, and consequently different regions in the inverted image have substantially different resolution, with important consequences, especially if a quantitative characterization is sought (Day-Lewis et al., 2005). The main resulting limitation is the difficulty in estimating the horizontal extent of anomalies, which are smeared laterally by regularization, and consequently it is basically impossible to assess horizontal and vertical anisotropy by looking at ERT data alone (Binley et al., 2002a).
Cross-Hole Ground Penetrating Radar
Cross-hole GPR can be acquired in different configurations. For multiple-offset gather (MOG) surveys, the receiver is moved to different locations in one borehole while the transmitter remains fixed. The transmitter is then moved and the process repeated. Following collection of all data in this mode and determination of the travel time for each wave path-line, it is possible to derive a tomogram of velocity within the plane of the borehole pair. In contrast, a zero-offset profile (ZOP) configuration is obtained by keeping transmitter and receiver at equal depth. In this manner, only measurements at successive depths are collected. By systematically lowering or raising the pair of antennas in the two boreholes, it is possible to build a one-dimensional profile of average interborehole travel time throughout the entire borehole length. While MOG has the greatest information content, and can represent the actual (albeit two-dimensional) distribution of dielectric properties via tomographic inversion, acquisition is slow and cumbersome. Zero-offset profile is fast to acquire and easy to interpret, but only provides one-dimensional information.
The boreholes drilled at the Gorgonzola site are suitable for cross-hole GPR acquisition, in spite of the presence of electrodes and wires deployed for the permanent ERT installation in the same boreholes. The only obvious limitation arises when the vertical offset between transmitter and receiver antennas exceeds an angle of 45°. This problem prevented us from using vertical radar profiles (e.g., Cassiani et al., 2004), and is also in partial contrast with the need for good multiangle coverage in MOG data. At present we are unable to attribute the problem to the existence of cables in the boreholes, as similar problems were noticed in boreholes with no cables or electrodes (Peterson, 2001; Alumbaugh et al., 2002). A cause of this limitation in cross-hole GPR may be related to first-arrival energy traveling between antenna tips rather than antenna centers (Irving and Knight, 2005).
At the Gorgonzola site, both MOG and ZOP surveys were acquired before, during, and after water injection. A PulseEkko 100 system (Sensors & Software, Mississauga, ON, Canada) was used with 100-MHz borehole antennas. In ZOP surveys, antennas were lowered with 0.25-m vertical spacing. For MOGs, a coarser 0.5-m spacing was selected to minimize the acquisition time. The acquisition parameters (number of stacks and time sampling) were selected to minimize the acquisition time for each shot and still preserve the data quality. In spite of all such efforts, the MOG acquisition time still amounted to about 1 h for each survey.
The GPR data analysis was limited to travel-time picking and inversion for velocity distribution to obtain images of dielectric constant converted to volumetric water content by means of the Topp et al. (1980) relationship. Evidence that this relationship is applicable at the Gorgonzola site comes from the comparison of porosity laboratory measurements on core material with in situ measurements below the water table using GPR ZOP profiles (see below).
The zero-offset profiles were analyzed under the assumption of straight-ray horizontal propagation. While this approach is known to contribute possible errors due to critically refracted waves, the proposed corrections (e.g., Rucker and Ferré, 2004) rely on the hypothesis that layers with sharp velocity changes can be identified. We cannot safely assume that such layering with sharp contrasts is present at the Gorgonzola site. The straight-ray assumption was also adopted for the tomographic inversion of MOG data, which was performed using the code Migratom (Jackson and Tweeton, 1994). For both ZOP and MOG data, time-lapse analysis was based purely on differences between the individual inverted images.
Monitoring of Natural Conditions
Starting in January 2005, monitoring of in situ volumetric water content conditions was performed on three borehole pairs (A–B, A–C, and B–C) using cross-hole ERT and GPR in both ZOP and MOG modes. The results of GPR ZOP monitoring from January to June 2005, converted to volumetric water content estimations using the Topp et al. (1980) relationship, are shown in Fig. 5
for the three pairs of boreholes. The data coming from the three independent data sets are fully consistent, confirming both data quality and the layered structure of the site. The cemented sand and gravel layer at 12- to 14-m depth is clearly visible, having a constantly higher water content than the surrounding formations. This fact is probably due to the higher residual water content of the cement, and possibly also to a lower permeability.

View larger version (19K):
[in this window]
[in a new window]
|
FIG. 5. Moisture content profiles estimated at the Gorgonzola test site from January to April 2005, derived from zero-offset profile radar and the constitutive relationship by Topp et al. (1980). The moisture content curves clearly show the falling water table from January to April. The residual moisture content values ( r) obtained in the laboratory on samples from the retrieved core (to the right) are shown on the profile relevant to Boreholes B to C (b.g.l. = below ground level).
|
|
Very little variation in volumetric water content is noticeable with time, if an exception is made of the changing water table, which is imaged very accurately (in the centimeter range, compared with direct water table measurements). We can therefore conclude that the site undergoes very small changes in saturation conditions under natural infiltration, which is probably strongly limited by the allochtonous soil layer on top. The quasi-steady-state condition observed in Fig. 5 may, in principle, be due to either (i) a dynamic equilibrium between infiltration and drainage through the system, or (ii) drainage to equilibrium, where the observed volumetric water content conditions correspond to water retained against gravity. The latter case is supported by direct measurements of pressure–saturation relationships on the core material, shown also in Fig. 5: the direct measurements confirm that the observed stationary water content profile is probably the drained profile. Note that under quasi-steady-state conditions, and even worse under drained conditions, it is impossible to identify uniquely the unsaturated flow parameters of the vadose zone (Cassiani and Binley, 2005); to give rise to transient, more informative conditions, a water injection test was planned and performed at the site.
 |
The Water Injection Test
|
|---|
Shortly after the trench installation in July 2005, a first preliminary water injection experiment was performed: 3.5 m3 of tap water was injected during a period of about 2 h by pouring water directly into the trench. The tap water electrical conductivity (450 µS/cm) was the same as the conductivity of the local groundwater (that is the source of water for the municipality). The evolution of the injected water plume in space and time was monitored using ERT and GPR measurements for a total of 93 h, from 5 July (the date of water injection) to 11 July 2005. Even though the monitoring technique proved successful, it became clear that the rapid dynamics of the system and the relatively small water volume injected did not allow an appropriate imaging of the water slug evolution using the frequency of acquisition adopted in that test. Also, the slug was totally dispersed 2 d after the injection and the vertical location of the center of mass was difficult to trace after 1 d, after reaching a depth of about 6 m. In spite of these limitations, a first attempt at calibrating a three-dimensional infiltration model based on Femwater (Lin et al., 1997) led to an estimate of saturated hydraulic conductivity in the range of 0.21 to 0.42 m/h, i.e., about one-third of the horizontal hydraulic conductivity measured by in situ testing. The calibration was based on the vertical movement of the center of mass of the injected slug (Binley et al., 2002a), as imaged by GPR ZOPs and ERT.
A second test with more frequent time sampling was designed and performed in January 2006 using the same injection trench. Also, a much larger water volume was injected, amounting to 23 m3 in 11 h, to ensure that the signal would be traceable well into depth. The timing of water injection and the monitoring scheme are shown in Fig. 6
. A further advantage of this second test was that the water table was close to the deepest point in its yearly oscillation (19 m below ground level), thus allowing characterization of the entire unsaturated zone thickness.

View larger version (29K):
[in this window]
[in a new window]
|
FIG. 6. Measured water injection flow rate at the Gorgonzola site during the test in January 2006, and corresponding monitoring plan using cross-hole ground penetrating radar (zero-offset profile [ZOP] and multiple-offset gather [MOG]) and cross-hole electrical resistance tomography (ERT). The photo shows the trench and the flow control device during the injection experiment.
|
|
Before the beginning of water injection, ERT, GPR MOG, and GPR ZOP background measurements were made. Figure 7
shows the background images for ERT and GPR MOG, both converted into estimates of volumetric water content using Archie's law and the Topp relationship, respectively. A strong borehole effect is noticeable in the ERT image, partly distorting the overall image pattern and probably affecting also the inverted resistivity values. We believe this effect was caused by the conductive backfilling of the boreholes. For this reason, here we mainly focus on the results obtained using GPR, even though the difference inversion adopted for ERT time-lapse data seems to successfully remove such borehole effects (see below). The MOG background image confirms the presence of a layer having lower electrical conductivity and higher dielectric constant at a depth of 12 to 16 m. This is a layer wetter than the surrounding material, and corresponds to the layer of cemented gravel and sand detected by the drilling core in Borehole D.

View larger version (53K):
[in this window]
[in a new window]
|
FIG. 7. Moisture content background images in January 2006 (before water injection), as derived from electrical resistance tomography (ERT) and multiple-offset gather (MOG) ground penetrating radar cross-hole imaging using a laboratory-calibrated Archie's law and the relationship of Topp et al. (1980).
|
|
Given the need to acquire frequent data, especially during water injection and in the first days afterward, ZOP surveys were preferred over the lengthier MOG acquisitions (see Fig. 6). At selected times, however, both modes of cross-hole GPR were acquired, to be compared with each other and with the results of ERT imaging. In Fig. 8
and 9
, the corresponding inversion results, converted to volumetric water content change, are shown at 3 and 45 h after the start of injection. The data at 3 h are fully consistent with each other both in terms of the extent of the anomalous injected slug and of the value of volumetric water content increase above background. The results are also fairly similar at 45 h, when the signal is obviously weaker.

View larger version (32K):
[in this window]
[in a new window]
|
FIG. 8. Moisture content changes with respect to background (Fig. 7) as imaged by zero-offset profile (ZOP) and multiple-offset gather (MOG) ground penetrating radar and electrical resistance tomography (ERT) at 3 h after injection start (b.g.l. = below ground level).
|
|

View larger version (31K):
[in this window]
[in a new window]
|
FIG. 9. Moisture content changes with respect to background (Fig. 7) as imaged by zero-offset profile (ZOP) and multiple-offset gather (MOG) ground penetrating radar and electrical resistance tomography (ERT) at 45 h after injection start (b.g.l. = below ground level).
|
|
The time evolution of the injected water is traced with evidence by both ERT and ZOP surveys, which compare very well in terms of both the center-of-mass depth and the values of volumetric water content change (Fig. 10
).
The complete sequence of volumetric water content change vertical profiles, as derived from ZOP GPR, is shown in Fig. 11
. Three phases can be identified in the water slug motion: (i) during the 11 h of water injection, the total mass "seen" by GPR steadily increased, while some water pooled on top of the cemented layer at a depth of about 12 m (Fig. 11A); (ii) in the early drainage phase (14–31 h), water was released progressively from the shallower layers and penetrated into the cemented layer (Fig. 11B); and (iii) in the late drainage phase, water was released from the cemented layer and reached the water table at 19-m depth (Fig. 11C).
The nature of the observed water slug movement shows a clear need to take into account vertical heterogeneity, in terms of both hydraulic conductivity and residual volumetric water content. This conclusion is consistent with the results of natural-condition monitoring discussed above.
 |
Calibration of the Unsaturated Flow Model
|
|---|
The geophysical data acquired during water injection, converted into volumetric water content, were used to calibrate a variably saturated flow model. The approach is very similar to the one successfully applied by Binley et al. (2002a) and is based on tracking the downward movement of the center of the injected water volume. As in Binley et al. (2002a), only the saturated hydraulic conductivity (Ks) was chosen as the key fitting parameter to be calibrated, because Ks has a dominant effect on gravitational sinking of the water mass.
The adopted variably saturated flow model was Femwater (Lin et al., 1997), which had been already successfully used in similar studies (Binley et al., 2002a; Winship et al., 2006). Femwater solves Richards' equation using a three-dimensional finite element formulation. The vadose zone underlying the Gorgonzola site was discretized using a vertical spacing of 0.05 m, with about 100,000 prism elements and 108,000 nodes. The three-dimensional mesh represents only one-quarter of the site, thanks to symmetry around the trench. Injection was simulated as a second-type boundary condition at the trench bottom, honoring the measured flow rates (Fig. 6). No-flow boundary conditions were applied on the model sides, while a water table condition was set at the model bottom. Initial conditions in the model were set to correspond to gravity-drained conditions as observed by long-term monitoring of natural conditions (Fig. 5). Since ZOP data were collected every 0.25 m in depth, the model was subdivided vertically in layers of the same thickness, each having a different residual volumetric water content (from ZOP data) and potentially a different saturated hydraulic conductivity. Porosity was set for the entire model to be equal to 0.28, consistent with the observed volumetric water content below the water table (see Fig. 5) and with direct laboratory measurements on core material. A van Genuchten parameterization was adopted for unsaturated flow properties, choosing parameters representative of a coarse sand and gravel: n = 2.2 and
= 2 m–1, adopted for the entire model. This van Genuchten parameterization adequately reproduces the pressure–saturation curves measured in the laboratory on core material (see Fig. 4). Only saturated hydraulic conductivity was modified by trial and error to match field data. Figure 12
shows the results of calibration in terms of the depth of the center of mass of the injected water. The ZOP, MOG, and ERT data are all consistent in terms of center-of-mass sinking.

View larger version (18K):
[in this window]
[in a new window]
|
FIG. 12. Results of three-dimensional flow model calibration based on electrical resistance tomography (ERT), multiple-offset gather (MOG), and zero-offset profile (ZOP) estimates of the center of mass depth as a function of time (Ks = saturated hydraulic conductivity; b.g.l. = below ground level).
|
|
When the laboratory-determined Ks data were used in the model, the center of mass sank to a depth of only 4 m: the low-conductivity region in the upper 6 m in the model, consistent with lab data (Fig. 3), trapped the water slug and released it only during a very long time. This demonstrations that there is a severe inconsistency between hydraulic conductivity measured in the lab for this upper layer and field evidence—probably the core material was contaminated with the fine material of the shallowest soil washed down during drilling.
On the other hand, if the entire model is given a homogeneous hydraulic conductivity equal to the value estimated in the preliminary injection test (0.42 m/h), the sinking is far too fast, even though the field data are honored in the first few hours, i.e., in the range of time and depth explored in the preliminary test.
Based on the assumption that Ks laboratory data are biased to low values for the upper 6 m, we modified Ks in the model to have a value equal to the mean of laboratory data below 6 m, i.e., 0.21 m/h (which accidentally happens to be half of the value estimated from the preliminary injection experiment). The corresponding curve of simulated center-of-mass sinking fits the field data (Fig. 12). A homogeneous model having the same mean Ks fits the data as well (Fig. 12). This fact is evidence that the calibration procedure is not capable of distinguishing details in the Ks profile, but is basically sensitive only to its mean value. Other simulated center-of-mass sinking curves vary with continuity as a function of the mean Ks adopted (not shown), as already observed by Binley et al. (2002a).
The individual volumetric water content profiles derived from ZOPs were compared with their analogs extracted from the simulation results, from the model that best fits the data. For each depth in the ZOP profile, the volume that was actually sampled by the GPR energy corresponds to the Fresnel volume (
erven
and Soares, 1992; Galadedara et al., 2003; Winship et al., 2006), so that three-dimensional volumetric water content distribution in the model is averaged across this volume for each ZOP depth and time. The Fresnel volume is approximated as an ellipsoid, with the foci corresponding to the transmitter and antenna locations. The width of the Fresnel zone is a function of the signal dominant wavelength, i.e., of the central frequency and mean GPR velocity. For the relevant equations, see
erven
and Soares (1992). The corresponding results are shown together with the ZOP field data in Fig. 13
. It is quite obvious that mass balance is not honored in the comparison between simulated and measured ZOP data: the Fresnel volume traced through the simulation model contains too much (at least twice as much) water compared with the field data. Note that this mass balance problem does not affect the goodness of fit of the center-of-mass motion.
 |
Discussion and Conclusions
|
|---|
The water injection experiment conducted at the Gorgonzola test site proves that infiltration in permeable Quaternary sediments can be monitored accurately by both borehole ERT and GPR in time-lapse mode, provided that frequent measurements are made to trace the rapid infiltration process. Note that while the adopted methodology is similar to that proposed by Binley et al. (2002a), the hydraulic conductivity of the Quaternary sediments considered is at least two orders of magnitude larger than observed in UK Permo-Triassic sandstones. As a consequence, the infiltration process took place in the time scale of hours rather than many days. Particularly useful in this respect are GPR ZOP surveys that can be conducted in a few minutes. Of course, ZOPs can only provide one-dimensional integrated profiles of volumetric water content estimates. Use of the cross-hole data to calibrate unsaturated flow models proved to be an effective tool to derive unsaturated hydraulic parameters under in situ conditions. It is extremely difficult to derive similar estimates in the unsaturated zone from direct measurements.
The presented three-dimensional model calibration is based on fitting the depth of the injected water center of mass as a function of time. The procedure proved to be very sensitive to small differences in (isotropic) saturated hydraulic conductivity as adopted in the calibrated model. A few limitations of the method have emerged, however.
First, the study shows that the center-of-mass calibration is indeed sensitive only to the mean saturated hydraulic conductivity, while moderate spatial (vertical in this case) variations of hydraulic conductivity are not detectable by this technique.
Second, an evident mass balance problem emerged from this study. This in itself is not new to this type of noninvasive monitoring of hydrologic processes. Similar mass balance misfits were observed, mostly using ERT (Binley et al., 2002a; Singha and Gorelick, 2005), because of resolution limitations at the center of the image away from the boreholes (Day-Lewis and Lane, 2004; Day-Lewis et al., 2005). In some cases, cross-hole GPR has also suffered from mass balance problems (Looms et al., 2008), while in other cases GPR has proven to adequately capture the injected water mass (Winship et al., 2006). A lack of mass balance is clearly disturbing. In the case of this study, we believe that the mass balance problem may be due to horizontal to vertical anisotropy of the medium. In all simulations presented, the hydraulic conductivity was assumed isotropic. A similar approach was also used by Binley et al. (2002a) and Winship et al. (2006), since horizontal spreading in the ERT or MOG images cannot, in general, be measured reliably because of the smearing characteristics of the inversion procedures, albeit resolution patterns are very different for GPR MOG and ERT (Day-Lewis et al., 2005). In our study, however, we have some evidence that anisotropy may be detectable because MOG and ERT images (Fig. 8 and 9, for example) are so similar that resolution issues probably are not prevalent, since MOG and ERT should smear the images in different ways. The horizontal spread (for instance in Fig. 8) is substantial and larger than in the simulation results (not shown), indicating a certain anisotropy in the hydraulic conductivity field that, in fact, brings water away from the center of the bulb. If anisotropy was indeed present, this might explain why the Fresnel volume averaging applied to isotropic simulation results shows too much water mass—in reality, part of this mass has escaped laterally. Future work is needed to calibrate an anisotropic flow model. Note that in this case the estimation of vertical hydraulic conductivity would probably change with respect to the results presented here.
 |
ACKNOWLEDGMENTS
|
|---|
We wish to acknowledge financial support from the Consorzio Industriale di Gorgonzola/Pessano con Bornago, that made this study possible. We also thank Nicoletta Fusi for help with laboratory analysis, and Andreas Kemna for support in using the software CRTomo.
 |
REFERENCES
|
|---|
- Annan, A.P. 2005. GPR methods for hydrogeological studies. p. 185–213. In Y. Rubin and S.S. Hubbard (ed.) Hydrogeophysics. Water Sci. Technol. Library Ser. 50. Springer, New York.
- Alumbaugh, D., P.Y. Chang, L. Paprocki, I.R. Brainard, R.J. Glass, and C.A. Rautman. 2002. Estimating moisture contents in the vadose zone using cross-borehole ground penetrating radar: A study of accuracy and repeatability. Water Resour. Res. 38(12):1309, doi:10.1029/2001WR000754.[CrossRef]
- Archie, G.E. 1942. The electrical resistivity log as an aid in determining some reservoir characteristics. Trans. AIME 146:54–67.
- ASTM International. 2004. Annual book of ASTM standards. Vol. 4.08. ASTM Int., West Conshohocken, PA.
- Bing, Z., and S.A. Greenhalgh. 2000. Cross-hole resistivity tomography using different electrode configurations. Geophys. Prospecting 48:887–912.[CrossRef]
- Binley, A., and K. Beven. 2003. Vadose zone flow model uncertainty as conditioned on geophysical data. Ground Water 41:119–127.[CrossRef][Web of Science][Medline]
- Binley, A.M., G. Cassiani, R. Middleton, and P. Winship. 2002a. Vadose zone flow model parameterisation using cross-borehole radar and resistivity imaging. J. Hydrol. 267:147–159.[CrossRef]
- Binley, A.M., and A. Kemna. 2005. DC resistivity and induced polarization methods. p. 129–156. In Y. Rubin and S.S. Hubbard (ed.) Hydrogeophysics. Water Sci. Technol. Library Ser. 50. Springer, New York.
- Binley, A.M., P. Winship, R. Middleton, M. Pokar, and J. West. 2001. High resolution characterization of vadose zone dynamics using cross-borehole radar. Water Resour. Res. 37:2639–2652.[CrossRef]
- Binley, A.M., P. Winship, L.J. West, M. Pokar, and R. Middleton. 2002b. Seasonal variation of moisture content in unsaturated sandstone inferred from borehole radar and resistivity profiles. J. Hydrol. 267:160–172.[CrossRef]
- Cassiani, G., and A.M. Binley. 2005. Modeling unsaturated flow in a layered formation under quasi-steady state conditions using geophysical data constraints. Adv. Water Resour. 28:467–477.[Medline]
- Cassiani, G., A.M. Binley, and T.P.A. Ferré. 2006a. Unsaturated zone processes. p. 75–116. In H. Vereecken et al. (ed.) Applied hydrogeophysics. NATO Sci. Ser. 71. Springer-Verlag, Berlin.
- Cassiani, G., V. Bruno, A. Villa, N. Fusi, and A.M. Binley. 2006b. A saline tracer test monitored via time-lapse surface electrical resistivity tomography. J. Appl. Geophys. 59:244–259.[CrossRef]
- Cassiani, G., C. Strobbia, and L. Gallotti. 2004. Vertical radar profiles for the characterization of deep vadose zones. Vadose Zone J. 3:1093–1115.[Abstract/Free Full Text]
erven
, S., and J.E. Soares. 1992. Fresnel volume ray tracing. Geophysics 57:902–915.[CrossRef][Web of Science]- Chang, P., D. Alumbaugh, J. Brainard, and L. Hall. 2004. The application of ground penetrating radar attenuation tomography in a vadose zone infiltration experiment. J. Contam. Hydrol. 71:67–87.[CrossRef][Medline]
- Daily, W., A. Ramirez, A. Binley, and D. LaBrecque. 2005. Electrical resistance tomography: Theory and practice. p. 525–550. In D.K. Butler (ed.) Near-surface geophysics. SEG Invest. Geophys. 13. Soc. of Explor. Geophys., Tulsa, OK.
- Daily, W.D., A.L. Ramirez, D.J. LaBrecque, and J. Nitao. 1992. Electrical resistivity tomography of vadose water movement. Water Resour. Res. 28:1429–1442.[CrossRef]
- Day-Lewis, F., K. Singha, and A. Binley. 2005. On the limitations of applying petrophysical models to tomograms: A comparison of correlation loss for cross-hole electrical-resistivity and radar tomography. J. Geophys. Res. 110(B8):B08206, doi:10.1029/2004JB003569.[CrossRef]
- Day-Lewis, F.D., and J.W. Lane, Jr. 2004. Assessing the resolution-dependent utility of tomograms for geostatistics. Geophys. Res. Lett. 31:L07503, doi:10.1029/2004GL019617.[CrossRef]
- French, H.K., C. Hardbattle, A. Binley, P. Winship, and L. Jakobsen. 2002. Monitoring snowmelt induced unsaturated flow and transport using electrical resistivity tomography. J. Hydrol. 267:273–284.[CrossRef]
- Galadedara, L.W., G.W. Parkin, J.D. Redman, and A.L. Endres. 2003. Assessment of soil moisture content measured by borehole GPR and TDR under transient irrigation and drainage. J. Environ. Eng. Geophys. 8:77–86.[Medline]
- Hubbard, S.S., J.E. Peterson, Jr., E.L. Majer, P.T. Zawislanski, K.H. Williams, J. Roberts, and F. Wobber. 1997. Estimation of permeable pathways and water content using tomographic radar data. Leading Edge 16:1623–1628.[Abstract]
- Huisman, J.A., S.S. Hubbard, J.D. Redman, and A.P. Annan. 2003. Measuring soil water content with ground penetrating radar: A review. Vadose Zone J. 2:477–491.
- Hvorslev, M.J. 1951. Time lag and soil permeability in ground water levels and pressures. Bull. 36. U.S. Army Eng. Waterways Exp. Stn., Vicksburg, MS.
- Irving, J.D., and R.J. Knight. 2005. Effects of antennas on velocity estimates obtained from crosshole GPR data. Geophysics 70:K39–K42.[CrossRef]
- Jackson, M.J., and D.R. Tweeton. 1994. MIGRATOM: Geophysical tomography using wavefront migration and fuzzy constraints. Rep. of Invest. 9497. U.S. Bureau of Mines, Minneapolis, MN.
- Kemna, A. 2000. Tomographic inversion of complex resistivity: Theory and application. Ph.D. diss. Univ. of Bochum, Germany.
- Kemna, A., J. Vanderborght, B. Kulessa, and H. Vereecken. 2002. Imaging and characterisation of subsurface solute transport using electrical resistivity tomography (ERT) and equivalent transport models. J. Hydrol. 267:125–146.[CrossRef]
- LaBrecque, D.J., M. Miletto, W.D. Daily, A.L. Ramirez, and E. Owen. 1996. The effects of noise on Occam's inversion of resistivity tomography data. Geophysics 61:538–548.[CrossRef][Web of Science]
- LaBrecque, D.J., and X. Yang. 2000. Difference inversion of ERT data: A fast inversion method for 3-D in-situ monitoring. J. Environ. Eng. Geophy. 6(2):83–89.
- Lin, H.J., D.R. Richards, C.A. Talbot, G.T. Yeh, J. Cheng, and H. Cheng. 1997. FEMWATER: A three-dimensional finite element computer model for simulating density-dependent flow and transport in variably saturated media. Tech. Rep. CHL-97-12. U.S. Army Eng. Waterways Exp. Stn., Vicksburg, MS.
- Looms, M.C., K.H. Jensen, A. Binley, and L. Nielsen. 2008. Monitoring unsaturated flow and transport using cross-borehole geophysical methods. Vadose Zone J. 7:227–237 (this issue).[Abstract/Free Full Text]
- Peterson, J.E., Jr. 2001. Pre-inversion corrections and analysis of radar tomographic data. J. Environ. Eng. Geophys. 6:1–18.
- Roth, K., R. Schulin, H. Fluhler, and W. Hattinger. 1990. Calibration of time domain reflectometry for water content measurements using a composite dielectric approach. Water Resour. Res. 26:2267–2273.[CrossRef]
- Rucker, D.F., and T.P.A. Ferré. 2004. Correcting water content measurement errors associated with critically refracted first arrivals on zero offset profiling borehole ground penetrating radar profiles. Vadose Zone J. 3:278–287.[Abstract/Free Full Text]
- Schmalholz, J., H. Stoffregen, A. Kemna, and U. Yaramanci. 2004. Imaging of water content distributions inside a lysimeter using GPR tomography. Vadose Zone J. 3:1106–1115.[Abstract/Free Full Text]
- Singha, K., and S.M. Gorelick. 2005. Saline tracer visualized with three-dimensional electrical resistivity tomography: Field-scale spatial moment analysis. Water Resour. Res. 41:W05023, doi:10.1029/2004WR003460.[CrossRef]
- Slater, L., M.D. Zaidman, A.M. Binley, and L.J. West. 1997. Electrical imaging of saline tracer migration for the investigation of unsaturated zone transport mechanisms. Hydrol. Earth System Sci. 1:291–302.
- Topp, G.C., J.L. Davis, and A.P. Annan. 1980. Electromagnetic determination of soil water content: Measurements in coaxial transmission lines. Water Resour. Res. 16:574–582.[CrossRef]
- Winship, P., A. Binley, and D. Gomez. 2006. Flow and transport in the unsaturated Sherwood Sandstone: Characterization using cross-borehole geophysical methods. p. 219–231. In R.D. Barker and J.H. Tellam (ed.) Fluid flow and solute movement in sandstones: The onshore UK Permo-Triassic red bed sequence. Spec. Publ. 263. Geol. Soc., London.
This article has been cited by other articles:

|
 |

|
 |
 
S. Lambot, E. Slob, J. Rhebergen, O. Lopera, K. Z. Jadoon, and H. Vereecken
Remote Estimation of the Hydraulic Properties of a Sand Using Full-Waveform Integrated Hydrogeophysical Inversion of Time-Lapse, Off-Ground GPR Data
Vadose Zone J.,
August 11, 2009;
8(3):
743 - 754.
[Abstract]
[Full Text]
[PDF]
|
 |
|

|
 |

|
 |
 
S. Lambot, A. Binley, E. Slob, and S. Hubbard
Ground Penetrating Radar in Hydrogeophysics
Vadose Zone J.,
February 25, 2008;
7(1):
137 - 139.
[Full Text]
[PDF]
|
 |
|