**1. Introduction**

Modern archaeological surveys require non-invasive techniques that can be employed on wide areas in order to study the shallow subsoil, highlighting the presence of physical anomalies that could be related to buried archaeological remains. These techniques should be able to identify the physical properties of the soils by means of measuring systems positioned on a mobile vehicle in order to reduce the time of acquisition compared to other more traditional methods [1].

Electrical resistivity prospection specifically proved to be an e ffective method to characterize the subsoil and materials in various fields of study, such as applied geology [2–4], hydrogeology [5,6], engineering [7–9], environmental sciences [10,11], and agriculture [12–14]. This method is widely used for archaeological surveys: many past research has shown that the variations of electrical resistivity are often related to walls, floorings, paving and burials [15–22]. These structures are generally characterized by di fferent electrical properties compared to the geological layers that, over time, have buried them.

Electrical methods based on fixed arrays permit the creation of resistivity models of the medium through a set of electrodes arranged along linear profiles or along two-dimensional configurations [23]. In the case of two-dimensional Electrical Resistivity Tomographies (ERT) the surveys provide a vertical section of the medium. The distance between the electrodes is constant. Each measurement is related to a specific distance between current and potential electrodes. Loops of electrodes with various geometries are used to perform tri-dimensional electrical resistivity tomographies in order to obtain a 3D model of the investigated medium [24–34].

Different systems for measuring and mapping electrical resistivity have been developed in order to overcome the main limitations of the fixed-array systems [35]. In particular, multi-polar mobile systems allow the performance of continuous resistivity profiling without a fixed device with an effective reduction in the time of acquisition. The first solution proposed a quadrupolar mobile device, provided with wheels equipped with electrodes for measuring the electrical resistivity of the soil at a constant theoretical depth [36]. The depth of investigation increased by varying the space between the wheels using an extendable support. Different solutions based on multi-polar mobile systems were developed and described in Panissod et al. [37]. These authors proposed an instrument based on eight electrodes named MUlti-pole Continuous Electrical Profiling (MUCEP), positioned along a V-shaped array (also called "Vol-de-canards" array) and located on a mobile measurement system. This equipment was used to explore shallow layers up to 3 meters in depth. In the same context, the comparison among resistivity prospections, obtained by the use of electrodes arrays towed by a vehicle, showed both advantages and limitations based on Pole–Pole arrays with a single depth of investigation with regards to the multi-pole systems. Multi-polar configuration provides a three-dimensional model of the subsoil with a significant reduction in acquisition time [38]. Further studies [35] have improved these multi-pole systems in order to reduce the noise effects on the raw data. In particular, the authors propose the equipment known as Automated Resistivity Profiling (ARP). This acquisition instrument is composed of a transversal current dipole and of three transversal receiving dipoles, set at an increasing distance from the transmitting one. This electrode array is designed to obtain simultaneously measurements of apparent resistivity simultaneously at three different subsoil depths. Papadopoulos et al. [39,40] described a three-dimensional inversion of resistivity data acquired with the ARP system in order to obtain a reliable numerical model capable of describing the spatial distribution of the electrical resistivity in the subsoil of the investigated area. The inversion process described by the authors is based on the least squares method and establishes smoothness constraints in order to consider the instability of the model and the non-uniqueness of solutions [41].

Several studies described multi-pole acquisition systems with the application of precision agriculture and soil characterization [42–45], archaeological prospection [46], and non-destructive evaluation of pavement and concrete infrastructures [47]. The ARP system was used by some authors [48] to evaluate the reduction of soil characterization costs for viticulture: the authors showed how the ARP system represents an effective method to reduce the cost of soil analysis, collecting a minor number of soil samples for direct tests. Many results of archaeological surveys carried out with the ARP system are extensively described in [1].

Integrated prospection with complementary geophysical methods permits to better define the physical model of the investigated subsoil [49]. One of the most integrated approaches in archaeology is the combination of electrical and Ground Penetrating Radar (GPR) methods [19,50].

In this paper, we present the results of surveys carried out by means of the ARP method at the archaeological site of Mont'e Prama (Sardinia, Italy). Geophysical investigations were designed to support future excavation activities and to define the spatial extension of the site in the northern side of an area in which archaeological digs took place in the 1970s. The analysis of electrical resistivity anomalies, both due to the buried archaeological remains and the geomorphological and geo-pedological spatial variability, permitted the evaluation of the effectiveness of this method for the site of study. Resistivity data were compared with the results achieved by extensive high-resolution surveys performed with multi-channel GPR at the same site. The archaeological interpretation of the experimental data was supported by the integration of different complementary geophysical methods and GPS data.

#### **2. Materials and Methods**

#### *2.1. The Archaeological Site of Mont'e Prama*

The study area of *Mont'e Prama* is one of the most important archaeological areas of the western Mediterranean [51]. The site is located in the Sinis peninsula (western-central Sardinia, Italy) as shown in Figure 1. The first archaeological remains were accidentally discovered in 1974 during soils ploughing when a local farmer found a fragment of a Nuragic sculpture. In 1975, 1977 and 1979 various campaigns of archaeological excavation were performed and a necropolis with characteristic monumental sculptures from the early Iron Age was studied. Many baetyls, models of "nuraghe" and numerous fragments of sculptures (usually called "Giants") were collected [51,52]. The last archaeological trials of these campaigns did not find new evidence and excavations were suspended.

**Figure 1.** (**a**) Location of the Sinis Peninsula in Sardinia (Italy); (**b**) geographic localization of the archaeological site of Mont'e Prama (red point); (**c**) detailed aerial image of the site of study where the surveyed area is white-bordered and highlighted.

Recently, new archaeological and geophysical investigations were carried out from 2013 to 2015 in the framework of a joint project between the Universities of Cagliari and Sassari and the Archaeological Superintendency of Cagliari and Oristano, in which context other significant discoveries were made [53–56]. Overall, 28 statues representing figures of boxers, archers, warriors, measuring up to 2.20 m in height, 16 models of "nuraghe" and baetyls were pieced together from over 5000 fragments of limestone. Archaeologists identified more than 50 shaft graves. Figure 1a,b indicates the geographic location of the site of study and the area with main archaeological digs carried out in 2015 (Figure 1c). After this campaign, the full extent of the archaeological site appears significantly greater than that supposed by the first researchers in 1970s.

## *2.2. Geological Setting*

The geology of the Sinis peninsula, the wide region including the Mont'e Prama site, consists of both sedimentary and volcanic rocks [57]. In particular, the geographic sector closest to prospection area is mainly characterized by sedimentary lithologies. These formations range from messinian fenestrae-rich limestones and breccias of supratidal and intertidal environment with typical benthic fossil content to evaporitic fine-grained limestones, sub-litoral marls and bioclastic limestones (*Calcari Laminati del Sinis*, *Capo S. Marco* formation). Quaternary conglomerates, sands and mud deposits, organized in river terraces and alluvial cones, are recognized in the upper part of the sequence.

In the context of 2014 archaeological digs, the stratigraphic sequence in correspondence to the excavations was reconstructed [58]: it is characterized by the presence of geological units related to the sedimentary complex of the Sinis peninsula. In particular, the site presents a thin layer of silt-sandy soils (Pleistocene and Olocene) with a maximum thickness between 40 cm and 50 cm. The soil layer covers a yellowish thick carbonate crust of pedogenetic origin. The hard crust is probably derived from salt precipitation as a result of capillarity action in past climatic conditions and reaches 1 m of maximum thickness. This horizon is observed with good spatial continuity in all over Sinis peninsula and covers previous messinian sediments, mainly composed by carbonate sandstones and marlstones. The area of geophysical surveys is expected to present mean characteristics similar to these results with a certain spatial variability.

#### *2.3. Methods and Instrumental Layout*

## 2.3.1. ARP Survey

The electrical surveys were carried out using the multi-polar ARP03 device (ARP ®, Automated Resistivity Profiling, Géocarta SA, Paris, France), a measuring system made of four pairs of rolling electrodes arranged in a characteristic V-shaped configuration. Figure 2a shows a picture of the ARP system. The electrodes located on the first pair of wheels represent the current dipole, through which an electrical current of 10 mA is flowing. The other three couples of rolling electrodes are three potential dipoles designed to measure the electrical potential with increasing depth (Figure 2b). The longitudinal distance between the potential dipole electrodes and the current electrodes corresponds to 0.50 m, 1.0 m and 1.70 m, respectively. Moreover, the three potential dipoles constitute of pairs of wheels with axial distances of 0.5 m, 1.0 m and 2.0 m, respectively. Therefore, the resistivity measurements are related to three di fferent depths of investigation, conventionally indicated as 0.5 m, 1.0 m and 1.70 m.

The sensor allows the acquisition of multiple measurements of apparent resistivity through the three pairs of rolling electrodes with a sampling frequency of 22.7 Hz. A radar doppler system is used to measure the distance. Depending on the vehicle speed, raw data can be acquired up to a 10 cm sampling interval.

Moreover, the positioning system, equipped with a di fferential correction DGPS, enables georeferencing of the acquired data and showing raw resistivity data in real time.

The speed of data acquisition is mainly due to the morphological features of the investigated area (slope, irregular topography, shallow obstacles, plan shape and irregularities, etc.). Usually, the employment of a quad bike to tow the carriage allows a maximum speed of 15 km/h.

The acquisition on the ground is typically performed along parallel profiles with variable length. For high-resolution surveys, such as archaeological investigations, the acquisition pattern is composed of profiles positioned ata1m distance.

 **Figure 2.** (**a**) Image of the Automated Resistivity Profiling (ARP) system towed by a quad; (**b**) scheme of ARP instrumental layout.

**(b)** 

The ARP survey at the site was carried out on November 21st, 2014, on a surface of about 2.28 hectares (Figure 1c). Two series of measurements were done:


Data were collected in autumn and on a day characterized by optimal weather conditions in order to perform the measurements with humidity conditions of the outermost shallow layers able to guarantee low contact resistances and high-quality measurements.

The acquisition allowed the detection of 277,319 punctual measurements of apparent electrical resistivity along a 27.8 km long route and with a velocity of 7.53 km/h along profiles oriented in an easterly direction. Table 1 summarizes the main acquisition parameters.


**Table 1.** Summary of the main parameters of the ARP prospection.

The spatial distribution of the measurement points is shown in Figure 3a. From this representation, it is possible to observe two distinctive points densities between the two main datasets, partially due to the different acquisition velocities, which are mainly related to the availability of undisturbed wide paths for the towing vehicle but mostly due to the overlap configurations; in fact, the smaller area was surveyed along shorter and orthogonal profiles, even because of superficial archaeological features partially prohibiting vehicle motion, which produced a denser points distribution. While in the northern free area it was possible to perform the acquisition with optimal sampling parameters, for the southern most difficult area, the requirement of spatial coverage and reliability of the local dataset entailed a little higher redundancy of data (rough density of measurements about 29 point/m2).

**Figure 3.** (**a**) Planimetric spatial distribution of apparent resistivity data. (**b**) Elevation contour map of the digital elevation model (DTM) obtained by integrated acquisition of GPS data over surveyed areas (dashed bordered areas).

The GPS data continuously and automatically acquired during the measurements at 1 Hz sample rate, integrated with isolated points external to the ARP acquisition area, also allowed the creation of a detailed digital terrain model (DTM) of the surveyed area, which provides useful information for the analysis and interpretation of ARP data. Based on DTM results, the site presents a regular topography with a weak slope angle (mean dip angle 4.3◦) indicatively oriented from west (top) to east (bottom) for a difference in height of about 16.5 m (Figure 3b).

In order to delete both negative resistivity values and outliers, raw data were filtered using a bi-dimensional median filter operating on a mobile spatial window with a 2 meter radius centered on each measurement point. The filtering algorithm removes resistivity values that are higher than the reference threshold function of the median of values measured within an area centered on every measuring point. Once filtered, the experimental data are interpolated using 2D cubic spline functions over plan surface to permit to obtain a regular grid that can be used to support a successive processing phase or the interpretation of the geophysical anomalies.

#### 2.3.2. GPR Survey

Ground Penetrating Radar methods study both electric and magnetic properties of materials by means of the interaction of microwave signals with the investigated bodies [59]. Electromagnetic waves, transmitted using di fferent superficial probes, can be refracted, reflected and di ffracted at the electromagnetic discontinuities and are absorbed passing through the medium. These e ffects are a function of material properties (e.g., electrical permittivity, magnetic permeability) and signal features (e.g., wavelength, waveform, etc.). These techniques have been widely used in recent decades in engineering, geology and archaeology [60].

In order to validate and interpret ARP data, a comparison with some GPR profiles is proposed. The Ground Penetrating Radar survey was performed in the context of the same research project on a wide area of study, even including the sector investigated by means of the ARP method [31–34]. The GPR prospection was done using a 15-channels (16 antennas) GPR system (STREAM-X, by IDS), with a 200 MHz central frequency that allows the acquisition of synchronous recordings of 15 parallel radargrams with a fast acquisition speed (about 10–15 km/h). Each radar profile is spaced 0.12 m from the others. Therefore, by means of this instrumental setup, it was possible to collect a big dataset characterized by a high spatial density (horizontal: 12 cm perpendicularly to the moving direction and 9 cm along it: vertical: less than 1 cm). A di fferential GPS antenna was used for positioning the measurements with a horizontal shift less than 5 cm. This equipment with respect to standard single/dual channels GPR systems includes di fferent advantages such as time e fficiency, fixed distance between the fifteen profiles (12 cm), absolute parallelism between all synchronous profiles and 3D highlighting of buried targets.

GPR data were processed using Reflexw software (by Sandmeier ®). Processing flow includes background removal, de-wow filtering, gain, bandpass filtering (100–300 MHz) and time–depth conversion.
