Journal of Applied Geophysics 47 Ž2001. 241–249 www.elsevier.comrlocaterjappgeo
Estimation of groundwater level by GPR in an area with multiple ambiguous reflections Yuichi Nakashima, Hui Zhou, Motoyuki Sato ) Center for Northeast Asian Studies, Tohoku UniÕersity, Sendai, 980-8576, Japan Accepted 8 May 2001
Abstract Ground-penetrating radar ŽGPR. was applied to estimate the groundwater level in Ulaanbaatar City, Mongolia. Multiple reflectors occur at different levels. Any one of these was a possible water table reflector. In order to identify the groundwater reflection amongst a set of ambiguous reflectors and to increase the signal-to-noise ratio, the common-midpoint ŽCMP. method was used. The operation frequency of the radar system was 100 MHz. The start position of a survey line was 10 m away from a building that influenced GPR data. However, a reflection from the building can be suppressed significantly by two-dimensional f–k filtering. By noise reduction, velocity analysis, normal moveout correction, stacking and migration using estimated velocities, the subsurface images are obtained. The signal-to-noise ratio of the stacked CMP profile is much higher than that of the common-offset profile. The usual detectable range of GPR is several meters, however, we can measure more than 10 m in depth because of the improvement of the signal-to-noise ratio by CMP and due to the very dry soil. Also, we can estimate the vertical dielectric constant distribution from the interval velocities obtained from the CMP. Then we can locate a groundwater reflection within a multitude of ambiguous reflectors, and find that the depth of the groundwater level is approximately 8.0 m. The experiment results show that GPR data collected by the CMP method not only improves the quality of the radar profile but also enables us to estimate the subsurface vertical profile of the dielectric constant. q 2001 Elsevier Science B.V. All rights reserved. Keywords: Ground-penetrating radar; Common midpoint method; Velocity analysis; Vertical dielectric constant estimation; Ground water survey
1. Introduction Ground-penetrating radar ŽGPR. has been widely applied to investigate subsurface structures or buried objects. GPR records electromagnetic waves reflected from boundaries of subsurface material having different dielectric constants. In general, GPR data are collected with the common-offset method, in )
Corresponding author. Fax: q81-22-217-6075. E-mail address:
[email protected] ŽM. Sato..
which the distance between a transmitting and receiving antenna is fixed ŽDavis and Annan, 1989; Arcone 1996; Lorenzo et al., 1998.. In low attenuation media, this method can produce high-resolution images of shallow subsurface. However, the dielectric constant cannot be calculated using the common-offset data. The common-midpoint ŽCMP. method was developed as a fundamental data processing method in the 1950s, and made public by Mayne Ž1962.. Recent studies have shown that the CMP data have some
0926-9851r01r$ - see front matter q 2001 Elsevier Science B.V. All rights reserved. PII: S 0 9 2 6 - 9 8 5 1 Ž 0 1 . 0 0 0 6 8 - 4
242
Y. Nakashima et al.r Journal of Applied Geophysics 47 (2001) 241–249
advantages over common-offset data ŽFisher et al., 1992; Grasmueck, 1994.. The first advantage is the improvement of the signal-to-noise ratio. The CMP method records the reflections from the same reflection point several times as the distance between a transmitting and a receiving antenna increases. Therefore, stacking of the signals increases the energy of the reflected signals, because the reflected waves are coherent, and the energy of noise is decreased due to its randomness or residual normal moveout. The coherent noise, like reflections from objects on the ground surface, can be reduced by two-dimensional f–k filtering before stacking, since it has very different velocity. The second advantage of the CMP is the possibility of normal-moveout ŽNMO. velocity estimation. The interval velocities can be obtained from NMO velocities using the Dix Ž1955. formula, and then the dielectric constant, which is directly related to the water content, can be estimated by using the interval velocities. The GPR profile acquired by the common-offset method shows the reflection from geological anomalies, however, it gives only the travel time from the anomaly, and we cannot obtain the information about the material of the anomaly. When there is a medium, like a groundwater aquifer, it will be identified by its NMO velocity, because of its high dielectric contrast compared to dry soil. Robert et al. Ž1996. show the effect of velocity analysis on the stacked radar image and particularly demonstrated that the NMO velocity can be used to determine the subsurface water content. The final advantage is the improvement in the accuracy of spatial positioning by applying imaging techniques to GPR data. Using velocity information obtained by velocity analysis, time or depth migration can be applied to the stacked profile. After migration, the quality of the image of a complex structure can be improved significantly. In order to verify that CMP improves GPR, we apply the processing to GPR data acquired in the south part of the city of Ulaanbaatar, Mongolia. The objective of the survey was to detect the groundwater level. We selected this site because we can use a large area of flat ground, which is especially suitable for CMP data acquisition. In addition, the soil is dry and easy to apply GPR. The dielectric constant of the groundwater layer is higher than that of other layers.
Therefore, it is possible to identify a groundwater reflection within a multitude of ambiguous reflectors by inverting the velocity information to obtain the dielectric constant. Arcone et al. Ž1998. show that the dielectric constant distribution of the subpermafrost groundwater can be obtained by inverting velocity information. In this paper, we discuss the acquisition, processing and interpretation of CMP data at the experiment site.
2. Experiment Ulaanbaatar, the capital city of Mongolia and the center of industry and commerce of the country, has a high concentration of population in the nation. Recently, the demand for construction of buildings in the city is increasing. Although it is a very dry region, the seasonal change of local groundwater level is a serious problem in the city. When buildings are constructed, the groundwater sometimes flows into the basement or foundation excavation. Therefore, it is important to know the depth of the water table prior to excavation. Japan International Cooperation Agency ŽJICA. has conducted a project, AThe Study on Water Supply System in Ulaanbaatar and SurroundingsB, during 1993–1995. The final report ŽJICA, 1995. summarizes that the ground water in the Ulaanbaatar area is supplied along Tuul River, and has an aquifer of 5–10-m depth. The rainfall occurs only in July and August and seasonal ground water level change is large. The experiment site is located at the south part of the city of Ulaanbaatar and between the city and the Tuul River. It is above the aquifer of the Tuul River. This site was grassy and flat as shown in Fig. 1Ža., and the soil was very dry on the day when the GPR survey was carried out in August 1999. GPR data were acquired along several survey lines. In this paper, we show one example from the GPR data. The length of the survey line is 100 m and a pumping well facility, which can influence GPR data, is located 10 m away from the beginning of the survey line as shown in Fig. 1Žb.. We applied three GPR survey methods, namely common-offset, common-source and CMP surveys. A common-offset profile can be collected very quickly because of its easy measurement style. Therefore, it is suitable for
Y. Nakashima et al.r Journal of Applied Geophysics 47 (2001) 241–249
243
Fig. 1. The experiment site Ža. and the survey line Žb.. The start position of the survey line is 10 m away from a pumping well facility.
determining the broad- or large-scale structure. In most GPR surveys, this common-offset method is used. A common-source gather is collected by fixing
the transmitting antenna position, and moving the receiving antenna continuously along a survey line and collect data with a spatial interval between two
244
Y. Nakashima et al.r Journal of Applied Geophysics 47 (2001) 241–249
receiving stations. These gathers can be transformed to CMP gathers according to the antenna positions. Common-source collection method is always used to acquire a large quantity of data. Sometimes, only several CMP gathers are needed to get velocities by velocity analysis. These CMP gathers are always obtained by putting the transmitting and receiving antenna symmetrically about a CMP position, stepping them out at constant increments from a central point and recording radar data with different offsets. The geometry of this observation method can be easily controlled. The velocity estimated by this CMP gather sometimes gives a better result than that from the CMP gathers transformed from commonsource gathers, because the former can easily collect more data at one center position than the latter. For the reasons mentioned above, we applied CMP measurements to estimate the dielectric constant, and used common-source measurements to determine the detailed underground structure. Common-offset gathers were collected every 0.10 m from 0 to 100 m along the survey line. The offset was fixed at 1.0 m. CMP gathers were collected at two center positions Ž14 and 15 m.. The smallest offset was 1.0 m, the largest offset was 30 m, and the offset increment was 0.40 m. Common-source
gathers were collected from 5.5 to 18 m along the same survey line. The nearest offset was 0.50 m, the farthest offset was 30 m, the interval of transmitting antenna was 0.50 m and the interval of receiving antenna was 0.10 m. A pulse radar system with 100 MHz antennas ˚ GeoScience, RAMAC GPR. was used in ŽMALA this experiment. We selected this frequency in order to survey relatively deep regions. Since this radar system has unshielded antennas, we have to remove reflections from objects on the ground surface. This radar system has an initial time uncertainty, and the initial time is corrected in the data processing. Since the initial time uncertainty can cause a considerable relative error in the final dielectric constant, we must be attentive to correct it.
3. Data processing A GPR survey has many similarities to a seismic survey. However, GPR suffers from strong attenuation and dispersion from the subsurface material, so not all seismic signal-processing schemes can be used in GPR. However, when the material is very
Fig. 2. The raw CMP gather whose center position is 15 m Ža. and after noise reduction Žb..
Y. Nakashima et al.r Journal of Applied Geophysics 47 (2001) 241–249
dry, the standard sequence of the seismic data processing can be followed in the GPR data processing ŽFisher et al., 1992.. The procedure applied here consists of noise reduction, velocity analysis, NMO correction, stacking and migration. 3.1. Noise reduction Fig. 2Ža. shows the raw CMP gather whose center position is 15 m. Note that all profiles are displayed after AGC processing in this paper. From Fig. 2Ža., we notice that there are two kinds of noise. One is high frequency noise which would be caused by the radar system. This noise can be removed by using a bandpass filter. A bandpass filter from 30 to 120 MHz was used in this processing. An arrow in Fig. 2Ža. indicates another kind of noise. The apparent velocity of this event is very high. In practice, the apparent velocities of subsurface-reflected waves in a CMP gather are much slower, therefore, these events would be reflections from objects on the ground surface ŽSun and Young, 1995.. We knew that there was a water facility made
245
Table 1 Estimated NMO velocities at 14 m Ža. and 15 m Žb. Two-way time Žns.
NMO velocity Žmrns.
(a) 63 80 120 190
0.148 0.148 0.142 0.126
(b) 55 80 125 180 225
0.146 0.148 0.139 0.124 0.120
of bricks located 10 m away from the beginning of the survey line. The event shown by the arrow is the reflection from the facility. This kind of noise can be easily suppressed by two-dimensional f–k filtering in a CMP gather. In f–k domain, this kind of noise is distributed around k s 0 because the apparent velocity of the noise is very high. An f–k filter is designed to remove small wave number components.
Fig. 3. Velocity spectra of the CMP profiles at 14 m Ža. and 15 m Žb..
246
Y. Nakashima et al.r Journal of Applied Geophysics 47 (2001) 241–249
Fig. 2Žb. shows the profile processed from Fig. 2Ža. after the noise reduction. The high frequency noise has been suppressed and the energy reflected from the facility on the ground has been greatly reduced. 3.2. Velocity analysis, NMO correction, stacking and migration The traveltime of a reflection in a CMP gather can be expressed by an equation that is a function of
an antenna offset and NMO velocities ŽYilmaz, 1987.. Therefore, a CMP data set can provide a velocity profile in the subsurface media by velocity analysis. We apply the energy-normalized cross-correlation sum velocity spectrum ŽYilmaz, 1987. to CMP gathers after noise reduction. Fig. 3Ža. and Žb. shows the velocity spectra of the CMP gather whose center positions are 14 and 15 m, respectively. In each figure, the normal-moveout velocity decreases with increasing depth. In general, the velocity of electromagnetic waves decreases rapidly with depth at shallow depths ŽDavis and
Fig. 4. The stacked profile Ža., the common-offset profile Žb., the migrated profile of the stacked profile Žc. and the migrated profile of the common-offset profile Žd..
Y. Nakashima et al.r Journal of Applied Geophysics 47 (2001) 241–249
Annan, 1989.. This is primarily a result of increasing water saturation with depth. From the velocity spectra, several pairs of two-way time and NMO velocities are picked and summarized in Table 1. NMO corrections calculated from the NMO velocity remove the offset dependence of the reflection traveltimes. The CMP gather after NMO correction becomes the data set with zero-offset traces. The traces are stacked to a single CMP trace after the NMO correction. Stacking increases the energy of the reflected signals because the reflected waves are spatially coherent, and decreases the energy of the noise due to its randomness or residual NMO. If the NMO correction is carried out with a correct NMO velocity, the stacked profiles have an improved signal-to-noise ratio. Theoretically, the signal-to-noise ratio of the stacked profile is improved 6n times over the random noise when the number of traces of a CMP gather is n. Finally, migration can be applied to the stacked profile using the NMO velocities. After migration, the accuracy of spatial positioning is improved. In this study, we used Kirchhoff migration ŽDillon, 1988.. Fig. 4Ža. shows the stacked profile from 5.5 to 18 m along the survey line. The raw data were collected with the common-source method. Common-source gathers are transformed to CMP gathers. The smallest offset is variable from 0.1 to 0.9 m, and the interval offset is 1.0 m. Then the signal processing methods outlined above are applied. A CMP gather is stacked using 14 traces. Fig. 4Žb. is a common offset GPR profile whose offset is fixed at 1.0 m and whose position is the same as in Fig. 4Ža.. The trace interval is 0.10 m in each figure. Comparing Fig. 4Ža. and Žb., we can see some reflections from the subsurface layers in Fig. 4Ža., for instance around 50, 100 and 150 ns, but we cannot see them clearly in Fig. 4Žb.. We can find that the former has higher signal-to-noise ratio than the latter. Theoretically, the signal-to-noise ratio of the stacked profile is improved 3.7 times over random noise because the number of traces of a CMP gather is 14. Fig. 4Žc. and Žd. shows Kirchhoff migration profiles using NMO velocities from 15 m. It is only the vertical axis that has been depth-converted using the average of NMO velocities from 15 m. The data is time-migrated not depth-migrated. Although the de-
247
tectable range of GPR is usually several meters in wet soil, we can detect the reflection from a layer more than 10 m deep. It is because the soil is very dry at this survey site and we get a high signal-tonoise ratio profile by CMP stacking.
4. Dielectric constant estimation The NMO velocity is approximately same as the RMS velocity assuming the medium from the surface to the boundary be homogeneous, when the subsurface consists of a multiple horizontal layers. Therefore, in order to estimate the dielectric constant of each layer, the RMS velocities have to be corrected to interval velocities. An interval velocity Õ i, n can be calculated using the Dix Ž1955. formula:
n i ,n s
)
2 2 t 0, n nrms , n y t 0, ny1 nrms , ny1
t 0, n y t 0, ny1
Ž 1.
where Õ i, n is the interval velocity of layer n; Õrms, n is the RMS velocity down to the bottom of layer n; Õrms, ny1 is the RMS velocity down to the bottom of layer n y 1; t 0, n is the two-way traveltime down to the bottom of layer n; and t ny 1 is the two-way traveltime down to the bottom of layer n y 1. The thickness of layer n can be calculated by: dn s
n i , n Ž t 0, n y t 0, ny1 . 2
.
Ž 2.
On the other hand, the relationship between interval velocity Õ i and dielectric constant ´ r can be expressed as:
nis
c
(´
Ž 3. r
where c is the electromagnetic velocity in free space. The dielectric constant of the groundwater layer is higher than those of other layers because the water content of the former is higher than that of the latter. Therefore, GPR is sensitive to the groundwater table. After the velocity analysis and dielectric constant estimation using the CMP gather, we can estimate the depth of the groundwater level.
Y. Nakashima et al.r Journal of Applied Geophysics 47 (2001) 241–249
248
Table 2 Estimated dielectric constant distributions at 14 m Ža. and 15 m Žb. Two-way time Žns.
Interval velocity Žmrns.
Depth Žm.
Dielectric constant
(a) 63 80 120 190
0.148 0.148 0.129 0.092
0–4.7 y5.9 y8.5 y12.0
4.1 4.1 5.4 10.6
(b) 55 80 125 180 225
0.146 0.152 0.121 0.080 0.102
0–4.0 y5.9 y8.7 y10.9 y17.7
4.2 4.0 6.0 14.0 8.6
Table 2Ža. and Žb. shows the vertical dielectric constant distributions at 14 and 15 m that are the center positions of the two CMP gathers. These are calculated with Eqs. Ž1–3. using the results of the velocity analysis shown in Table 1. By comparing Table 2Ža. and Žb., the dielectric constant distributions are quite consistent with each other. We find that the dielectric constant of the near surface layer is very small because the soil is very dry. However, the dielectric constant of a layer, which is deeper
than 8.0 m, is higher than those of any other layers. This value of the layer is over 10, and corresponds to a high water content layer. Therefore, we think that this layer is the groundwater table and the depth of the upper boundary is about 8.0 m. Fig. 5 shows the migrated profile that is same as Fig. 4Žc. and the estimated vertical dielectric constant at 15 m Ža vertical line.. Raw data of the former was collected with the common-source method and the latter was estimated with the CMP method. The vertical distribution of the dielectric constant corresponds to a multitude of reflected waves from the subsurface structures. Comparing the two results, we can identify the groundwater layer within a multitude of ambiguous reflecting layers and find the depth of the upper boundary shown with a white line is about 8.0 m. JICA curried out many measurements at this area from August 1993 to March 1995 in order to study the water supply system. JICA measured the geological samples near the experiment site. According to the JICA Ž1995. report, the surface is covered by a clay layer with a thickness of 0.50 m. The second layer between 0.50 and 6.0 m in depth is light gray sand and gravel with clay granule. The third layer between 6.0 and 8.0 m in depth is blue sand and gravel with clay. The fourth layer between 8.0 and
Fig. 5. The migrated profile of the stacked profile with the dielectric constant distribution.
Y. Nakashima et al.r Journal of Applied Geophysics 47 (2001) 241–249
30 m in depth is sand and rounded to subrounded gravel which may be a good aquifer. The estimated depth and dielectric constant of each layer from the GPR data correspond to the hydrogeologic data in the JICA report quite well. An important aspect is that the estimated depth of the upper boundary of the groundwater layer Ž8.0 m. is verified by the JICA report. Our results show that GPR has high potential to estimate dielectric constant distribution.
249
Acknowledgements During this field test, colleagues from the city of Ulaanbaatar and Mongolian Technical University including Prof. Batsukh, Mr. Byambaa and Nemer Buyankhishig helped us in better understanding the problems. We greatly acknowledge their support. The original planning of the radar survey was due to Prof. Gerel. We would like to thank the reviewers for carefully checking the manuscript and for their valuable suggestions.
5. Conclusions References GPR data collected by the CMP method improve the subsurface images over usual common-offset gathers. The CMP method improves the signal-tonoise ratio not only by stacking, but also by two-dimensional f–k filtering before stacking. In the GPR data, there is regular and irregular noise. The regular noise can be greatly suppressed by f–k filtering, and the random noise can be decreased by stacking. An important aspect is the ability to estimate velocities by velocity analysis using techniques from reflection seismology. The dielectric constants can be calculated from the estimated velocities. When there is a medium, like a groundwater layer, whose dielectric constant is quite different from the non-saturated media, it is possible to identify the medium within a multitude of ambiguous media and find the depth. We applied CMP method to a flat grassy plain in Ulaanbaatar. There was a pumping facility on the ground that influenced GPR data. The facility created noise like flat reflections. It is not so easy to remove this kind of noise from conventional GPR data. However, we can eliminate that from CMP gathers using f–k filtering. The stacked CMP profile confirmed the improvement of signal-to-noise ratio over common-offset profile, and we can measure more than 10 m in depth. Furthermore, we have successfully identified the groundwater layer at a multitude of ambiguous layers by the dielectric constant distribution with CMP method, and calculated the depth of the groundwater level. It is verified by JICA Ž1995..
Arcone, S.A., 1996. High resolution of glacial ice stratigraphy: a ground-penetrating radar study of Pegasus Runway McMurdo station, Antarctica. Geophysics 61, 1653–1663. Arcone, S.A., Lawson, D.E., Delaney, A.J., Strasser, J.C., Strasser, J.P., 1998. Ground-Penetrating Radar Reflection Profiling of Groundwater and Bedrock in an Area of Discontinuous Permafrost. Geophysics 63, 1573–1584. Davis, J.L., Annan, A.P., 1989. Ground-penetrating radar for high-resolution mapping of soil and rock stratigraphy. Geophys. Prospect. 37, 531–551. Dillon, P.B., 1988. Vertical seismic profile migration using the Kirchhoff integral. Geophysics 53, 786–799. Dix, C.H., 1955. Seismic velocities from surface measurements. Geophysics 20, 68–86. Fisher, E., McMechan, G.A., Annan, A.P., 1992. Acquisition and processing of wide-aperture ground-penetrating radar data. Geophys. Prospect. 57, 495–504. Grasmueck, M., 1994. Application of seismic processing technique to discontinuity mapping with GPR in crystalline rock of the Gotthard Massif, Switzerland. Proc. 5th Int. Conf. on GPR. University of Waterloo, pp. 1135–1149. JICA ŽJapan International Cooperation Agency., 1995. The Study on Water Supply System in Ulaanbaatar and Surroundings, Final Report. Lorenzo, H., Cuellar, V., Hernandez, M.C., 1998. GPR detection of archeological galleries and tunnels. Proc. 7th Int. Conf. on GPR. University of Kansas, pp. 661–666. Mayne, W.H., 1962. Common-reflection-point horizontal data stacking techniques. Geophysics 27, 927–938. Robert, J.G., Davis, P.L., Jung, M.L., Toksoz, ¨ M.N., 1996. Velocity variations and water content estimated from multi-offset, ground-penetrating radar. Geophysics 61, 683–695. Sun, J., Young, R.A., 1995. Recognizing surface scattering in ground-penetrating radar data. Geophysics 60, 1378–1385. Yilmaz, O., 1987. Seismic data processing. Soc. Exp. Geophys. 154–234.