Preprint
Article

Application of Experimental Configurations of Seismic and Electric Tomographic Techniques to the Investigation of Complex Geological Structures

Altmetrics

Downloads

244

Views

88

Comments

0

A peer-reviewed article of this preprint also exists.

Submitted:

26 July 2024

Posted:

29 July 2024

You are already at the latest version

Alerts
Abstract
The main purpose of this study is the determination of the optimal data acquisition and processing parameters of electric and seismic tomographic techniques, for the investigation of complex geological environments. Two different study areas, in central-east Peloponnese and SW Attica were selected, where a detailed geological mapping and surface geophysical survey were carried out. The applied geophysical survey included the application of electrical resistivity tomography (ERT) and seismic refraction tomography (SRT). The geoelectrical measurements were acquired with different arrays and electrode configurations. Moreover, various types of seismic sources were used at seventeen shot locations along the seismic arrays. For the processing of geoelectrical data, clustered datasets were created increasing the depth of investigation and the discriminatory capability. The seismic data processing included: (a) creation of synthetic models and seismic records to determine the effectiveness and capabilities of the technique, (b) spectral analysis of the seismic records to determine the optimal seismic source type and (c) inversion of the field data to create representative subsurface velocity models. The results of the two techniques successfully delineated the complex subsurface structure that characterizes these two geological environments. The application of the ERT combined with the SRT, are the two dominant, high-resolution techniques for the elucidation of complex subsurface structures.
Keywords: 
Subject: Environmental and Earth Sciences  -   Geophysics and Geology

1. Introduction

Subsurface investigation of complex geological environments had always been a major challenge for near-surface applied geophysics. A complex geological environment refers to an area where various geological factors interact with each other in complicated processes, resulting in a dynamic and rapidly changing environment. These environments usually consist by a variety of geological formations, usually succeeding each other over a relatively short distance. They are characterized by the presence of extensional and compressional tectonic structures (faults, folds, overthrusts, detachments), due to the multiple deformational phases that have undergone. In addition to the tectonic causes, the exposure of the geological formations to exogenous and endogenous processes, such as erosion, weathering, faulting, magmatism and hydrothermal fluid circulation, further contributes to the complexity of the geological environment. All the above processes pose significant challenges in terms of investigating and interpreting the factors that led to the creation of such environments.
One dimensional (1D) geophysical techniques, such as transient electromagnetics (TEM), borehole logging and vertical electrical sounding (VES), as well as geotechnical investigation techniques like boreholes, are incapable of delineating adequately the subsurface structure of such environments, due to the presence of strong lateral variations [1,2]. Moreover, the high-resolution seismic profiling technique often gets affected, regarding the data quality, by diffraction and scattering phenomena, multiple reflections and unsuccessful utilization of migration techniques, due to the strong near-surface velocity contrasts and lateral velocity variations in fault zones [3,4].
Among several geophysical methods, electrical resistivity tomography (ERT) is a high-resolution and effective technique, capable of producing 2D or/and 3D subsurface resistivity distributions [5,6,7], which can typically be corelated with spatial variations in the lithological composition and the presence of tectonic structures in the investigated subsurface section. It has been widely applied to gain information about the location and structural characteristics of near-surface fault zones, their effective depth and estimate the width of the affected zone [8,9,10,11]. Furthermore, the ERT has been applied in complex geological environments for the delineation of conductive [12,13] and weathered zones [14], as well as for geotechnical purposes regarding the investigation of subsurface air-filled voids, fractures and bedrock’s depth [15,16,17].
The seismic refraction method is a widely applied geophysical method for near-surface investigations. Conventional seismic refraction processing techniques fail to reveal the true subsurface structure in complex environments, since they utilize oversimplified geometry, by dividing the substratum into discrete layers of constant velocity. However, recent advances in seismic data acquisition and complex inversion algorithms have established the seismic refraction tomography (SRT) technique as a prime geophysical tool for subsurface investigation in such environments [17,18]. These algorithms allow for the inversion of the P-wave first arrivals, resolving both vertical and lateral velocity variations in the subsurface, due to the presence of localized anomalies such as fault zones, karstic voids, boulders and lateral lithological variations [19,20]. Several authors have recommended the application of seismic refraction tomography (SRT) technique, for subsurface investigation in complex geological environments, as a way to overcome the limitations of the seismic profiling techniques [21,22,23].
The combined application and evaluation of the ERT and SRT techniques is proved to provide a more comprehensive picture of the subsurface structure, particularly in complex geological environments, since they are sensitive to different physical properties variations [6,20,24]. This combination enhances the interpretation by allowing the cross-validation of the results and reduction of the ambiguities that might arise when using a single method, thereby improving the reliability and accuracy of the geological interpretation. It benefits also the identification and characterization of a wide range of subsurface structures, including lithological variations and structural discontinuities.
The main purpose of the present study was the determination of the optimal acquisition and processing parameters of ERT and SRT geophysical data, for the investigation of complex geological environments. To this end, two different study areas characterized by complex geological structure were chosen for investigation. The first site is located in the wider area of Ano Doliana Arcadia (Figure 1), where, according to the detailed geological mapping carried out at the survey site, it is characterized by a variety of formations with different lithological composition and intense tectonic activity, resulting in a particularly complex geological environment. The second site is located in the area of Plaka (Figure 3) in SE Attica, where an Upper Miocene granodiorite intrusion and the associated contact metamorphism and deformation of the surrounding formations, has further exacerbated the complexity of the local geological structure.

2. Geological Setting

2.1. The 1st Study Area Kleisoura Valley, Ano Doliana

The first study area of Kleisoura Valley is primarily composed of Alpine formations that belong to three distinct geotectonic units of the Outer Hellenides. These units have been tectonically juxtaposed during the Alpine orogenetic cycle, involving processes associated with the initial stage of subduction and the late-orogenetic extensional stage. For the needs of the present study, a detailed geological mapping at a scale of 1:7500 was carried out, at the area where the near-surface geophysical survey was to be applied. In this way, we managed to identify the geological environment and surface structure of the study area and collect all the necessary information to evaluate the geophysical results. In Figure 1, the detailed geological map constructed for the study area is presented, where the location to which the geophysical survey was applied is noted by the red mark.
The study area includes the following Alpine units from bottom to top: the high-pressure metamorphosed Phyllite-Quartzite Unit, the very-low grade to non-metamorphosed neritic Tripolitza Unit and the uppermost pelagic Pindos Unit.
The lower unit is the Phyllite-Quartzite series of Permo-Triassic age [25], which mainly consists of schists alternating with quartzites, metapelites and phyllites (Figure 2f). This series is characterized by abrupt changes in lithological facies both laterally and vertically. The main petrographic types in this series include greenschists with residual glaucophane and paragonite; schists with chloritoid and glaucophane; chloritic-paragonitic-muscovitic schists with garnet; mica-quartzites with hematite, ilmenite and tourmaline; and muscovitic-quartzose phyllites with chloritoid. According to the existing mineralogical paragenesis (presence or absence of Na-amphibole) and due to the wide variety of petrographic types, the Phyllite-Quartzite series reflect a broad range of metamorphic conditions (T=280-480C°/P=2-10kb), corresponding to the blueschist and subsequent greenschist metamorphism facies [26].
The Phyllite-Quartzite series are tectonically overlain by Tripolitza Unit, which is the lowest non-metamorphic unit observed in the study area, mainly comprising two formations: a) an Oligocene flysch (Figure 2b), which is characterized by alternating fine-grained thin-bedded brownish sandstones, thin interbeds of pelites and sparse limestones lenses or interbeds. b) massive to thick-bedded and occasionally medium-bedded, grey-black, often bituminous and intensely karstic limestones and dolomites of Paleocene-late Eocene age (Figure 2a,d). In many cases, the transition from the Eocene limestones to the flysch is either gradual, marked by a few meters thick light-colored marly limestones (Figure 2a,c), or abrupt through high-angle faults (Figure 2d). It should be noted that several isolated Tripolitza limestones were observed, west of the provincial road, within the flysch of the same unit (Figure 2a). In these locations, the limestones appear to be folded and in most cases are bounded at their contact with the flysch by faults. Moreover, along the provincial road, an abrupt transition from the limestones to the flysch of Tripolitza Unit was observed, characterized by the presence of thick breccia (non-cohesive to cohesive cataclasites) at the structural top of the limestones (Figure 2e).
Figure 1. Detailed geological map of the 1st study area (Kleisoura Valley, Ano Doliana).
Figure 1. Detailed geological map of the 1st study area (Kleisoura Valley, Ano Doliana).
Preprints 113417 g001
The structurally highest unit of the region is the Pindos Unit (or Pindos nappe). The Pindos nappe primarily occupies the highest topographic regions, both west and east, along the central axis of the study area. In the study area it includes formations of the upper part of the stratigraphic column (e.g., the “Arcadian nappe”) [26], primarily pelagic thin-bedded biomicritic Cretaceous limestones (Figure 2a), alternating with thin chert layers and nodules and thin red to yellow pelitic layers. For the most part, the Pindos limestones are seen overlying the Tripolitza flysch, but at the southern part, they overlie the Phyllite-Quartzite series.
Figure 2. (a) General view of the west flank of the Kleisoura Valley where the “Arcadian nappe” (Pindos limestones- Pl) overlies the Tripolitza Unit; (b) Close-up view of Tripolitza flysch (fl); (c-d) The transition from the Tripolitza limestones (Tl) to the flysch (fl) through transitional marly limestone (ml) beds (c), or high-angle fault (d); (e) Thick tectonic breccia located at the structurally highest levels of the Tripolitza limestones; (f) Phyllite-Quartzite series: boudinaged light colored quartzitic layers (Q) surrounded by gray-green phyllites (Ph).
Figure 2. (a) General view of the west flank of the Kleisoura Valley where the “Arcadian nappe” (Pindos limestones- Pl) overlies the Tripolitza Unit; (b) Close-up view of Tripolitza flysch (fl); (c-d) The transition from the Tripolitza limestones (Tl) to the flysch (fl) through transitional marly limestone (ml) beds (c), or high-angle fault (d); (e) Thick tectonic breccia located at the structurally highest levels of the Tripolitza limestones; (f) Phyllite-Quartzite series: boudinaged light colored quartzitic layers (Q) surrounded by gray-green phyllites (Ph).
Preprints 113417 g002
The post-alpine sediments are seen north of the study area and include polymict Pliocene conglomerates featuring clasts originating from all the units (including metamorphic rocks of the Phyllite-Quartzite series). Moreover, Quaternary river deposits are present, mainly found along the Kleisoura Valley stream and alluvial deposits of considerable thickness resulting from in situ weathering.
The thrusting of the Pindos Unit on the Tripolitza Unit took place during the Middle Eocene to Oligocene, when the Tripolitza Unit was accreted in the overriding plate. Subsequently, this was followed by a synorogenic exhumation stage in the subduction channel and later by an extensional late-orogenic stage. Eventually, the metamorphosed Phyllite-Quartzite series were exhumed to upper crustal levels through a crustal-scale detachment system.
The detachment is exposed at the northern part of Kleisoura Valley, where the Tripolitza flysch and limestone are juxtaposed through a low-angle normal fault (dip ~10°-20°), generally diping towards the south. Conversely, in the southern part, the detachment fault presents a similar low-angle geometry (10°-20°) but dips generally towards the north. This change in dip direction is attributed to later neotectonic high-angle normal faults that crosscut all the previous compressional and extensional structures.

2.2. The 2nd Study Area “Plaka”

In Figure 3, the detailed geological map, modified after [27], of the wider study area is presented, where the location to which the geophysical survey was applied is noted by the red mark.
Lavreotiki (south Attica) is located in the northwestern part of the Attic-Cycladic crystalline complex and its structure is primarily controlled by a significant detachment fault with a top-to-SW sense of shear associated with the West Cycladic Detachment System [28,29]. The footwall of this detachment system includes formations of the Kamariza Unit, while the hanging wall is characterized by a complex structure that encompasses the metamorphic Lavrion Unit, the non-metamorphic Berzekos Unit and scattered occurrences of Neogene sediments. Moreover, the metamorphic units have been intruded by late Miocene magmatic bodies.
The Kamariza Unit is the structurally lowermost unit of Attica comprising three formations: the generally coarse-grained white Lower Kamariza Marble, the lithologicaly heterogeneous brownish calcareous Kamariza Schists (Figure 4g) and the white to blue-gray, ultra-mylonitic Upper Kamariza Marble (Figure 4f), with its upper boundary identified as the surface of the detachment fault (Figure. 4a) [28,30,31]. Due to contact metamorphism, the Kamariza schists have become hornfels within a few hundred meters around the primary granodioritic intrusion [32].
The Lavrion Unit of Mesozoic age is tectonically overlying the Kamariza Unit and constitutes the intermediate tectonic unit. It comprises from bottom to top the Triassic-Jurassic Pounta Marble (Figure 4e), the Lavrion schists and the Mavrovouni Marble [28]. The Lavrion Unit is characterized by Eocene age blueschist-facies metamorphism, followed by greenschist metamorphism during the Oligocene.
Finally, the uppermost unit is the the non-metamorphic Berzekos Unit, which consists of radiolarites, sandstones and cherts and is rich in late Jurassic-early Cretaceous silicate fossils. On top of these formations, the late Cretaceous rudist-bearing limestones have been deposited unconformably. In most areas, these limestones have been infiltrated by hydrothermal fluids and have undergone ankeritization [32,33].
Figure 3. Detailed geological map of the 2nd study area (Plaka). Modified after [27].
Figure 3. Detailed geological map of the 2nd study area (Plaka). Modified after [27].
Preprints 113417 g003
The Neogene continental deposits primarily consist of marls, siltstones and conglomerates, with clasts originating from the surrounding non-metamorphic units [32,33]. However, these deposits are not observed in the study area.
The magmatic intrusions include the main Plaka granodiorite with limited surface exposure, intruded within the Kamariza Unit. This intrusion is associated with contact metamorphism and sulfide mineralization [32]. Additionally, numerous small-scale intrusions with compositions ranging from dacitic to granodioritic have been observed in the area [32,33,34,35]. The main granodiorite intrusion locally exhibits intense weathering, intensified by hydrothermal fluid circulation (Figure 4d), resulting in zones of varying mechanical and geophysical properties [27].
Figure 4. (a) Panoramic view of the Plaka study area where the geophysical survey was conducted. The trace of the detachment fault that juxtaposes the overlying Pounta Marble (PM) against the underlying Kamariza Schists (KS) is marked with the yellow line; (b) Closer view of the study area, where the mining debris (Md) is illustrated; (c) Contact metamorphosed marble lenses (Ml) in the Kamariza Schists (KM); (d) Strongly weathered granodiorite; (e) Pounta Marble (PM); (f) Ultra-mylonitic Upper Kamariza Marble (UM); (g) Kamariza Schists-hornfels (KS).
Figure 4. (a) Panoramic view of the Plaka study area where the geophysical survey was conducted. The trace of the detachment fault that juxtaposes the overlying Pounta Marble (PM) against the underlying Kamariza Schists (KS) is marked with the yellow line; (b) Closer view of the study area, where the mining debris (Md) is illustrated; (c) Contact metamorphosed marble lenses (Ml) in the Kamariza Schists (KM); (d) Strongly weathered granodiorite; (e) Pounta Marble (PM); (f) Ultra-mylonitic Upper Kamariza Marble (UM); (g) Kamariza Schists-hornfels (KS).
Preprints 113417 g004

3. Methodology

To investigate the geologically complex subsurface structure of the two study areas, we employed a multidisciplinary, near-surface geophysical approach, utilizing the electrical resistivity tomography (ERT), seismic refraction tomography (SRT) and ground penetrating radar (GPR) techniques. More specifically, the ERT and the SRT techniques were applied to both study areas (Figure 5), while the GPR technique was used exclusively in the study area of Plaka. Furthermore, all the geophysical measurements were topographically corrected using a GNSS system.
In Figure 6 and Figure 7, the acquisition layout of the geophysical techniques applied in the 1st and 2nd study areas is presented correspondingly.

3.1. Electrical Resistivity Tomography

3.1.1. Data Acquisition

Across the first study area of Kleisoura Valley, Ano Doliana, two ERT profiles of 470m (ERT 1 & ERT 2) and one of 235m (ERT 3) were conducted, using 48 electrodes with a spacing of 10m and 5m correspondingly (Figure 6). The two long profiles were executed along the same direction, with the ERT 2 profile shifted 40m towards the north, resulting in a combined profile of totally 510m length. The short profile was designed to coincide with the center of the ERT 2 section. The application of the ERT 3 profile, with the smaller electrode spacing (5m), aimed to increase the resolution of the technique in this intermediate part of the profile. This is expected to improve the delineation of the abrupt surface transition from the flysch to the limestones of Tripolitza Unit, verified by the presence of the tectonic breccia in this specific area.
Regarding the second study area of Plaka, two ERT sections of 470m and 475m were acquired (Figure 7), with 48 electrodes of 10m and 5m spacing respectively, resulting in a joined ERT profile. The roll along technique was implicated in two consecutive phases, to cover the same profile length for the 5m spacing section.
All the ERT measurements were acquired using the Wenner, Wenner-Schlumberger and Dipole-Dipole electrode configurations, due to their different sensitivity regarding the way of the subsurface distribution of the electrical resistivity. Using the Electre Pro software, the sequences for the different electrode configurations were created for each ERT profile, according to its geometrical characteristics. Subsequently, the sequences were imported to the Syscal Pro switch 48 (Iris Instruments) resistivity meter, to begin the automated data acquisition process. Moreover, a real time generation of the apparent resistivity pseudosection was established using the FieldView software, to directly inspect the measurements sequence in the field.

3.1.2. Processing

The geoelectrical data were processed using the Res2DInv software by Geotomo [5]. Before the main processing, quality control of the data was performed, which included the extermination of noisy and inconsistent points (bad points). For the solution of the inverse problem, the smoothness-constrained least-squares method was applied, which is a modification of the Gauss-Newton least-squares equation [36], while the finite element method [37] was implemented for forward model calculation.
Depending on how the electrical resistivity of the subsurface varies, two different optimization methods are provided by the software. In the case where the electrical resistivity varies in a smooth and gradual manner, the l2 norm smoothness-constrained optimization method is proposed [38]. Conversely, when the electrical resistivity distribution is characterized by sharp and strong lateral variations, the l1 norm smoothness-constrained optimization method is recommended [39,40,41].
The subsurface structure of the 1st study area is expected to be characterized by steep transitions between geological formations, both laterally and vertically. Laterally, due to the abrupt transition from the flysch to the limestone formations of Tripolitza Unit and vertically, due to the expected existence of the Phyllite-Quartzite Unit, underlying the flysch and limestone, through the detachment fault. For the above reasons, the l1 norm smoothness-constrained optimization method (robust data constraint) was selected for the inversion of the geoelectrical data of the 1st study area. On the other hand, a more gradual subsurface resistivity distribution was expected due to the irregular shape of the granodiorite intrusion and the contact metamorphism of the Kamariza Schists, in the 2nd study area of Plaka. Therefore, the l2 norm smoothness-constrained optimization method (standard data constraint) was considered the most adequate for the calculation of the model parameters in this area.
Every ERT profile was processed separately for each of the three different electrode configurations applied. With this procedure, a total of nine and six different subsurface resistivity distribution models were generated for the 1st and 2nd study areas respectively. The geoelectrical data obtained with the same electrode configuration were integrated into a unified clustered dataset, resulting in three ERT profiles (one for each electrode configuration) for each study area. The most representative ERT section for each study area was selected taking into consideration the data quality and RMS or absolute percentage error between the observed and calculated resistivities, as well as the overall characteristics of the subsurface resistivity distribution.

3.2. Seismic Refraction Tomography

3.2.1. Data Acquisition

In both study areas, the recording of the P-wave first arrival times was carried out along a 235m seismic line, using 48 geophones of 10Hz natural frequency, with 5m spacing. The Geometrics StrataView seismograph was used to record arrival times, to which the seismic receivers were connected via two multicore cables. The total recording time was 512ms and a sampling window of 0.250ms was used for signal digitization. Seismic waves were recorded from seventeen shot-points (Figure 6 and Figure 7), including four outshots on both ends of the seismic line, and nine shot-points in between the active spread with a 30m shot interval.
The accessibility of the 1st study area (Kleisoura Valley, Ano Doliana) and the fact that it was located outside a residential area, allowed us to test several types of seismic sources, including three impact seismic sources and two explosive seismic sources. The impact seismic sources were two seismic sledgehammers of 3.5kg and 6.5kg and an accelerated weight drop (Geodevice AWD-33PS model) with 20kg hammer weight. Subsequently, the electric seismic detonator and buffalo gun were used as explosive seismic sources.
In the 2nd study area of Plaka, only the 6.5kg seismic sledgehammer and the 20kg AWD-33PS seismic source were feasible to be used, with limitations. Since this area was located nearby a residential area, the use of explosive sources was prohibited.

3.2.2. Processing

Seismic data were processed using different modules of the Seismic Pro software (Geogiga tech. corp.), in three (3) main processing steps:
  • Creation of synthetic seismic data (2D forward modeling), to determine the extent to which the SRT technique can provide a velocity model that best represents the complex subsurface structure characterizing the two study areas. Furthermore, by comparing the synthetic seismic records with the field seismic records, it was possible to attribute certain characteristics observed in the records to specific subsurface structures.
  • Analysis and comparison of the seismic records acquired from the five different seismic sources to determine the type of seismic source that performed best in the specific geo-environment, providing the seismic records characterized by the highest S/N ratio.
  • First brakes picking and inversion of the seismic data, to calculate the subsurface velocity models in the two study areas.
The creation of the synthetic seismic records was conducted using the Modeling 2D module, which incorporates the shortest path method [42] for forward modeling calculation. Seismic rays trajectories and corresponding travel times were calculated based on a predefined velocity model, that was created taking into account the surface geological observations and the results of the ERT technique (Figure 8), given that it was the only available information regarding the subsurface structure for both study areas. In general, resistive formations were modeled as compacted formations, characterized by high seismic velocity values, while conductive formations were modeled as less compacted formations or unconsolidated sediments, with lower seismic velocity values.
An exception was the Kamariza Schists in the 2nd study area, which was characterized by low resistivity and high-velocity values. Seventeen synthetic seismic records were created for both study areas, simulating the data acquisition layout incorporated in the field. Subsequently, the first breaks were manually picked on the synthetic datasets and the synthetic traveltimes were then inverted using the smoothing-constrained regularized inversion approach of the DW tomo module. Afterwards, synthetic traveltimes obtained from the synthetic seismic records were superimposed on the field seismic records of the same shot location, to determine whether the synthetic traveltimes were consistent with the actual traveltimes. This processing sequence was repeated many times by altering the initial velocity model parameters and inversion settings, until a representative velocity model was obtained with minor deviations between the modeled and observed first break traveltimes.
The selection of the appropriate type of seismic source is crucial, as the characteristics of each source affect the quality of the data, the maximum depth of investigation and the maximum resolution that can be achieved [43,44]. According to several researchers [45,46], the selection of the appropriate type of seismic source should be based on the following criteria: Maximum depth of investigation, appropriate frequency content, S/N characteristics, investigation field limitations, availability and operation cost. In general, the most critical factors for selecting the proper seismic source type for seismic data acquisition are the energy and frequency content of the generated signal [44,47,48,49], as well as the repeatability of the source [50,51]. The signal characteristics emitted by a source are affected by many factors. The way the data are acquired (single or consecutive records), as well as the recording and processing parameters, are some of the factors that affect the signal’s characteristics. In addition, the surface and subsurface conditions in the area where the survey is carried out greatly influence the amplitude and frequency content of all types of recorded waves.
Seismic sources signal analysis and comparison was carried out only for the 1st study area, given that this area allowed us to utilize all five seismic sources for the data acquisition. Comparisons were made between the seismic records acquired at each shot-point, from the five seismic sources used. Both the recording and processing parameters were kept constant at each location. The amplitude spectral content diagrams were constructed for all seismic records, by applying the fast Fourier transform. For each seismic record, a specific time window was selected for the analysis to be performed, to analyze the part of the record characterized by the highest quality signal. The amplitude spectra of the seismic records are represented relative to the reference level (0 db), representing the maximum recorded amplitude. Negative db values indicate whether the amplitude of different frequencies is smaller than the reference level. In addition, the amplitude spectra diagrams of the different seismic sources, constructed at each shot-point, are presented normalized to the one characterized by the maximum energy.
The processing of the field seismic records was conducted using the DW tomo module. The first processing step was the creation of a group file list, which included the seventeen seismic records characterized by the best S/N ratio, with the clearest first breaks, thus facilitating the management and processing of the data. In both study areas, all the seismic records acquired with the AWD-33 PS seismic source were used, as well as the seismic records obtained with the 6.5kg sledgehammer, but only in the case of the 2nd study area. At this step, the geometrical characteristics and the topography of the seismic developments were also defined. Moreover, low-pass Butterworth filters with a cut-off frequency of 200Hz and a slope of 12db/oct, as well as adjustment of the signal gain were applied in order to remove the high-frequency noise and highlight the first breaks in the seismic records respectively.
Manual first break picking was efficiently accomplished in fourteen (14) and in ten (10) out of the seventeen (17) acquired seismic records for the 1st and 2nd study areas correspondingly. The high vegetation, the possible existence of underground mining galleries, the anthropogenic deposits (mining debris) of the old mines excavations located at the central part of the seismic section and the difficulty of carrying and utilizing the AWD-33PS seismic source in the outshot positions, resulted in the unsuccessful recording of the first arrivals in the 2nd study area and thus a great loss of information. Afterwards, a primary processing was carried out by applying the intercept-time method, in order to calculate a gradient initial velocity model. In the case of the 1st study area, P-wave velocity was calculated equal to 1070 m/s near the surface, gradually increasing to 2435 m/s until the depth of 17.3m and further increasing to 4240 m/s until the depth of 70m. Conversely, the initial velocity model calculated for the 2nd study area was characterized by 1360 m/s surface velocity, gradually increasing to 2550 m/s until the depth of 22m and further increasing to 5840 m/s until the depth of 60m.
For the calculation of the subsurface P-wave velocity models, the smoothing-constrained regularized inversion approach was applied to invert the seismic data. For model smoothing, 15m and 6.25m horizontal and vertical smoothing lengths were set respectively, while for the inversion limitation, a maximum change of 60% in average velocity per iteration was allowed. A total of 30 iterations were defined for the inversion process, with the possibility to select the desired result from any iteration cycle. For forward modeling solution, the shortest path method [42] was incorporated, with horizontal and vertical model cell dimensions equal to 2.5 x 1.25m, corresponding to the ½ and ¼ of the geophone spacing respectively.

3.3. Ground Penetrating Radar

The broader area of Plaka was submitted in extensive mining campaigns from the Classical period until the late 20th century. Thus, in the entire region, including the 2nd study area, several mining galleries are present. Due to the possible existence of such galleries throughout the selected geophysical line, the GPR method was implemented in an effort to locate any possible underground mining corridors [52]. The GPR method is a fast applied, low cost and high-resolution geophysical method, based on the recording of the electromagnetic (E/M) pulse reflections. The resolution and propagation depth of the E/M pulse is constrained by the central antenna frequency and the conductivity of the substratum. The reflectivity intensity is affected by the dielectric constant differences of the subsurface mediums [53].

3.3.1. Data Acquisition

GPR data were acquired along a 177m line, starting from the 14th geophone up to the end of the seismic line (Figure 7). A bistatic 100 MHz shielded antenna (Noggin 100-Sensors & Software) was used to achieve a satisfactory investigation depth (~15m), while the pulse transmission step (5cm trace interval) was controlled by an odometer wheel. In each step, several E/M pulses were stacked to improve the S/N ratio, by incorporating a feature of the GPR system, which automatically controls the number of stacks depending on the conductivity of the medium.

3.3.2. Processing

GPR data were processed using the Ekko_Project 6 software by applying a standard filtering procedure including Dewow, background removal, bandpass filter, SEC2 gain and f-k migration using the hyperbola fitting method.

4. Results and Discussion

4.1. Forward Modeling

In Figure 9, the seismic tomogram calculated after the inversion of the synthetic traveltimes is presented in respect to the ERT results, based on which the initial subsurface velocity model was created for the 1st study area (Kleisoura Valley, Ano Doliana). It is evident that all the structures simulated by the construction of the synthetic subsurface velocity model are highlighted. More specifically, both the high resistivity regions, simulated as high-velocity formations are adumbrated by the SRT technique, in 165-215m distance and at the second half of the profile. Moreover, the subvertical electrical discontinuity simulated as a strong lateral velocity variation, is delineated at 300m distance, as well as the low-velocity region in the middle part of the profile, with the characteristic U-shaped form of the iso-velocity curves. Therefore, we conclude that the SRT technique was able to render, to a fairly satisfactory degree, the complex subsurface structure considered for this case.
In Figure 10, we present three (3) indicative examples of overlaying the synthetic first arrival times onto the field seismic records, acquired with the AWD-33PS seismic source at three (3) different shot locations. The synthetic first arrivals are marked by the red symbols, while the position of the true first arrivals or their divergence from the synthetic ones is indicated with the blue arrows.
Regarding the outshot normal 1 shot-point (Figure 10a), located at 120m offset, a small divergence of 15ms on average is observed between the synthetic and true first arrivals. However, the general distribution of the synthetic first arrivals is consistent with that of the true ones. Furthermore, the simulation of the low-velocity region in between the high-velocity formations proves to be reasonable, as indicated by the increase in the slope of both the synthetic and true refracted arrival times, at 257.5-297.5m distance.
Similarly, a very small divergence of less than 5ms is observed in the case of the tomo shot 1 shot-point (Figure 10b). In general, the position of the synthetic first arrivals is identical to the position of the true arrivals in the overall record. Small divergencies are observed to a few seismic traces after the 300m distance, where the subvertical electrical discontinuity was investigated.
Finally, in the case where the shot-point is located at distances greater than 300m (Figure 10c), large divergencies are observed between the synthetic and true traveltimes, on both the direct and refracted arrivals. The higher slope of the true direct arrivals indicates a lower seismic velocity at shallow depths, for the seismic formation located at distances greater than 300m. While in the two previous cases, the synthetic and true refracted arrivals were almost identical for the same traces, the divergence observed on the direct arrivals indicates a vertical gradient in the seismic velocity of this formation.
The large divergencies observed between the direct and refracted arrivals in this case, are attributed to the strong lateral and vertical heterogeneity that characterizes the subsurface structure of the study area, which was not accounted for the initial velocity model constructed for modeling. On the contrary, the initial velocity model was simulated by a subsurface structure characterized by the presence of discrete, homogeneous and uniform seismic formations.

4.2. Source Comparsion

In Figure 11, an indicative example of the seismic records and their corresponding amplitude spectra plots are presented, for the outshot normal 2 shot-point (82m offset), acquired with the AWD-33PS (AWD 20kg), 6.5kg sledgehammer (SH 6.5kg), seismic detonator (SD) and buffalo gun (BG) seismic sources.
Regarding the seismic records acquired with the explosive seismic sources, from the 13th seismic trace onwards, i.e. at distances greater than 142m from the shot-point, high frequency noise predominates. As for the impact seismic sources records, these were acquired by four and six vertical stackings for the AWD 20kg and SH 6.5kg sources respectively. In the seismic record of the SH 6.5kg source, there is a lot of high-frequency noise contamination, from the 19th trace to the end of the record, which makes difficult the identification of the P-wave arrival times. Significant improvement is observed in the AWD 20kg seismic record, where P-wave arrival times are distinct and noise levels are low, even in the most distant seismic traces.
For signal analysis, the time window of 28-250ms, from the 2nd to the 13th seismic trace (Figure 11 - grey shaded box) was chosen. Within the selected window, P-wave refracted arrivals are discerned, at 37ms and 60ms for the 2nd and 13th seismic traces respectively, while after 75ms surface waves dominate. All the amplitude spectra plots are presented normalized to the AWD 20kg source plot, which is characterized by the highest energy amplitudes. The maximum energy amplitude of the BG, SD and SH 6.5kg sources is marked by the blue, yellow and green horizontal lines correspondingly and it is concentrated in a frequency band near 60Hz. In the case of the AWD 20kg source, the corresponding amplitude is increased by 30, 24 and 9 dB in comparison with the BG, SD and SH 6.5kg sources respectively.
Most of the body and surface waves energy is concentrated within a frequency bandwidth of 10-150 Hz, for the selected time window. Furthermore, there is much less damping of the high-frequency waves amplitude in the case of the AWD 20kg source compared to the others. This variation is due to the higher energy released by the AWD 20kg source, along with the fact that the damping of the seismic energy is proportional to the frequency and propagation distance of the seismic waves.

4.3. Field Data Results of the 1st Study Area Kleisoura Valley, Ano Doliana

The results of all the geophysical techniques applied in the 1st study area, along with their representative geological profile, created after the combined evaluation of the results, are presented in Figure 12. The Wenner-Schlumberger ERT profile was considered the most representative of the subsurface structure. Due to the wide range of resistivity values, the logarithmic scale was chosen to visualize the data distribution. The SRT profile was finalized after the 22nd iteration of the inversion process, achieving a divergence of less than 3ms between the calculated and observed first arrivals.
Regarding the subsurface resistivity distribution (Figure 12a), the lower resistivity values (<50 Ohm.m) appear to be constrained near the surface, between 190-485m distance, down to a depth of 3-12m. Very characteristic is the presence of a sharp, subvertical electrical discontinuity, located at 300m distance and at depths greater than 10m. This discontinuity separates two areas characterized by a strong contrast between their resistivity values, emphasizing the lateral heterogeneity of this study area. At distances greater than 300m (north of the electrical discontinuity) to the end of the profile, underlying the surface conductive zone, a highly resistive formation is identified, with resistivity values greater than 5600 Ohm.m. From the beginning of the ERT profile, up to 300m distance (south of the electrical discontinuity) and underlying the surface conductive zone, a less resistive region is investigated, with resistivity values lower than 3500 Ohm.m. Within this region, a further distinction can be made concerning the electrical resistivity distribution. Between 135-210m distance and for an average depth of 15-45m, a resistive formation is identified, characterized by resistivity values ranging from 1300-3500 Ohm.m. In the remaining part of the profile, the subsurface model is characterized as relatively resistive, with resistivity values between 100-1000 Ohm.m.
In Figure 12b, the final P-wave tomographic model obtained after the 22nd iteration of the inversion process, with 3ms RMS error, is presented. The subsurface velocity distribution is characterized by strong lateral variations and the presence of regions with high vertical gradient. Low P-wave velocity values of 800-1300 m/s are concentrated at the near-surface part of the profile, at distances greater than 180m to the end of the profile. Between 220-345m distance, an arched shape of the iso-velocity curves is observed, resulting in a strong lateral variation in the 12-40m depth range. From the start of the profile to 180m distance higher P-wave velocity values are investigated (1700-3000 m/s), which at greater distances underly the abovementioned low-velocity values. In particular, from the start of the profile to 135m distance, a high vertical velocity gradient is observed, where the VP rapidly increases from 1700 m/s at the surface to 3000 m/s, at the depth of 20m. Between 135-220m distance, a high-velocity region is present until the average depth of 45m, with VP = 2400-3000 m/s. On the central part of the profile (265-335m) and in the depth range of 40-65m, a region characterized by VP = 1900-2400 m/s is present, with the characteristic U-shaped form of the iso-velocity curves, resulting in a lateral velocity variation in respect with the adjacent regions at this depth. At distances greater than 345m to the end of the profile and at a depth range of 10-65m and 5-20m correspondingly, a wide variety of velocities is identified (1700-3000 m/s), characterized by a lower vertical velocity gradient in conjunction with the first part of the profile. Finally, a high-velocity (VP > 3000 m/s) formation has been investigated throughout the whole profile length and at the greatest investigation depths. At the southern and northern parts of the profile, the top of this formation is present at 20m depth, while towards the central part of the section it is investigated in progressively greater depths, reaching the maximum value of 65m at 315m distance.
The final geological interpretation of the geophysical results is presented in Figure 12c, which was derived from the combined evaluation of the ERT and SRT techniques and taking into consideration the geological observations and data collected during the detailed geological mapping that was carried out in the area.
The location where the subvertical electrical discontinuity was investigated by the ERT technique is in agreement with the location where the abrupt lateral transmission from the flysch to the limestones of the Tripolitza Unit was observed. Moreover, the tectonic breccia at this location led to evaluating this discontinuity as a normal south-dipping fault, which brings these two geological formations into contact. Another indication supporting this interpretation is the strong lateral variation observed in the distribution of the seismic velocities, in a 50m width vertical zone, symmetrically to the fault’s location. This zone is characterized by VP values of 1200-2400 m/s, indicating the degradation of the mechanical properties of the subsurface formations in this region, due to the fault’s impact. Furthermore, according to the SRT profile, the fault appears to develop down to 65m depth, as at greater depths a high-velocity (VP > 3000 m/s) formation is observed, which does not appear to have been affected by the fault.
From 190m to 485m distance, the low resistivity (< 50 Ohm.m) and low-velocity (800-1100 m/s) geophysical formation is interpreted as the incohesive soil sediments of the alluvial deposits. These deposits are characterized by 12m maximum thickness and overlie the fault zone at 300m distance. The high conductivity of these deposits is due to the phreatic aquifer that develops within them, while the low-velocity values are attributed to the fact that consist primarily of loose, non-cohesive soil materials, products of the weathering and erosion of the adjacent geological formations.
At distances greater than 300m (north of the fault’s location), the high resistive formation (> 5600 Ohm.m), which is also characterized by a wide variety of seismic velocities, VP = 1300-3000 m/s, is geologically interpreted as the Tripolitza Unit Eocene limestones. This formation develops underlying the alluvial deposits, except for the last 25m of the profile where it is exposed on the surface. The maximum thickness of this formation is located at 300m distance where it reaches 55m, while towards the end of the profile it continuously decreases, reaching the minimum value of 20m. According to field observations, the limestones appear to be highly karstified at this location, which justifies both the high resistivity values and lower seismic velocities that were investigated. Another factor that has significantly downgraded the formation’s seismic velocity (1300-1900 m/s) is the fault’s activity, in a zone of 25m width. The relatively intact part of the formation is observed at distances greater than 25m from the fault’s location to the end of the profile and at depths greater than 30m, where it is characterized by VP = 2200-3000 m/s.
The deepest investigated seismic formation, characterized by high VP values of 3000-4400 m/s, is interpreted as the bedrock of the study area, consisting of the metamorphic lithologies of the highly heterogeneous Phyllite-Quartzite formation. Quartzite is characterized by higher VP values than phyllite, which justifies the lateral variations in the formation’s velocity. This formation was investigated near the surface at the beginning and the end of the profile, at depths greater than 20m, while towards the central part of the profile (295-355m) it is observed at depths greater than 65m. This formation is in contact with all the overlying formations through the detachment fault, identified in the field at about 40m northern to the end of the geophysical line.
The resistive formation (1300-3500 Ohm.m) investigated between 135-210m distance and at a 15-45m average depth range, characterized also by high VP values of 2400-3000 m/s, is interpreted as a local occurrence of Tripolitza Unit limestones. At this location, the limestones are characterized by higher seismic velocity and lower resistivity values, than those north of the fault. This indicates that they have been affected to a much lesser extent by the karstification phenomenon, probably because they are surrounded by impermeable formations, which prevented the development of the karstification phenomenon.
The remaining part of the profile, characterized by the presence of a relatively resistive formation (100-1000 Ohm.m), with seismic velocities varying from 1200-3000 m/s, is interpreted as the flysch of the Tripolitza Unit. At distances greater than 185m, the flysch is underlying the alluvial deposits, up to 300m distance, where it comes into contact with the Tripolitza Unit limestones through the fault. Generally, the low-velocity values occur near the surface but with a high vertical velocity gradient. The lower velocity (1200-2200 m/s) and higher resistivity values (500 – 1000 Ohm.m) are observed near the fault’s location. The flysch is also a highly heterogeneous formation, consisting of sandstone and clay alternations, which justifies the lateral variation of both its resistivity and velocity values.

4.4. Field Data Results of the 2nd Study Area “Plaka”

The results of the ERT and SRT techniques applied along the geophysical line of the 2nd study area, as well as their geological interpretation are presented in Figure 13. The ERT profile considered to best represent the subsurface structure was that one obtained by the Wenner-Schlumberger dataset processing and is illustrated using a logarithmic scale for its resistivity distribution visualization. The subsurface velocity model was calculated after the completion of the 20th iteration of the inversion process, with a deviation between the calculated and observed arrival times of less than 3.2ms.
In the NW part of the ERT profile (Figure 13a), a highly resistive (> 1000 Ohm.m) formation is observed, from the beginning to 140m distance. The maximum investigated depth of this formation is 27m, located at 62m distance, while it seems to decrease towards the SE gradually. From 90-185m distance and underlying the aforementioned resistive formation, a conductive zone with resistivity values of 20-100 Ohm.m is delineated. At distances between 140-185m, the zone develops close to the surface and it is characterized by an average thickness of 22m. A second conductive zone, with the same resistivity values, is identified at 190-240m distance and at depths greater than 10m. The second zone is characterized by an average width of 60m and extends up to 85m depth, corresponding to the maximum investigation depth of the technique. These two conductive zones are separated by a resistive formation, with resistivity values ranging from 350-2500 Ohm.m. At distances less than 185m, this formation underlies the first investigated conductive zone, while between 185-240m distance develops close to the ground surface. The high resistivity values (1200-2500 Ohm.m) of this formation are confined near the surface at 185-205m distance, while the rest is characterized by lower resistivity values (350-1000 Ohm.m). At distances greater than 240m to the end of the ERT profile, the subsurface resistivity values range from 45 to 800 Ohm.m. From the surface and up to an average depth of 35m, lower resistivity values (45-200 Ohm.m) dominate, while at greater depths, higher resistivity values (250-800 Ohm.m) are observed.
Regarding the subsurface velocity distribution of the SRT profile (Figure 13b), the lower seismic velocities (800-1600 m/s) are constrained in the near-surface part of the section, where a slight increase in the velocities is observed at the 2nd half of the profile. From 120-290m distance, this low-velocity zone has an average thickness of 10m, while at greater distances an increase of the zone’s thickness occurs, reaching the maximum value of 35m at 330m distance. Underlying, a second velocity zone is delineated, characterized by a wide range of seismic velocities (1800-4700 m/s) and thickness varying from 10-45m. The majority of this zone is observed with seismic velocity ranging between 1800-2700 m/s. However, there are some local regions within which the highest seismic velocities of this zone are concentrated. They are located between 140-185m and 240-270m distance, in a depth range of 10-30m and are characterized by high VP values, ranging from 3700 - 4700 m/s. Finally, at the deepest part of the profile, a third velocity zone characterized by the highest velocity values is investigated. In this zone VP values greater than 3000 m/s are observed, while seismic velocities of 5800 m/s are also present. This zone is located at depths greater than 30m and 58m, located at 155m and 320m distance correspondingly, while between 205-240m distance, an updoming of the iso-velocity curves in this zone is observed.
The SRT profile of the 2nd study area was significantly limited compared to the corresponding model of the 1st study area, due to the unsuccessful recording and processing of the outshot seismic records. For that reason, the geological interpretation (Figure 13c) was mainly based on the ERT results, apart from the joint section of the two geophysical techniques that was based on their combined evaluation.
The highly resistive formation (> 1000 Ohm.m) that develops close to the near-surface, at the NW part of the ERT profile, is evaluated as the Pounta Marble which is consistent with the geological mapping of the area. According to surface geological observations, these marbles appear highly karstified in places, justifying the increased resistivity values. The maximum investigated thickness of this formation is 27m, located at 62m distance, gradually decreasing towards the SE up to 140m distance.
The two conductive zones (20-100 Ohm.m) investigated by the ERT technique, as well as part of the relatively resistive zone (350-1000 Ohm.m), located at depths greater than 17m, between the two conductive zones, are interpreted as the Kamariza Schists, which has been contact metamorphosed into hornfels. In several locations near the study area, surface occurrences of massive hornfels lenses of former marble have been observed within the Kamariza Schists (Figure 4c). These lenses may account for the increased resistivity values located between the two conductive zones. The first 10-15m of the formation’s thickness are characterized by velocities of VP = 800-1800 m/s, while at greater depths, a significant increase in the velocity values is observed, ranging from 2800-4800 m/s. The highest seismic velocity values are concentrated in a zone between 205-240m distance, at depths greater than 35m. At distances less than 140m, the Kamariza Schists underlies the Pounta Marble. The contact between these two formations is the detachment fault, which dips to the NW and is clearly delineated from the ERT profile (Figure 13a). According to geological observations, the Kamariza Schists is highly heterogeneous. In places its mineralogical composition consists primarily of phyllosilicate minerals (mainly chlorite and white micas), with intense weathering and hydrothermal alteration and in other places the carbonate minerals predominate. This heterogeneity accounts for the lateral variations observed in both the resistivity and seismic velocity investigated values.
The high resistivity (1150-3500 Ohm.m) and low-velocity (800-1000 m/s) values investigated between 190-210m distance, from the surface to the maximum depth of 7m, are attributed to the presence of anthropogenic deposits, consisting by the excavation waste materials from the old mining activity. During the data acquisition campaign, these deposits were also observed in the field at this exact location (Figure 4b).
Considering the subsurface distribution of the resistivity and P-wave velocity values at the profile’s segment, located at distances greater than 170m and overlying the Kamariza Schist formation, two different regions can be distinguished. The first one is located from the surface to an average depth of 12m, between 185-240m distance, while from 240-430m distance until the average depth of 25m. Within this region, low resistivity (50-350 Ohm.m) and low-velocity (800-1500 m/s) formations seem to dominate. Geologically, it is interpreted as the Plaka granodiorite formation, which, according to our field observations, is strongly weathered and hydrothermally altered [27]. Some of the circular forms could represent woolsack weathering. The thickness of this weathering-affected layer appears to decrease towards the SE, e.g., away from the ravine that runs through the granodiorite. The second region is located underlying the first one, extending to the maximum investigation depth of the two geophysical techniques. Within this region, higher resistivity (400-770 Ohm.m) and higher velocity (2000-5500 m/s) values are observed. These increased values can be justified by the presence of the Plaka granodiorite, which in this case, has been affected to a much lesser extent by the weathering processes. In general, the P-wave velocities of this formation appear to increase with respect to depth, indicating a compact, fresh granodiorite existing at greater depths. According to the distribution of the resistivity values at distances greater than 430m, this formation appears to develop close to the surface, which is consistent with the geological mapping and the geological observations carried out in the field.
In Figure 14, the processing results of the GPR line for the 2nd study area are presented. The method highlights the possible existence of air voids (GPR signal velocity ~0.3 m/ns) marked as red dashed circles, at 62m, 84m and 144m distance at a depth of approximately 7m. These air voids are interpreted as possible old mining galleries in the area.

5. Conclusions

The multidisciplinary approach carried out across the two geologically complex investigation fields, confirmed the advantages of the combined application and interpretation of the electrical resistivity tomography (ERT) and seismic refraction tomography (SRT) techniques. Beyond the technical assessment of optimum data acquisition and processing parameters for these geophysical techniques, this approach provided an initial understanding and connection between the surface geological observations with the near-surface geotectonic subsurface structure of the two study areas.
The present study demonstrated that the application of the ERT technique is indicated as the primary technique for subsurface investigation in such complex geological environments. Compared to the SRT technique, it provides greater flexibility in the field, making it less time-consuming and more cost-effective. The ability to acquire measurements with different electrode configurations and choose among various optimization methods regarding the inversion process, allows researchers to adopt the best possible approach for the ERT investigation, depending on the expected variations in resistivity. Additionally, the ERT technique does not require an initial model input for data processing, making it an ideal primary methodology, prior to the design and application of the SRT technique. Moreover, the results of the ERT technique, which may be the only available source of subsurface information, can be utilized for modelling of the seismic survey. This process allows for assessing the effectiveness of the SRT technique in revealing the considered subsurface structure and helps determine the distribution of the first arrivals, especially in cases where the seismic records are noisy.
The technical characteristics of the seismic source utilized for the SRT technique play a significant role in the determination of the maximum investigation depth, quality and resolution of the results. For great exploration depths, a seismic source of high energy and frequency content is considered essential. Its technical characteristics shall ensure easy operation and transportation in the field along the different shot locations and repeatability of the generated signal as well. In the present study, the accelerated weight-drop impact seismic source (AWD-33PS) provided the highest energy amplitude, as determined from the signal analysis and comparison procedure. However, the signal produced by the 6.5kg seismic sledgehammer was satisfactory enough, particularly when the shot position was located within the active spread length, where less energy is required compared to the outshot positions. Regarding the explosive seismic sources, the signal’s characteristics of the seismic detonator and buffalo gun were very similar, with lower amplitudes compared to the impact sources. Therefore, these type of explosive sources are recommended for shallow refraction surveys, with the deepest refractor located at a maximum depth of 10m. To increase the investigation depth, a greater quantity of explosive material should be used. The generation of a comprehensive and high-resolution subsurface velocity model requires the recording of arrival times from many different seismic source locations, which possesses significant challenges in data organization and management. Moreover, seismic data processing demands considerable involvement and expertise from scientific analysts.
In conclusion, a thorough understanding of the surficial geological conditions of the study area, combined with the application of the ERT and SRT geophysical techniques and the potential availability of geological reference data (detailed mapping and boreholes), constitutes the optimal combination for developing an accurate and scientifically substantiated subsurface scenario.

Author Contributions

Conceptualization, J.D.A., K.S., V.G. and D.M.; methodology, V.G., D.M., J.D.A. and I.-K.G.; software, V.G., I.-K.G. and S.D.; validation, V.G., D.M., J.D.A. and N.V.; formal analysis, V.G., D.M. and I.-K.G.; investigation, V.G., J.D.A., K.S., S.D. and I.-K.G; resources, J.D.A., K.S., N.V. and T.S.; data curation, V.G., D.M., I.-K.G. and S.D.; writing—original draft preparation, V.G., K.S. and I.-K.G.; writing—review and editing, J.D.A., S.D., N.V., D.M. and T.S.; visualization, V.G., K.S. and I.-K.G.; supervision, J.D.A., N.V., K.S. and S.D.; project administration, J.D.A. and K.S.; funding acquisition, J.D.A. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Special Account for Research Grants of the NKUA, grant number 19498 and 13504.

Data Availability Statement

Data available upon request.

Acknowledgments

The authors would like to acknowledge Mitsika G., Balomenou G., Flambouris E., Petroulias I. and Tzanaki S., for their contribution during the field measurements.

Conflicts of Interest

The authors declare no conflicts of interest.

References

  1. Quigley, P. T. Ground proving seismic refraction tomography (SRT) in laterally variable karstic limestone terrain. Doctoral dissertation, University of Florida, 2006.
  2. Akingboye, A.S.; Ogunyele, A.C. Insight into seismic refraction and electrical resistivity tomography techniques in subsurface investigations. Rudarsko-geološko-naftni zbornik 2019, 34, 93–111. [Google Scholar] [CrossRef]
  3. Mazzotti, A.; Stucchi, E.; Fradelizio, G.; Zanzi, L.; Scandone, P. Seismic exploration in complex terrains: A processing experience in the Southern Apennines. Geophysics 2000, 65, 1402–1417. [Google Scholar] [CrossRef]
  4. Improta, L.; Zollo, A.; Herrero, A.; Frattini, R.; Virieux, J.; Dell’Aversana, P. Seismic imaging of complex structures by non-linear traveltime inversion of dense wide-angle data: Application to a thrust belt. Geophysical Journal International 2002, 151, 264–278. [Google Scholar] [CrossRef]
  5. Loke, M. H. Tutorial: 2-D and 3-D electrical imaging surveys, 2004; 128p.
  6. Improta, L.; Ferranti, L.; De Martini, P. M.; Piscitelli, S.; Bruno, P. P.; Burrato, P.; Civico, R.; Giocoli, A.; Iorio, M.; D'Addezio, G.; Maschio, L. Detecting young, slow-slipping active faults by geologic and multidisciplinary high-resolution geophysical investigations: A case study from the Apennine seismic belt, Italy. Journal of Geophysical Research 2010, 115. [Google Scholar] [CrossRef]
  7. Giocoli, A.; Stabile, T.A.; Adurno, I.; Perrone, A.; Gallipoli, M.R.; Gueguen, E.; Norelli, E.; Piscitelli, S. Geological and geophysical characterization of the Southeastern side of the high agri valley (Southern Apennines, Italy). Natural Hazards and Earth System Sciences 2015, 15, 315–323. [Google Scholar] [CrossRef]
  8. Giano, S. I.; Lapenna, V.; Piscitelli, S.; Schiattarella, M. Electrical imaging and self-potential surveys to study the geological setting of the Quaternary, slope deposits in the Agri high valley (Southern Italy). ANNALS OF GEOPHYSICS 2000, 43, 409–419. [Google Scholar] [CrossRef]
  9. Nguyen, F.; Garambois, S.; Jongmans, D.; Pirard, E.; Loke, M.H. Image processing of 2D resistivity data for imaging faults. Journal of Applied Geophysics 2005, 57, 260–277. [Google Scholar] [CrossRef]
  10. Diaferia, I.; Barchi, M.; Loddo, M.; Schiavone, D.; Siniscalchi, A. Detailed imaging of tectonic structures by multiscale Earth resistivity tomographies: The Colfiorito normal faults (Central Italy). Geophysical Research Letters 2006, 33. [Google Scholar] [CrossRef]
  11. Giocoli, A.; Burrato, P.; Galli, P.; Lapenna, V.; Piscitelli, S.; Rizzo, E.; Romano, G.; Siniscalchi, A.; Magrí, C.; Vannoli, P. Using the ERT method in tectonically active areas: Hints from Southern Apennine (Italy). Advances in geosciences 2008, 19, 61–65. [Google Scholar] [CrossRef]
  12. Kumar, D.; Subba Rao, D.V.; Mondal, S.; Sridhar, K.; Rajesh, K.; Satyanarayanan, M. Gold - sulphide mineralization in ultramafic-mafic-granite complex of Jashpur, Bastar craton, central India: Evidences from geophysical studies. Journal of the Geological Society of India 2017, 90, 147–153. [Google Scholar] [CrossRef]
  13. Alexopoulos, J. D.; Dilalos, S.; Vassilakis, E. Adumbration of Amvrakia’s spring water pathways, based on detailed geophysical data (Kastraki-Meteora). In Advances in the Research of Aquatic Environment; Springer Berlin Heidelberg: Berlin, Germany, 2011; Volume 2, pp. 105–112. [Google Scholar]
  14. Muchingami, I.; Hlatywayo, D.J.; Nel, J.M.; Chuma, C. Electrical resistivity survey for groundwater investigations and shallow subsurface evaluation of the basaltic-greenstone formation of the urban Bulawayo aquifer. Physics and Chemistry of the Earth, 2012, 50-52, 44–51. [Google Scholar] [CrossRef]
  15. Gkosios, V.; Alexopoulos, J. D.; Giannopoulos, I. K.; Mitsika, G. S.; Dilalos, S.; Barbaresos, I.; Voulgaris, N. Determination of the subsurface geological regime and geotechnical characteristics at the area of Goudi (Athens, Greece) derived from geophysical measurements. In Bulletin of the Geological Society of Greece, Special Publication GSG2022-062. In Proceedings of the 16th International Congress of the Geological Society of Greece, Patras, Greece, 17-19 October 2022. [Google Scholar]
  16. Alexopoulos, J.D.; Gkosios, V.; Dilalos, S.; Giannopoulos, I. K.; Mitsika, G. S.; Barbaresos, I.; Voulgaris, N. Assessment of near-surface geophysical measurements for geotechnical purposes at the area of Goudi (Athens, Greece). In NSG2023 29th European Meeting of Environmental and Engineering Geophysics, 1, 1-5, Edinburgh, UK, 3-7 September 2023. [CrossRef]
  17. Alexopoulos, J.D.; Dilalos, S.; Voulgaris, N.; Gkosios, V.; Giannopoulos, I. K.; Kapetanidis, V.; Kaviris, G. The contribution of near-surface geophysics for the site characterization of seismological stations. Applied sciences 2023, 13, 4932–4932. [Google Scholar] [CrossRef]
  18. Zhang, J.; ten Brink, U.S.; Toksöz, M.N. Nonlinear refraction and reflection travel time tomography. Journal of Geophysical Research: Solid Earth 1998, 103, 29743–29757. [Google Scholar] [CrossRef]
  19. Cardarelli, E.; de Nardis, R. Seismic refraction, isotropic anisotropic seismic tomography on an ancient monument (Antonino and Faustina temple Ad 141). Geophysical Prospecting 2001, 49, 228–240. [Google Scholar] [CrossRef]
  20. Villani, F.; Tulliani, V.; Sapia, V.; Fierro, E.; Civico, R.; Pantosti, D. Shallow subsurface imaging of the Piano Di Pezza active normal fault (Central Italy) by high-resolution refraction and electrical resistivity tomography coupled with time-domain electromagnetic data. Geophysical Journal International 2015, 203, 1482–1494. [Google Scholar] [CrossRef]
  21. Improta, L.; Zollo, A.; Bruno, P.P.; Herrero, A.; Villani, F. High-resolution seismic tomography across the 1980 (Ms 6.9) Southern Italy earthquake fault scarp. Geophysical Research Letters 2003, 30. [Google Scholar] [CrossRef]
  22. Improta, L.; Bruno, P.P. Combining seismic reflection with multifold wide-aperture profiling: An effective strategy for high-resolution shallow imaging of active faults. Geophysical Research Letters 2007, 34. [Google Scholar] [CrossRef]
  23. Galone, L.; Villani, F.; Colica, E.; Pistillo, D.; Baccheschi, P.; Panzera, F.; Zaldívar, J. G.; D'Amico, S. Integrating near-surface geophysical methods and remote sensing techniques for reconstructing fault-bounded valleys (Mellieha Valley, Malta). Tectonophysics 2024, 875, 230263. [Google Scholar] [CrossRef]
  24. Alexopoulos, J.D.; Voulgaris, N.; Dilalos, S.; Gkosios, V.; Giannopoulos, I. K.; Mitsika, G.S.; Vassilakis, E.; Sakkas, V.; Kaviris, G. Near-surface geophysical characterization of lithologies in Corfu and Lefkada Towns (Ionian Islands, Greece). Geosciences 2022, 12, 446. [Google Scholar] [CrossRef]
  25. Zulauf, G.; W. Dörr; Marko, L.; Krahl, J. The Late Eo-Cimmerian Evolution of the External Hellenides: Constraints from Microfabrics and U–Pb Detrital Zircon Ages of Upper Triassic (Meta)Sediments (Crete, Greece). International journal of earth sciences 2018, 107, 2859–2894. [Google Scholar] [CrossRef]
  26. Lekkas, S.; Skourtsos, E. The nappe structure of the tectonic window of Doliana (Central Peloponnesus, Greece). Bulletin of the Geological Society of Greece 2004, 36, 1662. [Google Scholar] [CrossRef]
  27. Voudouris, P.; Melfos, V.; Spry, P.G.; Bonsall, T.; Tarkian, M.; Economou-Eliopoulos, M. Mineralogical and fluid inclusion constraints on the evolution of the Plaka intrusion-related ore system, Lavrion, Greece. Mineralogy and Petrology 2008, 93, 79–110. [Google Scholar] [CrossRef]
  28. Lekkas, S.; Skourtsos, E.; Soukis, K.; Kranis, H.; Lozios, S.; Alexopoulos, A.; Koutsovitis, P. Late Miocene detachment faulting and crustal extension in SE Attica (Greece). In Geophysical Research Abstracts EGU General Assembly, Vienna, Austria, 3-8 April 2011.
  29. Grasemann, B.; Schneider, D.; Stöckli, D.F.; Iglseder, C. Miocene bivergent crustal extension in the Aegean: Evidence from the Western Cyclades (Greece). Lithosphere 2012, 4, 23–39. [Google Scholar] [CrossRef]
  30. Berger, A.; Schneider, D.A.; Grasemann, B.; Stockli, D. Footwall mineralization during late Miocene extension along the west Cycladic detachment system, Lavrion, Greece. Terra nova 2012, 25, 181–191. [Google Scholar] [CrossRef]
  31. Coleman, M.; Dubosq, R.; Schneider, D.A.; Grasemann, B.; Soukis, K. Along-strike consistency of an extensional detachment system, West Cyclades, Greece. Terra Nova 2019, 31, 220–233. [Google Scholar] [CrossRef]
  32. Marinos, G.P.; Petrascheck, W.E. Lavrion, geological and geophysical research. Institute for Geology and Subsurface Research: Athens, Greece, 1956; 246 p.
  33. Photiades, A.; Carras, N. Stratigraphy and geological structure of the Lavrion Area (Attica, Greece). Bulletin of the Geological Society of Greece 2001, 34, 103. [Google Scholar] [CrossRef]
  34. Skarpelis, N.; Tsikouras, B.; Pe-Piper, G. The Miocene igneous rocks in the basal unit of Lavrion (SE Attica, Greece): Petrology and geodynamic implications. Geological magazine 2007, 145, 1–15. [Google Scholar] [CrossRef]
  35. Liati, A.; Skarpelis, N.; Pe-Piper, G. Late miocene magmatic activity in the Attic-Cycladic belt of the Aegean (Lavrion, SE Attica, Greece): Implications for the geodynamic evolution and timing of ore deposition. Geological Magazine 2009, 146, 732–742. [Google Scholar] [CrossRef]
  36. Loke, M.H.; Dahlin, T. A. Comparison of the Gauss–Newton and Quasi-Newton methods in resistivity imaging inversion. Journal of Applied Geophysics 2002, 49, 149–162. [Google Scholar] [CrossRef]
  37. Silvester, P. P.; Ferrari, R. L. Finite elements for electrical engineers, 3rd ed.; Cambridge University press: England, 1996; 494p. [Google Scholar] [CrossRef]
  38. deGroot -Hedlin, C.; Constable, S. Occam’s inversion to generate smooth, two-dimensional models from magnetotelluric data. GEOPHYSICS 1990, 55, 1613–1624. [Google Scholar] [CrossRef]
  39. Claerbout, J.F.; Muir, F. Robust modeling with erratic data. GEOPHYSICS 1973, 38, 826–844. [Google Scholar] [CrossRef]
  40. Olayinka, A.I.; Yaramanci, U. Assessment of the reliability of 2D inversion of apparent resistivity data. Geophysical Prospecting 2000, 48, 293–316. [Google Scholar] [CrossRef]
  41. Loke, M.H.; Acworth, I.; Dahlin, T. A comparison of smooth and blocky inversion methods in 2D electrical imaging surveys. Exploration Geophysics 2003, 34, 182–187. [Google Scholar] [CrossRef]
  42. Moser, T.J. Shortest path calculation of seismic rays. GEOPHYSICS 1991, 56, 59–67. [Google Scholar] [CrossRef]
  43. Knapp, R.W.; Steeples, D.W. High-resolution common-depth-point reflection profiling: Field acquisition parameter design. GEOPHYSICS 1986, 51, 283–294. [Google Scholar] [CrossRef]
  44. Feroci, M.; Orlando, L.; Balia, R.; Bosman, C.; Cardarelli, E.; Deidda, G. Some considerations on shallow seismic reflection surveys. Journal of Applied Geophysics 2000, 45, 127–139. [Google Scholar] [CrossRef]
  45. Evans, B. A handbook for seismic data acquisition in exploration. Society of Exploration Geophysicists: Houston, USA, 1997. [CrossRef]
  46. Steeples, D.W. A review of shallow seismic methods. Annals of Geophysics 2000, 43. [Google Scholar] [CrossRef]
  47. Miller, R.D.; Pullan, S.E.; Waldner, J.S.; Haeni, F.P. Field comparison of shallow seismic sources. GEOPHYSICS 1986, 51, 2067–2092. [Google Scholar] [CrossRef]
  48. Doll, W.E.; Miller, R.D.; Xia, J. A noninvasive shallow seismic source comparison on the oak ridge reservation, Tennessee. GEOPHYSICS 1998, 63, 1318–1331. [Google Scholar] [CrossRef]
  49. Bühnemann, J.; Holliger, K. Comparison of high-frequency seismic sources at the Grimsel test site, Central Alps, Switzerland. GEOPHYSICS 1998, 63, 1363–1370. [Google Scholar] [CrossRef]
  50. Aritman, B.C. Repeatability study of seismic source signatures. GEOPHYSICS 2001, 66, 1811–1817. [Google Scholar] [CrossRef]
  51. Aas, A.; Sinha, S.K. A novel weight-drop seismic energy source for subsurface characterization. Journal of Applied Geophysics 2023, 208, 104887. [Google Scholar] [CrossRef]
  52. Caselle, C.; Bonetto, S.; Comina, C.; Stocco, S. GPR surveys for the prevention of karst risk in underground gypsum quarries. Tunnelling and Underground Space Technology 2020, 95, 103137. [Google Scholar] [CrossRef]
  53. Daniels, D.J. Ground Penetrating Radar, 2nd ed.; The Institution of Electrical Engineers: London, UK, 2004; 734p. [Google Scholar]
Figure 5. Snapshots from data acquisition on the Kleisroua Valley, Ano Doliana (left) and Plaka (right) study areas.
Figure 5. Snapshots from data acquisition on the Kleisroua Valley, Ano Doliana (left) and Plaka (right) study areas.
Preprints 113417 g005
Figure 6. Detailed geological map of the 1st study area (Kleisoura Valley, Ano Doliana), along with the acquisition layout of the geophysical techniques.
Figure 6. Detailed geological map of the 1st study area (Kleisoura Valley, Ano Doliana), along with the acquisition layout of the geophysical techniques.
Preprints 113417 g006
Figure 7. Detailed geological map of the 2nd study area (Plaka), along with the acquisition layout of the geophysical techniques. Modified after [27].
Figure 7. Detailed geological map of the 2nd study area (Plaka), along with the acquisition layout of the geophysical techniques. Modified after [27].
Preprints 113417 g007
Figure 8. (a) ERT results of the 1st study area, based on which the synthetic velocity model was created; (b) Synthetic velocity model with the seismic wavefronts generated at the Outshot Normal 1 position (120m offset), with the geophones located between 157.5 – 392.5m distance (red marks); (c) Creation of the synthetic seismic record for the Outshot Normal 1 source location. First arrival times are noted by the red arrows.
Figure 8. (a) ERT results of the 1st study area, based on which the synthetic velocity model was created; (b) Synthetic velocity model with the seismic wavefronts generated at the Outshot Normal 1 position (120m offset), with the geophones located between 157.5 – 392.5m distance (red marks); (c) Creation of the synthetic seismic record for the Outshot Normal 1 source location. First arrival times are noted by the red arrows.
Preprints 113417 g008
Figure 9. (a) ERT results of the 1st study area, in relation with the (b) seismic tomogram generated after the inversion of the 17 synthetic seismic records.
Figure 9. (a) ERT results of the 1st study area, in relation with the (b) seismic tomogram generated after the inversion of the 17 synthetic seismic records.
Preprints 113417 g009
Figure 10. Superimposing of the synthetic first arrival times (red marks) on the field seismic records for the (a) Outshot normal 1, (b) Tomo shot 1 and (c) Tomo shot 4 seismic source locations (Figure 6). The blow arrows and the blue dashed lines indicate the divergence between the synthetic and true first arrivals.
Figure 10. Superimposing of the synthetic first arrival times (red marks) on the field seismic records for the (a) Outshot normal 1, (b) Tomo shot 1 and (c) Tomo shot 4 seismic source locations (Figure 6). The blow arrows and the blue dashed lines indicate the divergence between the synthetic and true first arrivals.
Preprints 113417 g010
Figure 11. Seismic records along with their corresponding amplitude spectra acquired with four (4) different seismic sources at the Outshot Normal 2 source location (Figure 6). AWD:20kg – Accelerated Weight Drop; SH:6.5kg – Sledgehammer; SD - Seismic detonator; BG – Buffalo gun.
Figure 11. Seismic records along with their corresponding amplitude spectra acquired with four (4) different seismic sources at the Outshot Normal 2 source location (Figure 6). AWD:20kg – Accelerated Weight Drop; SH:6.5kg – Sledgehammer; SD - Seismic detonator; BG – Buffalo gun.
Preprints 113417 g011
Figure 12. (a) ERT and (b) SRT profiles and their (c) geological interpretation derived from their combined evaluation, for the 1st study area (Kleisoura Valley, Ano Doliana).
Figure 12. (a) ERT and (b) SRT profiles and their (c) geological interpretation derived from their combined evaluation, for the 1st study area (Kleisoura Valley, Ano Doliana).
Preprints 113417 g012
Figure 13. (a) ERT and (b) SRT profiles and their (c) geological interpretation derived from their combined evaluation, for the 2nd study area (Plaka).
Figure 13. (a) ERT and (b) SRT profiles and their (c) geological interpretation derived from their combined evaluation, for the 2nd study area (Plaka).
Preprints 113417 g013
Figure 14. GPR profile of the 2nd study area. Underground air voids (possible mining galleries) are noted by the red dashed circles.
Figure 14. GPR profile of the 2nd study area. Underground air voids (possible mining galleries) are noted by the red dashed circles.
Preprints 113417 g014
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
Copyright: This open access article is published under a Creative Commons CC BY 4.0 license, which permit the free download, distribution, and reuse, provided that the author and preprint are cited in any reuse.
Prerpints.org logo

Preprints.org is a free preprint server supported by MDPI in Basel, Switzerland.

Subscribe

© 2024 MDPI (Basel, Switzerland) unless otherwise stated