GPR to constrain ERT data inversion in cavity searching: Theoretical and practical applications in archeology

GPR to constrain ERT data inversion in cavity searching: Theoretical and practical applications in archeology

Journal of Applied Geophysics 89 (2013) 35–47 Contents lists available at SciVerse ScienceDirect Journal of Applied Geophysics journal homepage: www...

6MB Sizes 2 Downloads 40 Views

Journal of Applied Geophysics 89 (2013) 35–47

Contents lists available at SciVerse ScienceDirect

Journal of Applied Geophysics journal homepage: www.elsevier.com/locate/jappgeo

GPR to constrain ERT data inversion in cavity searching: Theoretical and practical applications in archeology Luciana Orlando Dept. DICEA, Sapienza University of Rome — Italy

a r t i c l e

i n f o

Article history: Received 7 March 2012 Accepted 5 November 2012 Available online 22 November 2012 Keywords: ERT GPR Cavity A priori constrain Archeology

a b s t r a c t I used theoretical forward models to show that a cavity embedded in a stratified sedimentary sequence can induce an equivalence problem in the ERT data inversion. Conductive top soil increases the misfit between the ground feature and the ERT model. The misfit depends on array and stratigraphy sequences. The latter induce an equivalence problem that manifests itself as wrong cavity depth positioning. The misfit is greater in the data acquired with Schlumberger array than with dipole–dipole. The ambiguity of ERT data inversion problems was tested in the detection of cavities linked to an 8th–6th century B.C. Sabine tomb, 3 m wide × 3 m long × 2 m high, excavated from a shaly gray volcanic ash (cinerite) layer covered by semi-lithoid tuff and top soil layers. In the real study I reduced the ambiguity in the inverse problem of ERT data using a priori information on geometry and resistivity of the cavity. The constrains were carried out from georadar data acquired with 80 and 200 MHz antenna. I demonstrate that this procedure has a practical application in cavity detection, and is a key to the reduction of the uncertainty inherent in the inversion process of ERT data. © 2012 Elsevier B.V. All rights reserved.

1. Introduction The refinement of methods devoted to investigating cavities is one of the objectives of applied science, as the knowledge obtained from such methods can reduce the hazards that cavities induce in environmental, industrial and ancient buildings, and can help to detect voids related to buried archeological structures. Geophysical methods are an important tool for the detection of near-surface cavities, and enhance the probability of detecting cavities. The most effective methods in cavity exploration are microgravity (Buttler, 1984; Colley, 1963; Rybakov et al., 2001), electrical resistivity tomography (Bishop et al., 1997; McDonald, 2003; van Schoor, 2002), georadar (Benson, 1995; Malagodi et al., 1996; Piro and Gabrielli, 2009), self potential (Quarto and Sciavone, 1996; Vichabian and Morgan, 2002), and seismics (Cook, 1965; Grandjean, 2006; Grandjean and Leparoux, 2004; Piro et al., 2001; Sheehan et al., 2005). The efficiency of geophysical methods depends on the contrast between the physical properties of the material filling the cavity and those of the host material: the greater the contrast the more effective the method. Success also depends on the environmental setting of the area, and on the dimension and depth of the target in relation to the resolution of the applied method. For example, GPR gives high resolution but has been found more suitable for resistive environments and for mapping shallow targets; Gravity is sensitive to the dimension, density and depth of voids; ERT has proven advantageous for the detection of cavities filled by material that provides a large resistivity contrast with respect to the

E-mail address: [email protected]. 0926-9851/$ – see front matter © 2012 Elsevier B.V. All rights reserved. http://dx.doi.org/10.1016/j.jappgeo.2012.11.006

host material, but it suffers resolution limitations (Ellis and Oldenburg, 1994). Indeed, the application of several geophysical methods is required to set up a model of a subsurface to correctly depict the actual ground setting, and in most cases the data of each method is first inverted individually before being interpreted jointly. An excursus on the detection of cavities induced by sinkholes, using seismic refraction tomography, MASW and the joint interpretation of georadar and electrical resistivity tomography data, can be found in Dobecki and Upchurch (2006). Also Kaufmann et al. (2011) interpreted gravity and resistivity data jointly, as did El-Quadi et al. (2005) with their data from geoelectrical resistivity tomography and georadar, and Mochales et al. (2008) with data from gravity, magnetic and radar methods used to detect cavities. Small shallow cavities are detected successfully by combining GPR and ERT data. Indeed, the joint interpretation of GPR and ERT can help to remove ambiguity from data interpretation as GPR enables the reconstruction of the cavity in terms of geometry, and ERT in terms of resistivity. In fact, GPR can depict the top and bottom of a cavity with high resolution, but it often fails to give information about the material forming the structure, thus the data interpretation tends to be based on anomaly features. Instead ERT, which has less resolution than GPR, gives information on the lithology of host material and that of the target. However the two methods are limited by the low penetration of GPR in conductive media, and by the resolution and ambiguity problems of ERT. My approach to the problem has been to reduce the ERT data inversion ambiguity using a priori information on the boundary of the cavity obtained from GPR data. The actual data refer to an archeological necropolis characterized by several scattered Sabine tombs (8th–6th century B.C.). I focus the analysis on a chamber with a dromos tomb filled by air.

36

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

The joint interpretation of the archeological, geological, ERT and GPR data from this site showed that the ERT data inversion, depending on the array configuration, can be affected by non-uniqueness of solution because the inversion gave the incorrect depth of the tomb, the anomaly being, for the most part, incorporated in the tuff layer above the tomb. In the studied case the ambiguity found in the data inversion was due to the tomb merging with a more conductive media and being covered by resistive and outcropping humus layers. It is well known that in particular conditions all geophysical methods produce uncertainty in data interpretation because of the ambiguity inherent in data inversion (Langer, 1933). This occurs because field data are inverted by the linearization of the forward model and the accepted final solution is given by the minimization of the misfit function between calculated and experimental data (Sharma and Verma, 2011). In most cases global optimization techniques search for the solution in a predefined search range, and the system is free from constraints. It is well known that the error function includes several minima, and each minima corresponds to a model that can fit the observed data but not depict the actual ground setting. Inevitably, when a non-unique solution in the data inversion occurs, a single model to interpret the observations is sought. This can create ambiguity in the geophysical interpretation and quite often provides misleading results. Recently, considerable effort has been put into the use of joint inversion of geophysical data to reduce uncertainty in the data inversion. Geophysical data that are sensitive to the same physical quantity, for example magnetotelluric and DC resistivity, can be simultaneously inverted by minimizing an objective function that includes the data misfit of the different data types. Joint magnetotelluric-resistivity inversion was applied by Vozoff and Jupp (1975), Constable et al. (1987), and Albouy et al. (2001). In the case of data sensitive to different physical quantities, the joint inversion becomes more difficult because the methods involve different parameters that may not be directly correlated. Moreover the resolution of the resulting models varies throughout the model, and is different for different geophysical techniques (Day-Lewis et al., 2005; Gallardo and Meju 2003). Thus, a fundamental problem is to establish a precise and robust indirect relationship among the different characteristics of the media. One way to stabilize the solution is to impose the ground's geometric or physical characteristics, which can be derived from direct and indirect methods. Indeed, Beres et al., 2001 used GPR data as input to constrain the gravity data inversion in the cavity detection, and de Nardis et al. (2005) and Cardarelli et al. (2010) applied joint inversion to resistivity and seismic data. Amidu and Dunbar (2008) reduced uncertainty by the inclusion of a priori information on water resistivity in the inversion scheme, and Scott et al. (2000) inverted the ERT data constraining the boundary between a gravel and underlying clay layer known from the seismic refraction survey. The simplest way to overcome the ambiguity problem of data inversion is to impose a priori constraints on the boundary and/or physical parameters of strata (Musil et al., 2003). In the discussed case, I show that the ambiguity in the solution of ERT data can be overcome constraining the inversion with GPR data. Since surface GPR reflection was effective in detecting cavity boundaries, I used this structural information, derived from GPR in the inversion of surface ERT data, to provide the correct position of the resistive bulk of a cavity linked to a tomb. The same procedure was followed by Doetsch et al. (2012) and Linde et al. (2006) who constrained 3-D electrical resistance tomography with GPR data for aquifer characterization. The ambiguity of the ERT data inversion, and the validation of the procedure, was tested by synthetic models that took into account equivalence problem and array type. 2. Theoretical model The non-uniqueness of the solution of ERT data inversion was analyzed theoretically by simulating the archeo-geological scenario, which was taken into account during the actual data analysis.

The steps were: 1) Computation of the apparent resistivity based on geometric models, constitutive parameters and acquisition geometries similar to the actual data. 2) Data inversion of theoretical apparent resistivity using the same code used in the data inversion of the actual data. The real case was re-sampled using two synthetic models (Figs. 1–6). The models differ in the outcropping conductive layer simulating the top soil layer. In the model of Figs. 1–3, the cavity is covered by a resistive layer while in Figs. 4–6 a top soil layer is also included. The empty cavity is simulated by a block 2.5 m wide×2 m high with a resistivity of 30,000 Ω·m. The cavity is merged into a layer with 30 Ω·m resistivity. In the model of Figs. 1a, 2a, and 3a the empty cavity is covered by a 1 m thick layer with 200 Ω·m resistivity, in the model of Figs. 4a, 5a, and 6a the empty cavity is covered by two layers, the deeper one 1.5 m thick with a resistivity of 200 Ω·m and the shallower 0.5 m thick with resistivity of 30 Ω·m. The layer, 200 Ω·m resistive, simulates a semi-lithoid tuff and that, 30 Ω·m, the top soil. The cavity is carved into a layer 30 Ω·m resistive. The constitutive parameters of the models were selected taking into account the sedimentary sequence of actual data. Because the quality of data inversion is also influenced by the array type (Orlando et al., 1986) and electrode distance, the simulations were carried out for the dipole–dipole and Schlumberger arrays with electrodes 0.5, 1 and 2 m apart. The theoretical model was calculated with RES2DMOD program and the data inversion with RES2DINV program both commercialized by Geotomo. The theoretical apparent resistivity is calculated using the finite–difference method (Dey and Morrison, 1979; Loke, 1994). A random noise of 5% was added to the theoretical apparent resistivity. The code uses the smoothness-constrained least-squares method (deGroot-Hedlin and Constable, 1990) to link the model parameters to the model response, and the complete Gauss–Newton technique to determine changes in model parameters (Loke and Barker, 1996). The theoretical models (Figs. 1–6) show that the lateral limits of the cavity are properly identified in both models regardless of the device and electrode distances. The position of the top of the void of the model of Figs. 1–3 is correctly detected by both the arrays and by all the electrode distances, while the position of the bottom is best located by the dipole–dipole array rather than by Schlumberger, which underestimates it. The 200 Ω·m resistivity layer covering the tomb is properly estimated by both arrays with an electrode distance of 0.5 m (Fig. 1b,c), and less by electrodes 1 m and 2 m apart (Figs. 2 and 3). This layer is almost completely undetected by the Schlumberger array with electrodes 2 m apart (Fig. 3). The 30 Ω·m resistivity layer is properly defined by both the arrays and electrode distances while the resistivity of the tomb is always underestimated. The top of the layer at 200 Ω·m, of the model of Figs. 3–6, is properly located by both arrays, while the bottom is well located with the dipole– dipole array than the Schlumberger array that overestimates the depth. The top position of the cavity is always underestimated and is usually located with the top of layer, 200 Ω·m resistive. Both electrode arrays, dipole–dipole and Schlumberger, with electrodes 0.5 m apart, show the tomb almost completely positioned within the layer 200 Ω·m resistive. The results achieved with electrodes 1 and 2 m apart are better than that with electrodes 0.5 m apart. The 30 Ω·m resistivity is correctly defined by both arrays and with all electrode distances. The resistivity of the layer covering the cavity (200 Ω·m) is well defined by both arrays with electrodes spaced 0.5 m apart, while it is overestimated with electrodes 1 m apart and underestimated with electrodes 2 m apart. Also in this model the resistivity of the cavity is always underestimated. In the models with electrodes 2 m apart (Fig. 6) a ghost layer with resistivity b30 Ω·m and thickness similar to the 200 Ω·m resistive layer is present at the bottom of the 200 Ω resistive layer. This phenomenon is probably due to the distances between the electrodes (2 m), which are too wide to resolve this underground feature.

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

37

Distance (m) Depth (m)

0

0

15

25

20 200 ohm.m

30000 ohm.m

2

30 ohm.m

a)

4

Depth (m)

10

5

4

2

0

6

8

10

12

14

16

18

10

22

24

14

16

18

20

22

24

1

b)

2

Calculated apparent resistivity 25

0

Depth (m)

50

100

2

200

400

4

800 1600 3200 ohm.m

6

8

10

12

2

c)

4

Calculated apparent resistivity 25

50

100

200

400

800 1600 3200 ohm.m

Fig. 1. Theoretical model 1: a) Geometries and resistivity of the model; b) calculated resistivity for a dipole–dipole array; and c) calculated resistivity for a Schlumberger array. The array is 25 m long with electrodes 0.5 m apart.

The model of Figs. 1–3 is properly reconstructed by both arrays with electrodes spaced 0.5 and 1 m, even though the dipole–dipole provides a more accurate reconstruction of the cavity with respect to Schlumberger. The model of Figs. 4–6 is not properly represented by both arrays in terms of depth of the top and bottom of cavity, which are always underestimated. The phenomenon is more pronounced for

the electrodes 0.5 m apart. (Fig. 4b,c) where the cavity is mostly located in the resistive layer. The greater part of the anomaly linked to the cavity is, in both cases, embedded in the resistive layer simulating the lithoid-tuff. The theoretical models show that under this stratigraphic condition a non-uniqueness problem occurs and the solution model is further

Depth (m)

Distance (m) 0

0

5

10

15

20

25

30

35

40

45

200 ohm.m 30000 ohm.m

2

30 ohm.m

a)

4

5

10

15

20

25

30

35

40

45

25

30

35

40

45

Depth (m)

0 2 4

b) Calculated apparent resistivity 12

Depth (m)

0

0

25

50

100 200 400 800 1600 3200 ohm.m

5

10

15

20

5

c) Calculated apparent resistivity 25 50

100

200

400

800 1600 3200

ohm.m

Fig. 2. a) Geometries and resistivity of the model; b) calculated resistivity for a dipole–dipole array; and c) calculated resistivity for a Schlumberger array. The array is 53 m long with electrodes 1 m apart.

38

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

Distance (m) Depth (m)

0

20

40

60

80

0 200 ohm.m

30000 ohm.m 30 ohm.m

a)

5

20

40

60

80

60

80

Depth (m)

0

5

10

b) Calculated apparent resistivity 12

25

50

100

0

Depth (m)

200

400

800 1600 ohm.m

40

20

5

c)

10

Calculated apparent resistivity 25

50

100 200 400 800 1600 3200 ohm.m

Fig. 3. a) Geometries and resistivity of the model; b) calculated resistivity for a dipole–dipole array; and c) calculated resistivity for a Schlumberger array. The array is 100 m long with electrodes 2 m apart.

from the actual one when an outcropping conductive layer (Figs. 4–6) is present. The theoretical models have shown that this problem can only be partly solved by electrode array. This kind of ambiguity in ERT data inversion can be eliminated by imposing a priori constraints on the cavity's boundary and resistivity in the data inversion. In this study a priori constraints on the boundary of the tomb in the ERT data inversion were carried out by direct observation and GPR data.

3. Geo-archeological setting of the area The actual data discussed here were collected in a necropolis containing several 8th–6th century B.C. Sabine tombs. The surveyed area, located about 20 km north of Rome (Italy) (Fig. 7a), was partially excavated at the beginning of the 20th century (1915-16), revealing fossa tombs and some chamber tombs with entrance corridors (“dromos”, Buranelli et al., 1997) (Fig. 7b). Currently, nothing of the ancient necropolis is visible except for dromos chamber tomb n. 809

Distance (m) Depth (m)

0

0

10

15

20

25

200 ohm.m

2 30000 ohm.m

30 ohm.m

a)

4

Depth (m)

5

2

0

4

6

8

10

12

14

16

15

20

22

24

12

14

16

18

20

22

24

1 2

b) Calculated apparent resistivity 25

50

2

100

4

400 ohm.m

200

6

8

10

Depth (m)

0 2

c)

4

Calculated apparent resistivity 12

25

50

100

200

400 ohm.m

Fig. 4. a) Geometries and resistivity of the model; b) calculated resistivity for a dipole–dipole array; and c) calculated resistivity for a Schlumberger array. The array is 25 m long with electrodes 0.5 m apart.

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

39

Distance (m) Depth (m)

0

0

5

10

15

20

25

30

35

40

45

200 ohm.m

2 30000 ohm.m

30 ohm.m

a)

4

5

10

15

20

25

30

35

40

45

25

30

35

40

45

Depth (m)

0 2 4

b) Calculated apparent resistivity 25

50

5

100

10

400 ohm.m

200

15

20

Depth (m)

0

5

c) Calculated apparent resistivity 25

50

100

200

400 ohm.m

Fig. 5. a) Geometries and resistivity of the model; b) calculated resistivity for a dipole–dipole array; and c) calculated resistivity for a Schlumberger array. The array is 53 m long with electrodes 1 m apart.

(Buranelli et al., 1997), which, until recently, was used as a farmer's cellar. I used this tomb to validate the data inversion procedure, and to match the old excavation map and the geophysical surveys. Unfortunately, the tomb n. 809 is the only reference point for the old excavation so the overlapping of the old and the new maps was far from accurate. The necropolis lies on a small hill geologically characterized by lithoid tuff, covered by a thin, top-soil layer of variable thickness (0.20 to 0.5 m) (Cammarano et al., 1997) and sedimentary lacustrine formations (clay

and silt). The hill top itself consists of top soil, stratified semi-lithoid tuff and shaly gray volcanic ash (cinerite). Most of the tombs excavated in the early 1900s were of medium size, they were usually “fossa” but there were some chamber tombs (about 3×3×2 m) with entrance corridors (dromos) a few meters long. The fossa tombs, located mainly at the top of the hill, generally belong to the older period of the necropolis (8th century B.C) while the more recent dromos chamber tombs (7th and beginning of the 6th century B.C.)

Depth (m)

Distance (m) 0

20

0 5

60

40

80

200 ohm.m 30000 ohm.m

a) 20

40

30 ohm.m

60

80

60

80

Depth (m)

0

5

10

b) Calculated apparent resistivity 12

25

50

100

200

400 ohm.m

40

20

Depth (m)

0

5

10

c) Calculated apparent resistivity 12

25

50

100

200 ohm.m

Fig. 6. a) Geometries and resistivity of the model; b) calculated resistivity for a dipole–dipole array; and c) calculated resistivity for a Schlumberger array. The array is 100 m long with electrodes 2 m apart.

40

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

b)

a) *

Study area

Fossa tomb

Dromos chamber tomb (n. 809) Fossa with loculus tomb

Topographic contour line

GPR survayed area

ERT line

Fig. 7. a) Study area located about 20 km north of Rome (Italy). b) Fossa tombs and dromos chamber tomb with entrance corridors (“dromos”, Buranelli et al., 1997) excavated at the beginning of the 20th century (1915-16). In gray the topographic counter lines. Location of GPR survey (gray box) and ERT lines (gray lines). The dotted black line number 1 locates the ERT and GPR profiles position of Figs. 9, 13 and 14.

are located on the slopes. The fossa tombs are trenched in the semilithoid tuff. Chamber Tomb n. 809, which lies in the cinerite layer and is covered by semi-lithoid tuff and a top soil layer (Fig. 8), is one of the most important structures found so far in the necropolis. The tombs are scattered at different depths, generally not exceeding 3–4 m.

4. Geophysical surveys Geophysical surveys were carried out above and near Tomb n. 809 (Fig. 7). I analyzed the GPR and ERT data that respectively depict the underground in terms of electromagnetic impedance and resistivity

a)

b) Plan

Section

Chamber Top soil

Semi-lithoid tuff

Tomb entrance

After Buranelli et al, 1997

Fig. 8. a) Photograph of the chamber tomb entrance, excavated in a clay layer and covered by semi-lithoid tuff and top soil. b) Plan and section of the chamber tomb 809 (After Buranelli et al., 1997).

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

a)

0

0

10

Distance (m) 20

41

30

40

200 MHz

Time (ns)

Topographic level

100 Bottom of semi -lithoid tuff

Top of tomb

b) 0 Time (ns)

80 MHz

Top of tomb Topographic level

100 Bottom of tomb

200 Fig. 9. GPR transect collected with 200 MHz (a) and 80 MHz (b) antennas. The location is indicated by the dotted black line in Fig. 7.

variations. I focused my attention on removing the ambiguity in the ERT data inversion. 4.1. GPR The GPR data were acquired in a square area 40 × 40 m wide (Fig. 7). The area was surveyed with 40 profiles spaced 1 m apart with a 80 MHz

antenna, and 81 profiles spaced 0.5 m apart with a 200 MHz antenna. The data were collected in continuous mode along parallel transects stretching NW–SE in a time window of 256 ns. The sampling spatial interval was 0.048 m and 0.024 m for 80 and 200 MHz respectively, and the sampling time interval was 0.5 ns for both. The GPR data processing involved time zero calibration, zero-phase band-pass filtering, and background subtraction. Signal amplitudes were

Dromos chamber tomb (n. 809)

Topographic contour line Fig. 10. Time slice at 0.8 m collected with 200 MHz antenna, superimposing on the topographic contour and the map of the ancient excavation.

42

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

Dromos chamber tomb (n. 809)

Topographic contour line Fig. 11. Time slice at 2.6 m collected with 80 MHz antenna, superimposing on the topographic contour and the map of the ancient excavation.

Limit of semi-lithoid tuff

Fig. 12. Photograph of the limit of the semi-lithoid tuff clearly visible in the dry season.

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

43

a) Top soil

Semi-lithoid tuff

Tomb entrance

Semi-lithoid tuff

Depth (m)

115

Tomb

b)

110

shaly volcanic ash

Model resitivity with topography Iteration 7 RMSerror = 2.09

Depth (m)

115

c)

110 Model resitivity with topography

Unit electrode spacing = 1m

Iteration 7 RMSerror = 0.58

5

10

20

15

25

30

35

40

Calculated apparent resistivity 12.5

25

50

100

200

400

800 ohm.m

Resistivity in ohm.m Fig. 13. Photograph of the entrance of the chamber tomb (a); 2D ERT profile collected on the top of the dromos chamber tomb positioned as the dotted black line in figure acquired with dipole–dipole (b) and Schlumberger (c) arrays respectively. The same color scale was used in the mapping.

recovered by means of spherical exponential gain, and topography correction was also included in the data. The profiles were used to built a radar cube that was used to extract time slices at different depths.

Fig. 9 shows the GPR profiles acquired with 200 and 80 MHz antennas on dromos chamber tomb n. 809 (Figs. 7,8). Both GPR antennas detect the top of the tomb but only the 80 MHz antenna detects the bottom. The 200 MHz antenna detects the lower limit of the semi-lithoid

Semi-lithoid tuff Tomb

Depth (m)

115

110

shaly volcanic ash

Model resitivity with topography

Iteration 7 RMSerror = 2.1

Depth (m)

115

110 Model resitivity with topography Iteration 7 RMSerror = 0.57

5

10

Unit electrode spacing = 1m

15

20

25

30

35

40

45

Calculated apparent resistivity 12.5 25

50

100 200 400 800 ohm.m

Fig. 14. 2D ERT profile of Fig. 13 inverted with a priori constrain on the geometries and resistivity with dipole–dipole (b) and Schlumberger (c) arrays, respectively. The same color scale was used in the mapping.

44

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

tuff with higher resolution than the 80 MHz, and although the penetration depth of the two antennas was very similar higher resolution was achieved with 200, compared to 80 MHz. The 200 MHz highlights that the bottom of the tuff is characterized by several diffraction hyperbolas and detects several small anomalies that, for the most part, I related to fossa tombs. The GPR correctly detected the depth and shape of chamber tomb n. 809. The two-travel-time in-depth conversion was performed by direct measurement of the top soil thickness. The depths of the top and bottom, and the width of the tomb, were used as a priori information in the ERT data inversion. The data were used to build two 3D data cubes for the 200 and 80 MHz, and to extract the time slices at different depths. The time slice at 0.8 of the 200 (Fig. 10), superimposing the topographic contour and the map of the ancient excavation, shows several scattered anomalies located at the top of hill in correspondence with the excavated tombs. Note that not all the GPR anomalies match the tombs of the ancient excavation, this is because of the low georeferencing accuracy of the ancient excavation map with respect to the GPR data. The time slice at 2.4 m of 80 MHz antenna (Fig. 11) shows that at that depth few anomalies are present, and the most significant are scattered on a contour line. The biggest anomaly is linked to the chamber of dromos tomb n. 809 at the right end of the aligned anomalies. These anomalies are aligned parallel to the lateral edge of the tuff visible in the dry season (Fig. 12). 4.2. ERT The area was surveyed with 13 profiles, 1.5 m spacing, using Schlumberger and dipole–dipole arrays, which partially overlaps the GPR survey (Fig. 7): 10 NW–SE profiles and 3 perpendicular to them.

The data were acquired with a Syscal Pro instrument equipped with 48 electrodes, spaced 1 m apart. The 2D and 3D data inversions were obtained using the same code as that used for the theoretical data inversion. The topography, arranged with model blocks and a moderately damped distorted grid, was incorporated in the model. The Schlumberger and dipole–dipole ERT profiles inverted without constrains of Fig. 13 is located on the chamber tomb n. 809 (dotted black line in Fig. 7), positioned as the GPR profile of Fig. 9. The ERT section compared with the photograph (Fig. 13a) shows that the semi-lithoid tuff covering the tomb has a resistivity greater than 150 Ω·m and is laterally truncated as shown by the photograph of Fig. 12. The semi-lithoid tuff layer is characterized by several lateral discontinuities that are better defined with dipole–dipole array than Schlumberger array. My interpretation being that they inferred fossa tombs entrenched in the tuff and nowadays filled with loose conductive sediments. The semi-lithoid tuff layer overlays a 15–50 Ω·m resistive formation that is linked to the shaly gray volcanic ash where the dromos chamber tomb is excavated. In the ERT profile the chamber tomb manifests itself with a resistive anomaly >400 Ω·m at 22–25 m of the ERT profile. That anomaly is partially included in the tuff layer in the profile acquired with Schlumberger spread (Fig. 13c). The lateral limits are correctly defined by both dipole– dipole and Schlumberger arrays (Fig. 13b,c) that give a width similar to the actual one. With both spreads, dipole–dipole and Schlumberger the depth of the tomb is undervalued. The area with resistivity lower than 15 Ω·m seem due to the inversion and not to a real resistivity variation of the cinerite layer. The phenomenon is similar to that observed in the theoretical model (Fig. 6). In this case the ERT does not give a model that well fit the actual ground setting for what concerns the tomb thick and depth for both

Resistivity in ohm.m

Dromos chamber tomb n. 809

Topographic contour line Fig. 15. Tomolayer at 0.5 m obtained from a 3D cube inversion of 10 profiles stretched NW–SE with a priori information on the boundaries and resistivity of tomb.

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

dipole–dipole and Schlumberger array, even though dipole–dipole array gives solution more close to the actual setting than the Schlumberger array; though the RMS misfit error between the observed and calculated apparent resistivity is around 3% for dipole–dipole and 0.6% for Schlumberger array. Clearly this particular stratigraphic sequence generates a problem of ambiguity in the solution, the same as shown by the theoretical models discussed in the theoretical section. Theoretical models have shown that the ambiguity of the solution is further highlighted by the presence of an outcropping conductive layer (top soil), which better fits the actual case study. To overcome ambiguity in ERT data inversion I constrained the data inversion using the tomb geometry obtained from direct observation and GPR data interpretation taking into the data acquired with both dipole–dipole and Schlumberger arrays.

4.3. 2-3D actual data inversion with a priori information To overcome the ambiguity problem that occurred in the ERT data inversion for cavity detection, the ERT data inversion was obtained imposing the a priori constraints on the boundary and resistivity of the tomb. The shape of the tomb was simulated by a rectangle 3×2 m. In this case the code requires the specification of the top–left and bottom– right corners of the rectangle and the resistivity. The coordinates of the rectangle were ascertained by direct observation and from GPR data. I used the Marquardt and Occam inversion that combines the damped least squares method with the smoothness-constrained method, which in this case gives better results (Loke, 1994). The resistivity in that region was imposed to be 2000 Ω·m. The damping factor weight for the

45

resistivity of that region was 2. Such a small value was used because the resistivity can have some degree of uncertainty in small regions. The data inversion was calculated including the topography. Fig. 14 shows the inverted model of the profile acquired on the tomb, depicting the ground as being very similar to the real ground setting. A priori conditions on the tomb geometry and resistivity in the data inversion provide a model very similar to the actual ground setting, positioning the tomb at the true depth. In this case the final inverted model was obtained with 7 iterations, and the RMS misfit error (2.1% for dipole–dipole and 0.57 for Schlumberger arrays) was quite similar to that obtained with no constraints in the data inversion of Fig. 13. Also in this cases zones with resistivity b 15 Ω·m below the tomb are probably due to data inversion and not the real resistivity variation of the cenerite layer. The 10 profiles stretched NW–SE (Fig. 7), acquired with Schlumberger array were used for a pseudo-3D data inversion where a priori information on the boundary of the tomb was imposed. In this case the coordinates of top–left–back corner and bottom–right–front corner of rectangle and the resistivity were imposed. The same code used for the 2D data was used. Also in this case the topography was included in the data inversion. The 3D-ERT cube was used to extract the resistive tomo-layers at different depths. Figs. 15 and 16 show the tomo-layer at 0.5 and 4.5 m from the higher topographic level of the 3D cube. The tomo-layer at 0.5 m depicts the limit of the semi-lithoid tuff that includes several resistive and conductive anomalies, probably inferring tombs. The inverted model constrained with a priori information positions the tomb (Fig. 16) correctly in the space, showing that the procedure adopted allows the easy removal of the ambiguity problem induced by the particular sedimentary sequence of the study area.

Resistivity in ohm.m

Dromos chamber tomb n. 809

Topographic contour line Fig. 16. Tomolayer at 4.25 m obtained from a 3D cube inversion of 10 profiles stretched NW–SE with a priori information on the boundaries and resistivity of tomb.

46

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

Comparing the ERT and GPR profiles (Fig. 17) we observe that we obtained a good overlap of the boundaries of tuff and the positioning of the tomb. Moreover the ERT allows the deposits and the material filling the tomb to be characterized.

5. Conclusion I used GPR and ERT data to reduce ambiguity in the ERT data inversion, and for a joint interpretation of data devoted to characterizing anomalies linked to a tomb. GPR was a useful tool to define the limits of the strata and geometry of the tomb, but it was unable to characterize the strata lithology, which was obtained by ERT. But in this particular sedimentary setting the ERT was affected by ambiguity in the tomb depth positioning with both dipole–dipole and Schlumberger array, even if the dipole–dipole array is more accurate than the Schlumberger array. The georadar reconstructed, in good detail, the main stratigraphic geometries of the ground, including the real dimensions and the position of the top and bottom of the tomb. The 200 MHz antenna detected with high resolution the bottom of the semi-lithoid tuff, the top of the tomb and several anomalies in the tuff that I linked to fossa tombs trenched in the tuff. The 80 MHz antenna, given its lower resolution compared to the 200 MHz, was not suitable for fossa tomb detection

Depth (m)

114

but it did allow the detection of the bottom of the chamber tomb that were used to constrain the ERT data inversion. I observed that the inversion of the ERT data, collected to detect a void that merged into a conductive covering media consisting of two layers (the deeper layer more resistive than the outcropping one), led to an ambiguous solution. In fact, the forward model located the tomb at a position shallower than the real one in the covering resistive layer. The theoretical models show that the ambiguity in the data inversion depends on the stratigraphic resistivity sequence of the ground that only in part can be overcome by array type. Indeed, in the model with top soil covering a more resistive layer the ambiguity problem is enhanced, and to overcome this ambiguity problem inherent in ERT data inversion I used a priori information obtained from GPR data concerning the boundary and the resistivity of the cavity. I found that the 2–3D ERT data inversion constrained with a priori information gives a more robust model of the ground setting and that the dipole–dipole can help in reduction of uncertainty. Acknowledgment I wish to thank Professor Luciana Drago for providing the archeological knowledge and access to the site, and Doctor Valeria Poscetti for her help in the geophysical surveys and implementation of the cartography of the area.

a)

112

108 106

Semi-lithoid tuff

Tomb

110

Model resitivity with topography Iteration 7 RMS error = 2.1 5

12.5 25

50

10

100

200

400

20

15

30

25

35

800 ohm.m

0

Time (ns)

b) Tomb

100

200 MHz 0

Time (ns)

c) 100

Tomb

200

80 MHz 10

20

30

40

Distance (m) Fig. 17. a) ERT profile; GPR profiles collected with 200 MHz (b) and 80 MHz (c) antennas. Data collected on the dromos chamber tomb.

L. Orlando / Journal of Applied Geophysics 89 (2013) 35–47

References Albouy, Y., Andrieux, P., Rakotondrasoa, G., Descloitres, M., Joint, J.L., Rasolomanana, E., 2001. Mapping coastal aquifers by joint inversion of DC and TEM soundings: three case histories. Ground Water 39 (1), 87–97. Amidu, S.A., Dunbar, J.A., 2008. An evaluation of the electric-resistivity method for water-reservoir salinity studies. Geophysics 73 (4), G39–G49. Benson, A.K., 1995. Applications of ground penetrating radar in assessing some geological hazards: examples of groundwater contamination, fault, cavities. Journal of Applied Geophysics 33, 177–193. Beres, M., Luetscher, M., Olivier, R., 2001. Integration of ground penetrating radar and microgravimetric methods to map shallow caves. Journal of Applied Geophysics 46, 249–262. Bishop, I., Styles, P., Emsley, S.J., Ferguson, N.S., 1997. The detection of cavities using the microgravity technique: case histories from mining and karstic environments. Geological Society, Engineering Geology Special Publication 12, 153–166. Buranelli, F., Drago, L., Paolini, L., 1997. La necropoli di Casale del Fosso. Le necropoli arcaiche di Veio, Giornata di studio in memoria di Massimo Pallottino (a cura di Gilda Bartoloni), Roma, pp. 63–83. Buttler, D.K., 1984. Microgravimetric and gravity gradient techniques for the detection of subsurface cavities. Geophysics 49, 084–1096. Cammarano, F., Mauriello, P., Piro, S., 1997. High-resolution geophysical prospecting with integrated methods. The Ancient Acropolis of Veio (Rome, Italy). Archaeological Prospection 4, 157–163. Cardarelli, E., Cercato, M., Cerreto, A., Di Filippo, G., 2010. Electrical resistivity and seismic refraction tomography to detect buried cavities. Geophysical Prospecting 58, 685–695. Colley, G.C., 1963. The detection of caves by gravity measurements. Geophysical Prospecting XI, 1–9. Constable, S.C., Parker, R.L., Constable, C.G., 1987. Occam's inversion: a practical algorithm for generating smooth models from electromagnetic sounding data. Geophysics 52 (3), 289–300. Cook, J.C., 1965. Seismic mapping of underground cavities using reflection amplitudes. Geophysics 30 (4), 527–538. Day-Lewis, F.D., Singha, K., Binley, A.M., 2005. Applying petrophysical models to radar travel time and electric resistivity tomograms: resolution-dependent limitations. Journal of Geophysical Research 110, B08206. de Nardis, R., Cardarelli, E., Dobroka, M., 2005. Quasi-2D hybrid joint inversion of seismic and geoelectric data. Geophysical Prospecting 53 (5), 705–715. deGroot-Hedlin, C., Constable, S., 1990. Occam's inversion to generate smooth, twodimensional models from magnetotelluric data. Geophysics 55, 1613. Dey, A., Morrison, H.F., 1979. Resistivity modelling for arbitrary shaped two-dimensional structures. Geophysical Prospecting 27, 1020–1036. Dobecki, T.L., Upchurch, S.B., 2006. Geophysical applications to detect sinkholes and ground subsidence. The Leading Edge 25, 336–341 http://dx.doi.org/10.1190/1.2184102. Doetsch, J., Linde, N., Pessognelli, M., Green, A.G., Günther, T., 2012. Constraining 3-D electrical resi stance tomography with GPR reflection data for improbe acquifer characterization. Journal of Applied Geohpysics 78, 68–76. Ellis, R.G., Oldenburg, D.W., 1994. Applied geophysical inversion. Geophysical Journal International 116 (1), 5–11. El-Quadi, G., Hafez, M., Abdalla, M.A., Ushijima, K., 2005. Imaging subsurface cavities geoelectric tomography and ground penetrating radar. Journal of Cave and Karst Studies 67 (3), 174–181. Gallardo, L.A., Meju, M., 2003. Characterization of heterogeneous near-surface materials by joint 2D inversion of dc resistivity and seismic data. Geophysical Research Letters 30 (13), 1–4. Grandjean, G., 2006. Imaging subsurface objects by seismic-wave tomography: numerical and experimental validation. Near Surface Geophysics 4, 279–287.

47

Grandjean, G., Leparoux, D., 2004. The potential of seismic methods for detecting cavities and buried objects: experimentation at a test site. Journal of Applied Geophysics 56, 93–106. Kaufmann, G., Romanov, D., Nielbock, R., 2011. Cave detection using multiple geophysical methods: Unicorn Cave, Harz Mountains, Germany. Geophysics 76 (3), B71–B77. Langer, R.E., 1933. An inversion problem in differential equation. Bulletin of the American Mathematical Society, Series 2 29, 814–820. Linde, N., Binley, A., Tryggvason, A., Pedersen, L., Revil, A., 2006. Improved hydrogeophysical characterization using joint inversion of cross-hole electric resistance and groundpenetrating radar traveltime data. Water Resources Research 42, W12404. Loke, M.H., 1994. The inversion of two-dimensional apparent resistivity data. Unpubl. Ph.D. thesis, Un. Of Birmingham (U.K.). Loke, M.H., Barker, R.D., 1996. Practical techniques for 3D resistivity surveys and data inversion. Geophysical Prospecting 44 (3), 499–523. Malagodi, S., Orlando, L., Piro, S., Rosso, F., 1996. Location of archaeological structures using GPR method: three-dimensional data acquisition and radar signal processing. Archaeological Prospection 3, 13–23. McDonald, Davies R., 2003. Detection of sinkholes using 2D electrical resistivity imaging. First Breack 21, 32–35. Mochales, T., Casas, A.M., Pueyo, E.L., Piryo, O., Roman, M.T., Pocovi, A., Soriano, M.A., Anson, D., 2008. Detection of under-ground cavities by combining gravity, magnetic and ground penetrating radar survey: a case study from Zaragoza area NE Spain. Environmental Geology 53, 1067–1077 http://dx.doi.org/10.1007/s00254-007-0733-7. Musil, M., Maurer, H.R., Green, A.G., 2003. Discrete tomography and joint inversion for loosely connected or unconnected physical properties: application to crosshole and georadar data sets. Geophysical Journal International 153 (2), 389–402. Orlando, L., Piro, S., Versino, L., 1986. Location of sub-surface geoelectric anomalies for archaeological work: a comparison between experimental arrays and interpretation using numerical methods. Geoexploration 24, 227–237. Piro, S., Gabrielli, R., 2009. Multimethodological approach to investigate chamber tombs in the Sabine Necropolis at Colle del Forno (CNR, Rome, Italy). Archaeological Prospection 16, 1–14. Piro, S., Tsourlos, P.I., Tsokas, G.N., 2001. Cavity detection employing advanced geophysical techniques: a case study. European Journal of Environmental and Engineering Geophysics 6, 3–31. Quarto, R., Sciavone, D., 1996. Detection of cavities by the self-potential method. First Breack 14 (11), 419–431. Rybakov, M., Goldshmidt, V., Fleischer, L., Rotstein, Y., 2001. Cave detection and 4-D monitoring: a microgravity case history near the Dead Sea. The Leading Edge (Society of Exploration Geophysicists) 20 (8), 896–900. Scott, J.B.T., Barker, R.D., Peacock, S., 2000. Combined seismic refraction and electrical imaging. Procs. 6th Meeting of the European Association for Environmental and Engineering Geophysics, 3–7 Sept. 1997, Bochum, Germany, EL05. Sharma, S.P., Verma, S.K., 2011. Solution of the inherent problem of the equivalence in direct current resistivity and electromagnetic methods through global optimization and joint inversion by successive refinement of model space. Geophysical Prospecting 59, 760–776. Sheehan, J.R., Doll, W.E., Mandell, W.E., 2005. An evaluation of methods and available software for seismic refraction tomography analysis. Journal of Environmental and Engineering Geophysics 10, 21–34. Van Schoor, M., 2002. Detection of sinkholes using 2D electric resistivity imaging. Journal of Applied Geophysics 50, 393–399. Vichabian, Y., Morgan, F.D., 2002. Self potential in cave detection. The Leading Edge 866–871. Vozoff, K., Jupp, D.L.B., 1975. Joint inversion of geophysical data. Geophysical Journal of the Royal Astronomical Society 42, 977–991.