1. Introduction and Objectives
Geophysical applications within the context of urban planning are increasingly recurrent in the literature due to the success of geophysical methods in providing rapid subsurface information, such as identifying buried utilities, predicting their depth and geometry. Knowing this information prior to any activity that aims to drill or excavate underground, such as engineering jobs that aim to the repair and/or install gas, water, sewage, electricity and telephone networks is of great importance. Accidental contact between drilling/excavation tools and buried urban structures can generate unplanned costs on the worksite, require more time for the initial goal to be achieved, or lead potentially to serious accidents, putting people's lives at risk.
One of the geophysical methods that receives great emphasis in this context is GPR (Ground Penetrating Radar), which is able to map underground structures using the principle of transmission and reflection of electromagnetic waves of high frequency (10 MHz - 2.6 GHz)   . In addition to identifying buried targets, GPR can provide a good estimate of target depth, which may be critical to avoid accidental rupture of buried urban artifacts, which can cause serious problems.
The accuracy of depth estimates using GPR depends on precise knowledge of the velocity of electromagnetic wave in the medium. Once the velocity is determined, the double time is converted to depth, and consequently, the depth of the object is determined. This velocity information usually comes from the adjustment of the hyperbola equation in radargram, where the identified diffraction hyperbolas can be associated with buried targets. Depth estimates can also be obtained by means of direct information from previously known targets, together with assumption that the velocity is constant for the entire GPR profile. This is the usual procedure in studies using GPR as evident in the literature, as can be verified in     among others.
The above approach may be appropriate in situations where there are few variations in water content and/or lithology of the host soil of the targets of interest, but may be inappropriate in others, where such variations are more pronounced and may generate erroneous estimates of depth. In this context, the concept of interval velocity developed by Dix  can be explored and applied as an alternative to the constant velocity technique. The Dix principle was initially developed for seismology  , where there are different geological layers and, consequently, different seismic wave velocities within each layer. This principle allows the conversion of the transit time of seismic waves in depth. The same idea is used in the present work; however, these layers are defined by the vertical distances between the tops of the targets that can be identified in radargrams, and not by geological layers, as used in the seismic method. Determining the layers in this manner is results in having more velocity information and transit times of the identified targets, as well as the strategy to smooth the effects of the variation of the dielectric constant in depth, and as a way, to improve as estimates depths of utilities buried underground.
In addition to improving the depth estimates of buried underground utilities, the interval velocity technique holds great potential for understanding wave velocity in the soil matrix, leading to more precise information about lithology, soil compaction and concentration of water, thus expanding the field of applications beyond the context of urban planning.
The present study employs a comparison between the interval velocity technique and the more conventional constant velocity technique to estimate the depth of metal drums. Numerical modeling studies were performed aiming to validate the studies with real data. The research was developed at the IAG/USP Geophysics Test Site  .
2. Study Area: IAG/USP Geophysical Test Site
The investigated area was the Shallow Geophysics Test Site located in front of the Institute of Astronomy, Geophysics and Atmospheric Sciences (IAG) of the University of São Paulo (Figure 1).
Figure 1. Location of study area and the lines that form the IAG test site. The targets buried in Line 4 (L4) are shown in blue, orange and black, and are distributed along the red line (Adapted image from Google Earth).
drums, representing the context of chemical waste discards. In Line 5 has empty metallic drums filled with water and brine, simulating contamination environments. Finally, in Lines 6 and 7 contain metallic pipes, electrical cables and PVC conductors, representing water, gas, electricity and telephone networks.
Precise information is recorded for each line: depth, target position (horizontal or vertical) and the types of fluids (fresh or salt water, bryne or air). This information can be used for testing geophysical methods, thus simulating real problems and conditions. Another important aspect of this test site is that it is located in an urban context, offering the presence of noise from electric power grids, automobile traffic and people, and the presence of electronic devices. This interference is routine for data acquisition, and permit the verification of noise for instrumentation uses at the site. Further details about the IAG/USP test site construction can be obtained in  .
The present work was carried out on Line 4, consisting of metallic drums (Figure 1), which simulate environmental problems related to the disposal of waste. In New Jersey  , the GPR was employed to find metallic drums in a controlled test site where these drums were associated to dumps of chemical rejects. Similarly, in Rome  GPR was applied to find and to characterize abandoned dumps. Metallic targets were selected because they were easily identified by GPR, thus avoiding ambiguities in the interpretation of the results.
The targets buried under Line 4 used previously by Porsani et al.  to characterize the GPR signal reflected by the metal drums using a 270 MHz antenna in conjunction with numerical modeling studies. Porsani and Sauck  identified artifices in the radargrams, generated by the reflection of the electromagnetic wave between the metal drums positioned close to each other. The researchers proposed a data processing routine to remove these reflections. Porsani et al.  presented a study about the effect of the orientation of the transmitting and receiving antennas (unshielded) to map the buried targets. To accomplish this, the antennas were placed both parallel and perpendicular to the GPR profile direction with the advantages and disadvantages of each assessed for mapping the desired targets.
We apply the interval velocity technique of Dix  to estimate with greater precision the depth of the metallic drums buried under Line 4, as an alternative to the conventional constant velocity technique applied for the same purpose. A good estimate of the velocity of the electromagnetic wave depends on the degree of knowledge of the lithology and the water content in the soil.
The geology of this area is formed by a clay sandy embankment of dark red color, approximately 3 m thick. Just below this stratum there is a predominance of sandy clay sediments from the Resende and São Paulo formations, with a thickness that can reach up to 53 m, passing to the granitic gneiss basement. This information was verified through systematic geophysical studies in this area, together with information from three monitoring wells located at the site  .
3.1. Numerical Modeling
Numerical modeling studies consist of computational simulation of the propagation of electromagnetic waves in a soil model that hosts metallic drums, reproducing the controlled environment where the real data was acquired. The synthetic data are used with the same objective of the real data, specifically, a comparative study between constant and interval velocity techniques to estimate target depths, as well as to validate the studies done with real data.
Comparing the velocity analysis techniques in a computational environment is owes to the fact that hyperbolas present in synthetic radargrams are very well defined, being the characteristic responses of targets present in the subsurface where the depths of each target are known. This study in a controlled environment with real and synthetic data allows greater control of interpreter error by reading double time and adjusting the hyperbola equation in each target identified in the radargrams. In this way, the analysis of accuracy between the techniques can be done in an appropriate way.
The proposed numerical model is a computational simulation of the propagation of electromagnetic waves in a profile 30 m long by 2.5 m deep, containing empty metallic drums arranged horizontally and vertically (Figure 2). The targets simulate metallic drums in cross-section, appearing as circles and rectangles, where the radius of each circle corresponds to the top radius and base of the cylinder representing the drum, and the width and length of the rectangle are equivalent to the diameter and height of that cylinder, respectively. Parameters for the metallic drums were set as follows: electrical conductivity of σ = 109 S/m, magnetic permeability of μ = 1 and dielectric constant of ε = 1  .
The chosen model is formed by a sandy clay soil characteristic of the study area with the dielectric constant varying in depth, where this gradient of dielectric constant is justified by differences in lithology and soil compaction. Besides, the rainy season in the data acquisition phase, also have influences in the variation of water content in the first meters of soil. The dielectric constant chosen varied from 12 to 18 at 2.5 m depth (Figure 2). These values were obtained by the hyperbolic adjustment technique applied for each target in the studied profile. More details about these studies are in 4.2 item of this paper.
Although these variations also influence the electrical conductivity and the magnetic permeability of the medium, constant values were attributed to these parameters, being σ = 10−3 S/m e μ = 1, because the influence of the dielectric constant is more significant as regards the velocity of the electromagnetic wave in the medium in relation to the electrical conductivity and magnetic permeability (Equation (1)).
Information on the physical properties associated with the local geology was taken from Porsani et al.  and information on the arrangement of metal drums installed in the IAG/USP test site was obtained in Porsani et al.  . The choice of this model aims to reproduce the real environment where GPR data were acquired, and thus compare the real data with synthetic data.
Figure 2 shows the model used in the numerical simulation with a frequency of 270 MHz resulting in a synthetic radargram. The frequency chosen in the simulation is the same as that used in the field work.
The numerical modeling was performed using ReflexWIN® processing program  . The software is based on the Time Domain Finite Differences (FDTD) method to resolve the Maxwell equations, following numerical stability criteria  . In the numerical simulation, a model involving the geological background and the buried metal drums mentioned above was constructed. The
Figure 2. 2D model of soil that hosts metal drums buried horizontally and vertically, with vertical variation of dielectric constant ranging from 12 to 18 in 2.5 m depth.
simulated responses were obtained by using the spacing between the traces of 0.02 m (2 cm) and the frequency of the 270 MHz electromagnetic wave. The frequency used in the numerical model is the same used in the acquisition of real data.
3.2. Real Data
To apply the constant velocity technique to estimate the depth of subsurface targets, it is first necessary to obtain the velocity of the medium in question. This estimation can be done based on dielectric constants of soil lithology that are associated with velocity values already known in the literature (Equation (1)), being a more susceptible method to uncertainties, because it is mean values for each material, not necessarily representing the medium where the work is done. Another way is to use employ Equation (2) to calculate the velocity of the medium using the double time information provided by the GPR referring to some known target present in the subsurface that can be identified in radargrams, The third way is based on hyperbolas present in radargrams (Figure 3), coming from underground targets. Once these hyperbolas are identified, one can adjust the mathematical equation of a hyperbola using Equation (3) to obtain the velocity of the medium, being the approach used in this work. There are also other ways to determine the velocity of the medium, such as the CMP and WARR techniques  , not discussed in the present text.
The relation between the velocity of the electromagnetic wave in the medium (v) and the dielectric constant (ε) is given by Equation (1):
where (c) is the speed of light in the free space.
The velocity can also be defined according to Equation (2). Since (h) the depth of the target in the subsurface and (tdouble) is the transit time (two-way travel time) of the electromagnetic wave in the medium, obtained through GPR, given by:
Figure 3. Relation between the position (X1, X2 and X3) of the GPR profile in relation to the center of target projected on the surface and the depth (h) of the target. The distance between the GPR system and the target determines the travel time of the electromagnetic wave through the subsurface, which characterizes the hyperbolic shape.
Figure 3 shows the relationship between the position (x) of the GPR system along a profile with respect to the center projected on the surface of a target located at depth (h), (t) is the transit time of electromagnetic wave that travel from transmitting antenna, reflects on the target and then returns to the receiving antenna, thereby illustrating how the hyperbolas identified in radargrams are constructed.
Equation (3) shows that the hyperbola (Figure 3) is constructed from geometric relations between the position of the GPR system in relation to the center of target in subsurface and the velocity of the medium, which can be obtained by adjusting this equation in hyperbolas present in radargrams, being the approach used in this work.
Figure 4 shows the sketch of a radargram, where the vertical axis is the double transit time of the electromagnetic wave and the red hyperbola is the type of response for a point target located at depth (h). The velocity can be obtained by using the information from that target, which can be done through the previously known depth information along with the double time provided by GPR (Equation (2)), or by adjusting the mathematical equation of hyperbola Equation (3), the latter being the most usual way in GPR practices.
Since the velocity of the electromagnetic wave is predetermined, the constant velocity technique assumes that such velocity is constant for the entire GPR profile (Figure 4), the depth of the targets of interest can be obtained using the double time reading on the radargram for the target of interest by means of Equation (2).
For the GPR, it is considered that the wave velocity is related to the dielectric constant of the medium (Equation (1)), being dependent on the materials in the
Figure 4. Sketch of radargram containing hyperbolas from targets in subsurface. Since velocity (v) is determined, it is assumed that velocity (v) is the same and invariant for the remainder of the radargram.
Figure 5. Sketch of radargram containing hyperbolas at different depths. For each hyperbola, the velocity (v) is determined by the mathematical adjustment of the hyperbola, the double time through the direct reading of the top of each hyperbola, and finally the width of the layers (h), being defined by the vertical distance between the identified targets.
subsoil and the water content   . The limitation of Equation (1) owes to the assumption that the velocity of the medium is constant for every area. This procedure may not be a good strategy, especially if there are marked variations of lithology and/or water content in the study area, which may occur mainly in large areas.
As an alternative, the concept of interval velocity developed by Dix  , initially developed for seismic reflection investigations, can be used. This principle considers the different speeds with which the seismic waves cross the distances between the interfaces of the geological layers present in the area of investigation. The same concept is used in this work, but instead of geological layers, we consider the vertical distance between the surface and the top of each target. It is assumed that different degrees of soil compaction, variations in lithology and water content occur predominantly vertically. Figure 5 shows the schematic of how these layers are defined. Note that each target provides velocity information, which is associated with every ground packet. The deeper the target is located more variations can be present. The hyperbolas illustrated in Figure 5 are the same as discussed previously and associated with buried targets. The different colors and apertures of the hyperbolas show that speeds are different due to the higher water content in the materials or variation of lithology with the depth. In the vertical axis (t) is the double time, (v1, v2 and v3) are the velocities of the soil package above each target and (h1, h2 and h3) are the width of the layers defined by the vertical distances between the identified hyperbolas.
Once velocities and path times have been obtained, the layers are delimited, and then the widths of those layers are calculated by Equations (4)-(6). Note that layer width is determined from the shallower targets to the deeper ones, because in this way, the velocity variations are smoothed.
The depths (H) of each target are calculated by Equations (7)-(9):
The percentage error of the depth estimate is calculated by Equation (10). This error is how far the depth obtained is far from the real depth, serving as a parameter of control and comparison between the techniques approached in this work.
Note that to employ the interval velocity technique in radargrams to estimate the depth of the buried targets it is necessary to have the GPR reflection of at least two targets located at different depths. Another consideration is that the higher the number of targets located at different depths, the more layers can be defined, and consequently more velocity information is added.
3.3. Data Acquisition and Processing
The GPR data were acquired by means of the reflection profile technique with constant spacing, using the SIR4000 (GSSI) equipment with a shielded antenna of 270 MHz. A profile 30 m long with spacing between the traces of 0.02 m was acquired on Line 4 of the IAG/USP test site, where empty metallic drums were precisely buried, arranged horizontally and vertically, with depths varying between 0.5 and 2.0 m (Figure 2).
After the data acquisition step, the processing was done to increase the signal to noise ratio. Data processing included several steps with the objective of improving the radargram reflections. The first step was zero time correction, where the arrival of the air wave at the initial time of each trait is attributed, necessary so that the depth estimation of the targets is not impaired. Frequency filters (bandpass) were applied to eliminate high and low frequency noise, based on the spectrum of the signal obtained and the center frequency of the antenna used. Temporal gains were applied to highlight the anomalies of interest and compensate the losses by geometric scattering of the electromagnetic wave. Background removal was applied in order to minimize the effects of the background and evidence the hyperbolas of interest.
4. Analysis of Results
4.1. Synthetic Data
Figure 6 shows the radargram at a frequency of 270 MHz referring to the model shown in Figure 2. The velocity of the electromagnetic wave in the medium in question was determined for all the targets by means of hyperbola adjustment (Equation (3)).
In Figure 6 the diffraction hyperbolas of the targets are clearly identified, since it is a synthetic radargram. The hyperbolic adjustments for horizontally arranged metal targets (A - F targets) have perfect hyperbolas, and consequently their velocities are easily obtained by adjustments. On the other hand, targets G and H are metallic drums arranged vertically, generating diffraction hyperbolas flattened at the top, due to the positions of the ends of the drums in relation to the GPR profile, which generates greater difficulties in adjusting the hyperbolas and consequently obtaining the velocities for these targets.
Once velocities have been determined for each target by means of hyperbolic adjustments, the layers are defined considering the tops of the drums located at different depths and laterally proximate. Taking into account the position of the targets identified in the radargram (Figure 6), 5 models of 3 layers each were defined, composed by near-side targets and located at different depths.
The first model was elaborated using targets A, B and C, since these are located at different depths and close to each other (Figure 7). Once the layer model is established and the velocities are known, the widths of the layers (Equations (4)-(6)) and hence the depths (Equations (7)-(9)) of the targets involved in this model are calculated. The second model was constructed using targets A, C and D. Note that depths of targets A and C have already been obtained in the first model, i.e, in this second model depths A and C are obtained again as well as the depth of target D. The third model consists of targets D, E and H, where the depths of the targets E and H are obtained, since the depth of the target D was determined in the second model. The fourth model was constructed taking into account targets D, F and H, with the depth of only target F calculated, since the depth of the targets E and H had been calculated previously. In the fifth and final model, targets D, G and H were used to define the layers, with only the depth of target G determined as the depths of D and H had already been determined.
The depths of some targets are obtained more than once, since some targets were part of more than one model. An arithmetic mean is assigned to the depths obtained more than once, as the final depth of the target in question.
In Figure 7 the model of layers defined for the targets A, B and C is presented, where VA, VB and VC are the speeds determined by the hyperbola adjustment (Equation (3)) and h1, h2 and h3 are the widths of each layer defined by hyperbolic adjustments.
Table 1 presents the velocity values obtained by adjusting the hyperbola for each target, the actual depths estimated by the interval velocity technique and the errors associated with the estimates calculated using Equation (10).
Table 1 shows that the estimates of target depths using the interval velocity technique are close to the exact values. Since the model created (Figure 2) for the numerical simulation counts on the vertical variation of the dielectric constant with depth, consequently the deeper targets A and H (~2 m) are more sub-
Figure 6. Synthetic radargram referring to the proposed model and the hyperbolic adjustments with the respective velocities.
Figure 7. Synthetic radargram with the layer model defined for targets A, B and C.
Table 1. Velocity values obtained, real depths, depths estimated by the interval velocity technique and associated errors.
ject to velocity variations than the other targets. However, the discrepancy between the estimated and actual depth for these targets is less than or equal to 4 cm, being an indication that the interval velocity technique is adequate to bypass the effect of the vertical variation of dielectric constant with depth.
In an analogous way, in Table 2 the depth estimates of the determined targets by means of the constant velocity technique are presented. For this purpose, the velocity of target C was adopted as reference velocity for all the targets, and the depth estimated obtained using Equation (2). Target C was used as a reference because a greater number of targets are buried at the same depth.
In Table 2 the discrepancies between the estimated and actual depths reach 2 cm for targets B to G. Note that in Figure 6 the velocity values obtained between these targets are very close due to the fact that they are located at similar depths, thus justifying the small margin of error attributed to the depth estimates for
Table 2. Values of the real depths, and depths estimated by means of the constant velocity technique and associated errors.
these targets by the constant velocity technique. For the deeper targets A and H this discrepancy reaches 12 cm, which reinforces the thesis that the marked variation of dielectric constant poses difficulties in applying the constant velocity technique to estimate the depth of buried objects.
Similar to that observed in the analysis of the real data, when comparing the percentage results presented in Table 1 and Table 2, it is possible to verify that the interval velocity technique applied to estimate target depths is less susceptible to errors in relation to the constant speed technique. Table 2 shows also that the analysis of targets B to G indicates that the constant speed technique can be applied to estimate the depth of underground utilities with precision, since in the study area there are no significant variations of dielectric constant.
4.2. Real Data
Figure 8 shows the result of the GPR 270 MHz profile acquired along Line 4 of the IAG/USP test site. The diffraction hyperbolas for the buried metal drums were adjusted by means of Equation (3).
Observe in Figure 8 that the horizontal and vertical metal drums and the metal pipe reference guide (target D) have clear diffraction hyperbolas. Note also that the metal drums arranged vertically (G and H) do not have perfect hyperbolas due to their extremities. In this case, it is difficult to adjust the hyperbolas.
For target H it was possible to visualize the diffraction hyperbola and then make the hyperbolic adjustment to obtain the velocity. The same does not occur for target G, where the diffraction hyperbola is confused with the geological background, making it difficult to visualize the hyperbola, and consequently to make the hyperbolic adjustment. To overcome this problem, it was assumed that the target velocity L is equal to target F, since they are located close together and at similar depths.
The layer models defined for the actual data are the same as the templates defined for the synthetic data. Once the layers have been fixed and the velocities established the widths of the layers and finally the depths of the targets are cal-
Figure 8. Radargram of the profile of metallic drums and respective velocity adjustments by means of the hyperbola equation.
Figure 9. Layer model defined for A, B and C targets.
Table 3. Velocity values obtained, real depths, depths estimated by the interval velocity technique and associated errors.
culated. Figure 9 shows the defined layers for targets A, B and C, where VA, VB and VC are the velocities determined by hyperbola adjustment (Equation (3)) and h1, h2 and h3 are the defined layer widths.
After the layers are established and the target depths determined, Equation (10) allows the comparison between the estimated and actual depths. Table 3 shows the results obtained using the interval velocity technique applied in real data.
Table 3 shows that the discrepancies between the estimated depths and the actual depths are less than or equal to 4 cm for targets B to G. For targets A and H the difference between the estimated depth and the actual depth reaches 7 cm. Note that targets A and M are deeper (~2 m) and, accordingly, are more susceptible to velocity changes.
For the use of the constant velocity technique, the velocity of target C was set
Table 4. Actual depth values, estimated depths using the constant velocity technique and associated errors.
as the velocity of the medium, following the same procedure as applied to synthetic data. Table 4 shows the depth values estimated by means of the constant velocity technique applied in synthetic data.
In Table 4, the depths estimated for targets B to G show differences of up to 13 cm in relation to the real depths. The estimated depth of target A is 25 cm off the real depth, and the estimated depth of target H is 19 cm off the real depth. These discrepancies are indicative that the soil may have heterogeneities of lithology or variations of water content.
Comparing the percentage results presented in Table 3 and Table 4, the difference between the estimates of target depths through constant velocity and interval techniques is significant, indicating that the approach using the interval velocity technique is less susceptible to errors in the depth estimates.
Estimates employing the depth interval velocity technique appear more suitable and less susceptible to errors, since they are conducted using proximate targets and more velocity information compared to constant velocity technique. Additionally, the depth estimation by the constant velocity technique is more susceptible to errors, as shown in the analysis of real GPR data. These errors can mean potential hazards associated with installing and/or repairing underground utilities.
The studies of the real and synthetic data permit quantifying the vertical variation of depth associated with the variation of the dielectric constant, demonstrating the limiting effects of the constant velocity technique to obtain the depth of buried utility structures and showing that the use of the interval velocity technique can be applied satisfactorily to overcome this problem.
Since the dielectric constant depends on several parameters, such as soil compaction, lithology and water concentration, the results attained in this study indicate that application of the interval velocity technique where there are several buried utilities provides more detailed information about the velocity of the medium. This method has positive implications for both utility depth estimates, knowledge of soil matrices and can be used beyond urban planning contexts, such as for archaeological excavations in urban contexts.
BP is grateful CNPq-Conselho Nacional de Desenvolvimento Científico e Tecnológico for providing a research scholarship (134647/2015-7). JLP also is grateful CNPq (grants: 301692/2013-0 and 406653/2013-5) and FAPESP-Fun- dação de Amparo à Pesquisa no Estado de São Paulo for providing ﬁnancial support for the construction of the IAG/USP Geophysics Test Site (2002/07509-1), both are Brazilian research agencies. IAG/USP is acknowledged for providing infrastructure support. We thank Ernande, Marcelo and colleagues for helping in geophysical data acquisition.