Geophysical Characterization of Aquifers in Southeast Spain Using ERT, TDEM, and Vertical Seismic Reflection

Abstract

We assess the effectiveness of complementary geophysical techniques to characterize a Jurassic dolomite confined aquifer at Loma de Ubeda, Spain. This aquifer, which is penetrated by wells in the 100–600-m depth range, is confined by Triassic clays (bottom) and Miocene marls (top). The Jurassic dolomite is characterized by prominent seismic reflectors of high amplitude. Thus, it is readily differentiated from the low-amplitude reflectors of the confining clay-rich Triassic and Miocene materials. Electrical resistivity tomography (ERT) allowed us to detail the characteristics of the aquifer up to a maximum depth of 220 m. Lateral changes in facies and small faults have been identified using ERT. Time-domain electromagnetic (TDEM) is an excellent complement to the two above-mentioned techniques in order to widen the analyzed depth range. We acquire TDEM data with different configurations at multiple study sites while simultaneously varying measurement parameters. In doing so and by comparing the effectiveness of these different configurations, we expand the use of TDEM for aquifer characterization.

1. Introduction

Geophysical prospecting techniques for groundwater exploration have evolved considerably in recent years. The automatic data acquisition systems, the availability of powerful computers, and the development of two- and three-dimensional modeling software have greatly improved the resolution of complex geological models [1,2,3,4].

A great diversity of geophysical methods (gravimetric, magnetic, electrical, electromagnetic, and seismic) can contribute to the groundwater exploration. Electrical, electromagnetic, and seismic methods are the most efficient techniques in this kind of setting and they are more cost-effective.

Electrical methods are particularly suitable for groundwater research because hydrogeological parameters, such as porosity and permeability, can be correlated to electrical resistivity values. In addition, they are cheap and easy to implement techniques [5]. Electromagnetic sounding techniques, on the other hand, have not been used as extensively as electrical methods for groundwater research because the equipment is more expensive and the interpretation methods are more complicated. The use has focused on detecting the interface between fresh groundwater and saline seawater intrusion into coastal aquifers [6]. These two techniques have been used recently in an unconfined sandstone aquifer in semiarid southwestern Niger [7]. This study is a good example of the efficiency of using complementary geophysical techniques.

Herein, we analyze the effectiveness of electrical (electrical resistivity tomography, ERT), electromagnetic (time-domain electromagnetic, TDEM), and reflection seismic techniques to characterize dolomite aquifers. Our study area (Figure 1) not only has detailed surface geological information, but also numerous pumping boreholes and seismic lines, which were employed to verify and correlate the new geophysical data collected in this research.

The economy of Loma de Úbeda area in Southern Spain is based almost exclusively on the cultivation of olive groves. Groundwater exploitation has increased olive production and ensured the survival of crops during the severe droughts that affect this region.

Until the late 1980s, only a Miocene shallow aquifer was exploited, with pumping flows on the order of 1 L/s. In the early 1990s, the Loma de Úbeda Jurassic dolomite aquifer began to be exploited at depths between 100 m and 600 m, with pumping rates exceeding 40 L/s in some cases [8]. Numerous deep boreholes were drilled, many of them without the mandatory permits from the Water Authorities. Estimated irrigation in 2005 exceeded 35 hm³/year over an irrigation area of 24,050 ha of olive groves [9]. The aquifer has an average annual deficit close to 13.3 hm³ and it is clearly overexploited [10]. We applied different geophysical techniques to make comparisons with existing boreholes and surface geologic data. This allowed us to verify the effectiveness of each technique to investigate this aquifer. We emphasized TDEM and we used different configurations to find the optimal parameters. The different techniques used here can be applied to investigate other carbonate aquifers showing similar behavioral.

2. Description of the Study Region

2.1. Geological Setting

Two geological units can be differentiated in the study area: the Paleozoic basement and the post-Hercynian sedimentary cover. The first unit is dominated by intensely folded phyllites intruded by a granitic batholith (Figure 1: 1,2). Subhorizontal to gently dipping, the post-Hercynian sedimentary cover unconformably overlies the basement. It is composed of Triassic, Jurassic, and Neogene formations (Figure 1). Sparse Quaternary alluvial sediments are also present.

The Triassic strata belongs to the Chiclana de Segura Formation [11,12], which extensively outcrops at the north of the Guadalimar River (Figure 1: 3). This formation varies in thickness (50–400 m) and is essentially horizontal or dips gently southward. In these rocks, their reddish hues and the presence of shales and sandstones are characteristic (Figure 2). The lower part of the Triassic series is mostly composed by sandstones levels. Towards the upper-middle part, the number and thick of clayey facies increases, and gypsum and salt levels intercalations appear (Figure 2).

In the Guadalimar River valley, the Jurassic dolomites directly outcrop over the Triassic materials in apparent stratigraphic continuity (Figure 1: 4). They appear subhorizontal or with small dips, generally towards the south–southeast. These are strongly brecciated and dolomitized carbonates, thus causing an important secondary porosity that must affect the aquifer storability. The thickness, highly variable in the 0–70-m range, generally decreases towards the west. Although there are no biostratigraphic criteria, this lithological unit has been attributed to the Lower Jurassic following facies and stratigraphic position criteria [12,13].

Jurassic strata dip slightly towards the south–southeast and they folded gently on NE–SW axes [5]. On the other hand, two NE–SW and NNW–SSE fault directions in both the basement and the dolomites. The orientations of theses faults overall influenced large-scale regional patterns of early Tortonian sedimentation [14,15,16]. Jurassic strata pinch out towards the west–northwest. Several pumping boreholes near Rus-Canena document this trend (Figure 1). West of Ibros, Jurassic strata are absent and the Miocene strata rest directly over the Triassic ones. The disappearance of Jurassic dolomites also marks the western limit of the aquifer [8].

Three Miocene units are differentiated (Figure 1: 6): at the base, the Early Tortonian, which crops out locally, filling a small trough and pinching out laterally (Figure 2). These are gravel-rich sediments with dolomite fragments and siliceous sands, with little or no cementation, alternating with marls [14]. Late Tortonian strata either directly overlie the Early Tortonian ones or they overlie the Jurassic dolomites. They average 200 m in thickness, attaining a maximum thickness of 500 m. These strata consist chiefly of marls and marly limestone, but there are also minor calcarenite intercalations [13,17]. Messinian strata with more abundant calcarenites overlie the Late Tortonian strata [13,17] (Figure 2).

2.2. Hydrogeological Context

The Jurassic dolomites constitute an important hydrogeological unit underlying nearly 800 km² (Figure 1: 4). This aquifer behaves as unconfined in the northwestern portion of the study area (Figure 2). Towards the south, it is confined by Miocene marls, while it is limited by dipping faults of the Prebetic units to the east, by the olistostromic units of the Guadalquivir River depression to the south and west it pinches out [8,18]. In the western sector, the general groundwater flow path is NNW to SSE. Fracture systems compartmentalize the aquifer, changing piezometry levels and groundwater flow paths from place to place [18].

Diffuse recharge occurs by direct infiltration from precipitation whereas preferential recharge from the Guadalimar River streamflow infiltration occurs in those northern sectors where riverbed intersects the carbonate formations (Figure 1). It is also worth considering the existence of lateral transferences from eastern Prebetic units (Figure 1: 5).

The exploitation of the confined sector of the aquifer began in the early 1990s. The mechanical drilling methods before that data did not exceed 100 m depth from shallow Miocene strata only. The Jurassic dolomite aquifer has been exploited with pumping rates between 15 and 40 L/s only since 10 years ago thanks to deeper boreholes with depths in the 200–600-m range. Today, the pumping rate is much lower given the overexploitation that the aquifer unit suffers.

Different groundwater chemical facies can be differentiated. Shallow boreholes (dolomites pumped at 200–300 m) generally yield slightly mineralized waters. Salinity increases considerably with depth, and deeper chemical facies pass from magnesium-calcium sulfate, sodium bicarbonate-sulfate, sodium bicarbonate, to sodium chloride [18]. Chemical reduction processes in the deeper confined parts [18] and mixing of groundwater from Triassic, Jurassic, and Miocene aquifers have been described [9,19].

3. Methods

3.1. Seismic Reflection

Chevron Oil Company of Spain acquired seismic profiles in southeastern Spain in 1983 and 1984. These surveys used arrays of 18 geophones (two chains of nine geophones per trace) at intervals of 20 m. The processing flow followed a standard processing sequence, including migration. The information was processed by the ‘Compagnie Générale de Géophysique’. The integrated information of the different seismic lines, together with the direct data from the pumping boreholes, allowed a first reconstruction of the 3D geometry of the Jurassic dolomite unit [8].

This study focuses on seismic profile S84-68 (Figure 2). A correction speed of 2000 m/s was applied to the two-way travel time/depth transformation in Miocene strata. This is the wave velocity that is marked in the processing of the line. The thicknesses that result applying this wave velocity are correlated in pumping boreholes.

3.2. Electrical Resistivity Tomography

This geophysical prospecting technique consists of determining the distribution of the electrical resistivity of the subsoil from a very large number of measurements collected from the ground surface. The electrical resistivity quantifies how a material resists or conducts electric current [20,21]. The different electrical behavior of geological materials allows us to obtain 2D resistivity models, making ERT one of the most effective non-destructive tools to study and characterize subsurface discontinuities [20,22]. In recent years, this technique has attained widespread success in stratigraphic, hydrogeological, and environmental studies [3,23,24,25,26,27,28,29].

ERT involves the installation of numerous electrodes along a line (profile), with a given separation that determines the resolution and prospecting depth to be reached. A smaller separation between electrodes increases resolution, whereas a larger separation increases prospecting depth [22]. To make one measurement, only four electrodes are needed. Two of them act as ‘current electrodes’ and the other two act as ‘potential electrodes’. The way in which the electrodes are selected is termed ‘electrode configuration’. In stratigraphic studies, the Wenner–Schlumberger configuration is usually selected [30,31]. This configuration has good behavior and stability against resistivity changes, both vertical and horizontal ones, so it is useful for the investigation of horizontal or slightly inclined layers that can present lateral changes of facies and/or verticalized structures, as in our case [30,31].

Electrodes are connected to the measuring equipment (resistivity-meter), and through a predefined sequence, the groups of electrodes are selected. For each electrode quadrupole, a voltage and an intensity measurement are made. With these two readings, the ‘apparent’ resistivity of the ground is calculated, which is attributed to a certain geometric point in the subsoil.

We used a multi-channel, multi-electrode DC resistivity meter system manufactured by Deutsche Montan Technologie (RESECS model). The apparent resistivity values measured in the field were inverted to obtain electrical resistivity models by using the RES2DINV software [2]. This software applies a least-squares method with a damped smoothing, modified with the quasi-Newton optimization method. This inversion method constructs a subsoil model using rectangular prisms and determines the resistivity values for each of them, minimizing the difference between the observed and calculated apparent resistivity values [2,32]. We produced a NW–SE electrical tomography profile, with a total length of 1110 m using 112 electrodes spaced at 10 m (Figure 2).

3.3. Time-Domain Electromagnetic

The operating principle of the time-domain electromagnetic (TDEM) method is to circulate an electric current through a transmitter coil (usually square in shape) for short time intervals. When the current flow is abruptly interrupted, a magnetic field is produced that induces, according to Faraday’s law, a variable electric current in the subsurface, which in turn generates a transient secondary magnetic field. These currents flow in closed paths and migrate at depth, decreasing in intensity over time. Changes in the secondary magnetic field over time induce a transient voltage in the receiving coil. The shape of the decay of this voltage provides information on the conductivity distribution of the subsoil, which can be used to characterize it [21,33,34,35]. TDEM has been typically applied to finding mineral deposits [36,37,38], to investigate groundwater bodies [35,39], and to analyze sedimentary basins [40,41]. In recent years, this method has also been used in environmental studies [42,43] and to characterize marine intrusion [44,45].

The maximum prospecting depth depends on the decay time of the signal from the current cut, the current intensity, the loop size, the signal frequency, and the conductivity of the subsoil [33,34,35]. Although there are various measurement configurations, the most common one is to place the (smaller) receiving coil in the center of the (much larger) transmitting coil (central-loop configuration). Another option is to use the same coil to carry out both functions (single-loop configuration).

We used a TEM measuring system developed and produced by ELGEO Research & Production Company (AIE-2 model). It is a device with a maximum output power of 200 W and a current intensity of up to 10 A. The TDEM receiver is based on a 16-bit analogue-to-digital converter and a signal processor that provides an analogue-to-digital conversion immune to input voltage noise and real-time signal pre-processing. It offers a range of measurements between 5 µs and 10 µs, with an input voltage of 5–20 V with signal compensation.

In the field campaign, a total of 10 TDEM measurement stations were acquired along a north–south profile. The points located in Figure 2 refer to the center of the measuring stations for square loops of 200 m on each side. At each of these points, measurements were collected with central-loop and single-loop configurations, and the measurement parameters (intensity, time, voltage) were varied, which allowed us to compare the effectiveness of the different setups. To reach prospecting depths on the order of 400 m, transmission loops of 200 m × 200 m are needed. A 6-mm² section copper wire was used. This (larger than usual) wire decreased the resistance and thereby allowed it to obtain a higher effective current intensity. When the central-loop configuration was used, the receiving coil had a dimension of 20 m × 20 m, using a 5-turn cable and a ×10 amplifier. For the visualization and editing of the different curves, the TEMBIN software was used. The modeling and inversion processes were performed with ZondTEM1D and ZondTEM2D software http://zond-geo.com/english/zond-software/electromagnetic-sounding/zondtem1d/.

4. Results and Discussion

4.1. Seismic Reflection

The sedimentary cover sequence at the top of the seismic profiles produces reflectors of moderate amplitude, which we interpreted as Miocene sandstones and marls. The seismic response of the Jurassic dolomites consists of prominent reflectors of high amplitude and continuity. This different response allows us to trace both the upper and lower boundaries of the Jurassic aquifer. Under the Jurassic unit, locally discontinuous reflectors of low amplitude are associated to the Late Triassic unit, which contains clays and evaporates (gypsum and halite). At the base of the Triassic unit, continuous levels of greater amplitude are detected, which we have correlated with the basal, thick Triassic sandstones levels. The acoustic material is the Paleozoic basement, which is characterized by chaotic features and disorganized seismic facies.

Different criteria to calculate the stratigraphic units’ thickness were used. A criterion, namely the measured thicknesses of stratigraphic units, was measured in the field. The other criterion was velocity values included with the profiles. These two criteria generally coincided. Finally, surveys carried out in the area provide some direct assessment of local hydrogeological conditions.

The seismic profile selected in this study (Figure 2 and Figure 3A) shows Jurassic dolomites dip approximately 20° S–SE. The unit outcrops in the northern sector and occupies deeper positions in the southernmost sectors. The resolution of the profile is very low near the land surface, so the dolomites are not readily identifiable to the north (Figure 3A).

Using a correction speed of 2000 m/s, the dolomites can be identified at 200 m depth near the pumping boreholes, as pointed out in Figure 3A. This depth is confirmed from borehole drilling information. Under the dolomites, higher-amplitude reflectors in the lower part of the Triassic unit correspond to sandstones. Some of the faults detected with ERT techniques affecting the sedimentary cover are also identified by seismic methods.

4.2. Electrical Resistivity Tomography

The ERT profile was acquired at the northernmost end (Figure 2), just where the dolomite aquifer is too shallow to be discerned in the seismic profile (Figure 3A). Figure 3B shows the ERT profile, which reaches a depth of 220 m in the maximum penetration zone. From geoelectrical records, three resistivity ranges are differentiated. The first is characterized by high resistivity values in the 100–500 Ωm range and corresponds to the Jurassic dolomite aquifer. These geoelectrical facies outcrops at the northern profile end and dip gently towards the south–southeast. The ERT results trace this facies up to 500 m where it suddenly ends, and this could be explained by the presence of a normal fault. In the local outcrops, normal faults affecting the Jurassic dolomites have been detected. In the same ERT profile, normal faults affecting this unit can be deduced. Field observations and geological maps (Figure 2) corroborate these interpretations. This Jurassic aquifer crops out in the northern sector and it is also found under the Miocene strata in the southern sector.

Over the dolomite aquifer, the ERT profile also shows a discontinuous unit of variable thickness characterized by resistivity values in the 40–80 Ωm range. This electrical behavior of this second unit may correspond to Early Tortonian sandy and gravel-rich sediments with dolostone fragments. The presence of pre-Late Tortonian faults determines the accommodation space for these syntectonic facies (Figure 2), as interpreted in previous studies [14]. In the northernmost sector, some of these facies having similar resistivity values could also correspond to weathered dolomites, conferring some uncertainty in the interpretation.

A third unit characterized by resistivity values in the 5–20 Ωm range can be identified at shallow depth from 400 m onwards. Based on its electrical behavior, this unit is interpreted to be Tortonian marls and marly limestones, which have great thickness in the southern sector. As can be deduced from the geological map (Figure 2), the Messinian facies are absent in the ERT profile.

The resistivity of the different lithologies is related to different physical parameters of the geological materials, such as texture, pore fluids, density, etc. [21]. As deduced, if knowledge of the subsoil is scarce, results are more uncertainty. In contrast, if the information available is abundant, these techniques can have multiple applications, for instance to identify lateral changes of facies, thickness variations of stratigraphic units, water-table positions, or degree of saline contamination of an aquifer.

4.3. Time-Domain Electromagnetic

For each of the measurement stations of the TDEM profile, the configuration parameters of loops, current injection, voltage, or measurement time were varied. Figure 4A shows the different induced voltage decay curves measured at station 1 (S1) for a single-loop configuration and obtained using the TEMBIN software.

In each of the tests, the ramp time is directly proportional to the current used [46]. In our study, according to the technical specifications of the equipment, the ramp time was 12 µs or 30 µs, depending on whether the injected current was 1 or 8 A (for loops of 200 m × 200 m). In station S1, when comparing the curves ‘C_0.bem’ and ‘C_1.bem’ (Figure 4A), currents of low amperage led the plateau at the top of the curve to be lower, and this increases with increasing amperage. Therefore, the information of the surface part decreases with the amplitude of the current. However, the capacity for in-depth investigation increases with the magnetic moment (I ∗ A). According to Spies [34]:

$d~0.55{\left(\frac{I\ast A}{\sigma \ast Nr}\right)}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$5$}\right.}$

where: d is depth in m, I is current of the transmitter in Amps, A the effective area of the transmitter in m², σ is electrical conductivity of the medium in S/m, and Nr is noise level in V/m².

By keeping the surface of the transmitter loop constant (Tx = 200 m × 200 m), the magnetic moment and the penetration capability remain proportional to the current, as deduced when curves ‘C_0.bem’ and ‘C_1.bem’ are compared (Figure 4A). Therefore, high currents are used even if information is lost near the top of the profile. Such information from very shallow depths is not relevant to the hydrogeological interpretation.

Another important aspect to consider is the voltage range. We compare the response with the measurements at 0.1, 1, and 10 V at station S1 (Figure 4A). The superposition of the different curves allows us to deduce that greater stability is achieved in the depth of the induced voltage curve using low voltages as 0.1 V in the ‘C_3.bem’ curve of Figure 4A.

Another factor analyzed was the measurement time window; longer time windows produce larger depths of investigation [34]. In our case, different time-offs at station S1 (20, 100, and 500 ms) were used. When the ‘C_2.bem’ and ‘C_4.bem’ curves are compared, times greater than 100 ms are observed, the noise level increases considerably, and the voltage values are extremely low, so time intervals in this range should be ignored due to their low reliability.

All the aforementioned curves were obtained with a single loop. The readings were repeated at station S1 (Figure 4B) for the measurements with the different set-ups but using a central-loop device. As in the previous case, the data for research at deeper levels show that the magnetic moment must be increased (if the area remains constant, it will be proportional to the increase in current) and the voltage must be decreased. The measurement time was held constant at 100 ms.

In all the stations indicated in the profile (Figure 2), the previously mentioned parameters were used. Figure 5 shows the fit between the induced voltage curves measured in the field (dots) and those generated (curves) by the resistivity model (thick red line). The figure also shows the apparent resistivity obtained as a function of time in the first and last of these stations (S-1 and S-10). For each of these curves, the fit between the field and the model (by the ZONDTEM1D software) lines is observed.

The curves observed in the field at station S1 and the modeled ones (Figure 2) are represented and compared when we use a central-loop device (‘Center 1’ in Figure 5A) and a single-loop device (‘C_3’ in Figure 5B). First, the good fit of the curves stands out in both cases, with a root mean square (RMS) error of 2.7% and 3.3%, respectively. After the data inversion in both cases, similar results in the resistivity/depth curves are observed. In both cases, at shallow depths, there is an increase in resistivity values that correlates with the Jurassic dolomite aquifer, unit that outcrop in nearby sectors (Figure 2). Under this unit, a set of more conductive facies appears that would correspond to the Triassic clays. Within the Triassic unit, generally in the upper part, the increase in resistivity can be associated with the increase in evaporitic layers (Figure 5A,B). Below 400 m, the sharp increase in resistivity is associated with phyllites of the Paleozoic basement. This interpretation is corroborated by the information of the thicknesses observed in the nearby outcrops (Figure 2B) and to the information provided by the seismic profile in Figure 3A. Therefore, TDEM curves could determine the thickness of a geological unit accurately. However, TDEM curves do not allow to establish the correct resistivity value for each level, as other authors have already pointed out. In this study, Figure 5 represents resistivity values after inversion; their abrupt increases or decreases allow deducing the presence of different lithologies. Nevertheless, the modeled resistivity values do not correspond to the real resistivity of the proposed lithologies.

It should be noticed that the shape of the curve, with changes in the tendency of increase or decrease in resistivity, allows detecting lithological changes. However, absolute values should not be taken into account: the dolomites offer average resistivity values in the 200–400 Ωm range (Figure 5), that can be lower than 100 Ωm in some stations (S1, S8, S9, and S10) and greater than 10⁵ Ωm in others (S4). On the other hand, despite the overall similarity in the models of Figure 5A,B, a significant difference must be highlighted: the slope change of the decay curve of the induced potential. In the case of the central-loop device, the slope of the curve is lower, reaching depths of approximately 600 m with induced potential data on the order of 10 µV. With much lower voltage values, the curve continues to show a highly stable trend. However, using the single-loop device, the slope of the potential decay curve is greater. Thus, from approximately 150 m depth, the induced voltage values would be less than 10 µV.

The same comparison is performed at station S10 (see its location in Figure 2). Again, the good fit of the curves stands out, having an RMS error of 0.6% with a central loop (‘qcenter_8’ in Figure 5C) and 1.6% with a single loop (‘q_2’ in Figure 5D). With the central-loop device, three units are differentiated in depth: the shallowest (around 100 m depth) is characterized by a drop in resistivity values, although the local rises appearance would be associated to Messinian marls with sandstones interbedded. Between 100 and 300 m, conductive facies appear associated to Late Tortonian marls. Under these facies, the increase in resistivity values is associated to Jurassic dolomite aquifer. This interpretation also agrees with the information from direct data obtained in nearby outcrops (Figure 2B) and with the information provided by the seismic profile of Figure 3B.

In the case of the central-loop device, the potential curve reaches approximately a depth of 510 m with stable values (at approximately 150 m, the voltages are approximately 10 µV). In contrast, using the single-loop device, it offers values above 10 µV only up to 65 m, and the curve is no longer stable at a depth of 250. Therefore, at this station, the Jurassic aquifer is not reached with this configuration. From all this, it can be deduced that the central-loop device performs better than the single-loop device for deep aquifer research (Figure 5 and Figure 6).

For the rest of the stations (S2 to S9), measurements with both configurations were made. For brevity, only the curves obtained by the central loop are shown (Figure 6). In all of them, the position of the Jurassic aquifer is detected at variable depths. In Figure 3A, the position of each of the TDEM stations on the seismic profile is shown. This allows us to compare the depth of Jurassic dolomites obtained by seismic (Figure 3) and electromagnetic methods (Figure 5 and Figure 6).

Figure 7A shows the resistivity section obtained with the electromagnetic method, resulting from the joint inversion of data from the ten measurement stations through the ZondTem2D software. At the northern end of the section, the method is able to detect both the Jurassic dolomites and the contact with the underlying Paleozoic basement. In the southern sector, the Jurassic unit is also detected by a change in the trend of conductive to resistive materials. The method fits poorly in the central sector of the profile, where is unable to detect the decrease in resistivity values associated to the Triassic facies at deeper levels. The seismic reflection profile and the resistivity profile obtained by the TDEM method are represented together in Figure 7B. A good correlation can be observed between the two techniques for modeling the Jurassic unit. For the upper 200 m, at the northern end of the section, the TDEM resistivity profile and the profile obtained by the electrical resistivity tomography (ERT) method are also compared (Figure 7C). The superposition of the ERT and TDEM profiles allows deducing a good correlation in the lateral and vertical variations of the resistivity values obtained in both techniques.

5. Conclusions

In the seismic reflection profiles, abrupt changes in the amplitude of the wave associated to the upper and lower boundaries of the Jurassic dolomite aquifer are observed, making it easy to identify thickness variations at depth. Due to attenuation, such sudden changes in amplitude cannot be discerned in places where the aquifer is very shallow. At depths not exceeding 200 m, variations in resistivity values within ERT profiles allow identifying this aquifer level clearly. Simultaneously, the technique allows deductions of the lateral changes in facies, the thickness of the different units, structural dips, and the presence of faults.

The TDEM technique is versatile and can be used with different amplitudes, voltages, times, devices, and loop sizes, depending on the target penetration depths at which the research is scheduled. In our study, to reach more than 300 m depth, low voltages, high magnetic moments (increasing the amplitudes and the area of the loop), and central-loop-type devices were used. The thickness of the layer can be determined using TDEM curves, but not its real resistivity. In general, there is a good fit between seismic, electrical, and electromagnetic data.

Author Contributions

Conceptualization and investigation, J.R., J.M., R.M., S.S., and K.M.; software, V.T. and A.K.; writing and review, J.R., S.S., and M.C.H. All authors have read and agreed the published version of the manuscript.

Funding

This research received no external funding.

Acknowledgments

Careful reviews by R.M. Joeckel and by anonymous referees substantially improved the manuscript.

Conflicts of Interest

The authors declare no conflict of interest.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Figures

Figure 1. Regional hydrogeological map. The location of the study area is shown. Legend: (1) Paleozoic basement (Phyllites). (2) Paleozoic basement (Granite). (3) Triassic clays from Chiclana de Segura Formation (Impervious base). (4) Jurassic tabular dolomites (Aquifer). (5) Prebetic Jurassic dolomites (Aquifer). (6) Late Miocene marls with detrital levels. (7) Pliocene–Quaternary clays, marls, and conglomerates. (8) Undifferentiated stratigraphic contact. (9) Undifferentiated tectonic contact. (10) Groundwater flow paths.

Figure 2. (A) Geological mapping of the studied sector. The positions of the electrical resistivity tomography (ERT) line of the seismic reflection profile S84-68 and the time-domain electromagnetic (TDEM) are marked. For this profile, 10 measuring stations (S1 to S10) were set up. Legend: (a) Triassic clays and sandstones (Chiclana de Segura Formation), (b) Jurassic dolomites, (c) Miocene (Tortonian) marl and marly limestones, (d) Miocene (Messinian) sandstones and marls, (e) boreholes, and (f) normal faults. (B) Synthetic stratigraphic column. Legend as (1) phyllites, (2) rudites, (3) sandstones, (4) clays, (5) gypsum, (6) dolomites, (7) gravel-rich sediments, (8) marls, and (9) marls and marly limestones.

Figure 3. Seismic reflection profile and interpretation (see location in Figure 2). The positions of the pumping boreholes and the TDEM (S1 to S10) measuring stations are indicated. M.a.s.l. (metres above sea level) (A). Electrical resistivity tomography (ERT) profile (B).

Figure 4. Different induced voltage curves as a function of time, all collected at station S1. Different measurement parameters were used in each of the curves (variations in Current, Voltage, Time). (A) The curves represent a single-loop configuration so that the transmitter area (Tx) and the receiver area (Rx) coincide (200 m × 200 m = 40,000 m2). (B) The curves represent a central-loop configuration, where the same transmitter loop is maintained, but the receiver loop (Rx) has dimensions of 20 m × 20 m; in this case, a five-wire cable and an amplifier of ×10 were used, so the effective area of Rx in this case is 20 × 20 × 5 × 10 = 20,000 m2.

Figure 5. Induced voltage curves (blue line) and apparent resistivity (red line) as a function of time. For each of them, the fit between the field and the model (by the ZONDTEM1D software) lines is indicated. The interpretation after the inversion of the data is also included; red straight lines represent the resistivity as a function of the depth. In charts A and B, data collected at station S1 (see Figure 2) using two different devices are represented: central-loop (A) and single-loop (B). In charts (C,D), data collected at station S10 (see Figure 2) using central-loop (C) and single-loop (D) devices are also shown.

Figure 6. Representation of the induced voltage curves (blue line) and apparent resistivity (red line) as a function of time at stations S2 (A), S3 (B), S4 (C), S5 (D), S6 (E), S7 (F), S8 (G), and S9 (H) using a central-loop device. For each of them, the fit between the field line and the line of the ZONDTEM1D software model has been represented. The straight red lines represent the resistivity as a function of the depth, interpreted after data inversion. A good fit is deduced, with RMS error value ranges between 1.7% and 5.6%.

Figure 7. 2D section showing resistivity distribution as a result of smooth-model inversion implemented in the program ZondTEM2D (A). Seismic reflection profile and TDEM section correlation (B). Correlation of seismic reflection and electrical resistivity tomography (C).