<figure><img src="image_1.png" /></figure>

Article

# Advancing the Understanding of Complex Piezometric Information: A Methodological Approach Integrating Long-Term Piezometry, Surface Nuclear Magnetic Resonance, and Fracture Analysis Using Insights from the "Calcaires du Barrois" Series, France

Mathieu Bertrand<sup>1,2,*</sup>, Catherine Bertrand<sup>1</sup>, Naomi Mazzilli<sup>3,†</sup>, Sylvain Gigleux<sup>2,†</sup>, Sophie Denimal<sup>1</sup>, Rémi Valois<sup>3,‡</sup>, Lise-Marie Girod<sup>1,2</sup>, Guillaume Cinkus<sup>3</sup>, Valentine Busquet<sup>1,2</sup> and Konstantinos Chalikakis<sup>3,* †‡</sup>

<figure><img src="image_2.png" /></figure>

check for
updates

**Citation:** Bertrand, M.; Bertrand, C.; Mazzilli, N.; Gigleux, S.; Denimal, S.; Valois, R.; Girod, L.-M.; Cinkus, G.; Busquet, V.; Chalikakis, K. Advancing the Understanding of Complex Piezometric Information: A Methodological Approach Integrating Long-Term Piezometry, Surface Nuclear Magnetic Resonance, and Fracture Analysis Using Insights from the "Calcaires du Barrois" Series, France. *Water* **2024**, *16*, 1700. https://doi.org/10.3390/w16121700

Academic Editors: Maurizio Barbieri and Juan José Durán

Received: 8 April 2024

Revised: 4 June 2024

Accepted: 10 June 2024

Published: 14 June 2024

<figure><img src="image_3.png" /></figure>

**Copyright:** © 2024 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

**Abstract:** This study aims to analyze spatio-temporal piezometric data and integrate them with geological, geotechnical, and geophysical data to enhance their interpretation. The research focuses on a site located in the Meuse and Haute-Marne departments of France, which has been under investigation since 1994 as part of the surface facilities characterization for the Meuse-Haute-Marne underground laboratory and the CIGEO (Centre Industriel de stockage Géologique) Andra project. These investigations span different spatial and temporal scales. We observed the contribution of water masses associated with external forcing and identified two types of aquifer responses: a rapid response to rainfall events when fractures are well-connected, and minor reactivity at the matrix level. Additionally, we demonstrated that the matrix compartment can be finely characterized through a combined interpretation of piezometric response analysis, fracture analysis, and surface nuclear magnetic resonance (SNMR) soundings. The methodology developed in this project offers an improved understanding of karst piezometry and/or unsaturated zone extension, which is essential for comprehending flow dynamics and better constraining the functioning of karst aquifers. Furthermore, this site serves as an ideal workshop for studying flow in fractured media, providing valuable insights into hydrodynamic behavior in complex subsurface environments.

**Keywords:** hydrogeology; karst; piezometric analysis; fracturing; SNMR

## 1. Introduction

Karst aquifers generally have significant water resources and reserves, which are still relatively under-utilized and in some regions of the world represent the essential, if not the only, resource [1–3]. As a result, the search for new water resources (water for food, water shortage support, agriculture...) leads to an increasing interest in karst systems, either to optimize seasonal exploitation (active management) or to establish systems to protect quality (protection perimeters and protection works) or quantity (sizing of works in urbanized areas) [1]. Additionally, studying these systems is essential during geotechnical or civil engineering investigations [4].

If karst aquifers are first to be studied using a descriptive and comparative approach [5], then piezometric records are fundamental data for characterizing these systems. They are measured in the field, either at the source or by direct measurement of the water

<footnote n="1">¹ UMR6249 ChronoEnvironnement (CNRS-UFC), University of Franche-Comté, 25000 Besançon, France; catherine.bertrand@univ-fcomte.fr (C.B.); sophie.denimal@univ-fcomte.fr (S.D.); lise-marie.girod@mines-paristech.fr (L.-M.G.); vambusquet@orange.fr (V.B.)</footnote>

<footnote n="2">² Agence Nationale pour la Gestion des Déchets Radioactifs (Andra), RD960, BP9, 55290 Bure, France; sylvain.gigleux@andra.fr</footnote>

<footnote n="3">³ UMR 1114 EMMAH (AU-INRAE), Avignon Université, 84029 Avignon, France; naomi.mazzilli@univ-avignon.fr (N.M.); remi.valois@univ-avignon.fr (R.V.); guillaume.cinkus@univ-avignon.fr (G.C.)</footnote>

<footnote>* Correspondence: mathieu.bertrand@univ-fcomte.fr (M.B.); konstantinos.chalikakis@univ-avignon.fr (K.C.)</footnote>
column in the borehole [6]. It is an essential measure for understanding the behavior of an aquifer [7,8], for characterizing and evaluating its functioning and its hydrodynamic parameters (storage coefficient and hydraulic conductivity) [9–12], and also for detecting interactions between different structures during its operation [13]. These parameters can be extracted by analyzing signals acquired under active stress (usual hydraulic tests) or passive stress (tidal effect, [14]). These data are essential for the exploitation of water resources [15,16] as well as for civil engineering investigations [17,18].

However, in complex systems such as karst aquifers, the structure and the functioning of aquifers are difficult to understand due to the high heterogeneity degree and the different types of porosity that are involved in groundwater circulation. Indeed, the piezometric levels measured in boreholes are influenced by different types of water volumes coming from the matrix, fractures, or conduits [19–21]. Boreholes that are only a few meters apart can have piezometric levels that differ significantly, in part due to different dynamics between matrix and fractures and/or karst conduits. These differences in piezometric levels may also be due to varying degrees of saturation within the aquifer. Water circulation in this type of environment is thus constrained by hydrodynamic heterogeneities, which induce preferential circulation and, therefore, non-uniform water (and pressure) transfers [22,23], with a continuum between very fast, non-dispersive flows (piston or conduit flows) and very slow, highly dispersive flows (diffuse flows) [24–26]. Borehole logs, wall imaging techniques, and core analysis provide high-resolution fracturing information. Combined with hydraulic tests, they provide relevant information on the hydrodynamic functioning of the aquifer, such as fracturing, karstification degree, and connectivity. Nonetheless, the obtained information only reflects the well spatial scale, so there is a high uncertainty far from the wells. In these subsurface geological environments, flows are generally not linked to a particular structure or homogeneous reservoir but are distributed at all scales [27]. This complexity poses the problem of detecting, measuring (and being representative of), and monitoring these flows in 4D [19,28].

Non-intrusive techniques such as geophysical methods can complement and spatialize the geological and hydrogeological knowledge gained by analyzing the water levels in boreholes [29]. These methods are based on analyzing the propagation of physical variables (e.g., electromagnetic field, acoustic waves). This propagation is influenced by various factors, including the presence or absence of water in the environment being studied [30]. Among these methods, surface nuclear magnetic resonance (SNMR) provides a signal directly linked to the presence of groundwater. The initial amplitude of the SNMR signal is directly related to the presence of mobile water in the subsurface, while its relaxation time can be related to pore size [31,32]. This method thus provides direct information on groundwater hydrodynamic properties in karst environments [33].

When industrial activities are spread out above or within fractured and/or karstified aquifers, it is important to develop a comprehensive hydrodynamic functioning model to ensure the durability of the resources and minimize related risks. As an example, the French National Agency for the Management of Radioactive Waste (Andra) aims to securely dispose of high and low-intermediate level long-lived radioactive waste at the CIGEO (Centre Industriel de stockage Géologique) project (Meuse, France). The work support zone, and the ramp zone facilities, will be used for digging and building the underground structures located at a depth of approximately 500 m and for receiving and inspecting waste packages and preparing them for transfer to the underground facilities via a ramp these zones measuring around 550 hectares lies on the "Calcaires du Barrois" series. The basal part of the "Calcaires du Barrois" is the dominant formation at the outcrop, consisting of a complex hydrosystem with two horizons of differing permeability. The upper horizon is fractured, resulting from the decompression and cryoturbation phenomena of the last deglaciations, while the underlying horizon exhibits diffuse karstification. However, this karstification is not concentrated enough to produce large-scale dissolution features such as caves or metric drains. To properly size surface facilities and reservoir access shafts, a better understanding of the flow dynamics within
this formation is essential. Additionally, the hydrogeological characterization of the
“Calcaires du Barrois” is a key part of the project, aimed at establishing an effective
environmental monitoring network.

The hydrodynamic behavior of the “Calcaires du Barrois” is directly linked to the fracturing of the environment. Numerous geological and geotechnical investigations have been conducted over hundreds of boreholes at this site over the years, providing data on fracturing, and petrophysics. Piezometry and hydrodynamics investigation has been conducted according to dozens of piezometers. Nineteen boreholes were re-used as piezometers. These were selected from an area of 13 km² as they have provided good quality data series for more than 8 years. Because this density of investigation is not often available for studying karst systems [29], we designed this project to use these facilities to study the hydrogeologic complexities of this “fractured-rock karst site”.

We decided to bring together the data collected at the study site and evaluate the data
as a whole, complemented by SNMR. This investigation is directly related to the presence
of mobile water in the subsurface and provides direct information on the hydrodynamic
properties of mobile groundwater in karst environments. It is particularly suited to the
study area. We thus combined (1) the analysis of long piezometric records with (2) full
fracture analysis from borehole logs and (3) geophysical measurements to gain further
insight into the complex “Calcaires du Barrois” system. This approach is used to shed
more light on hydrodynamic responses and their origins in the complex “Calcaires du
Barrois” system. We also provide methodological insights that are easily transposable to
other studies aiming at water body characterization in karst environments.

## 2. Context

### 2.1. *Study Site*

The Meuse/Haute-Marne region is located on the south-eastern edge of the Paris Basin, between the valleys of the Marne and the Meuse rivers (France, Figure 1A) and, more precisely, in the outcrop zone of the Upper Jurassic formations ("Calcaires du Barrois" series), at the boundary of the presence of Cretaceous terrains. The "Calcaires du Barrois" series is a succession of limestone strata (Sublithographic limestone, Dommartin limestone, "Oolithe de Bure" formation, and Decayed limestones) and a marl layer that is approximately 10 m thick ("Pierre Chaline" formation), separating the Sublithographic limestone and the Dommartin limestone. The Kimmeridgian marls form the bedrock on which the sedimentary sequence rests (Figure 1A,B,D). This series dips northwest towards the center of the Paris Basin with an average dip of about 1 degree [34]. This shallow dip is associated with a plateau morphology: the "Calcaires du Barrois" series plateau. This plateau is cut by relatively narrow, incised valleys. It is bordered to the east by the Ormançon River and to the west by the Saulx River. Between these two rivers flow the Bureau and Orge rivers, tributaries of the Saulx.

The study site covers 13 km² and is located at the head of the Orge catchment. It is bounded by the Orge River to the west and the Bureau River to the east (Figure 1C). These two rivers are characterized by dry periods from April to October. The Bureau joins the Orge near Bure (Figure 1A). In this sector, only the Sublithographic limestone is exposed. In detail, this limestone is composed of five levels, which are distinguished by variations in facies (multi-centimetric limestone-mudstone alternation at the bottom of the series) and thickness (up to 15 m for level N5) (Figure 1B).

Between 2002 and 2023, rainfall was measured at the CMHM (Centre Meuse-Haute-Marne) site with a weather station (hourly and daily). The average annual rainfall (over 21 years) recorded at the CMHM was 946 mm. The annual precipitation at CMHM ranged between 750 mm (in 2004) and 1200 mm (in 2007).
<figure><img src="image_4.png" /><figcaption>**Figure 1.** (A)—Geological map of the study area. The white dashes represent the location of the future waste reception, control, and preparation area of the CIGEO site area called Descendrie area. (B)—Lithostratigraphic succession of the study area and location of the aquifers in this succession. (C)—Focus on the study area, showing the location of piezometers, geological boreholes, and SNMR soundings. The numbers correspond to the piezometer names in Appendix A. (D)—North-West (A)/South-East (B) cross-section with vertical exaggeration (10×).</figcaption></figure>
## 2.2. Geological and Hydrogeological Characteristics of the "Calcaires du Barrois" Series

The local hydrogeological system of the "Calcaires du Barrois" series consists of three more-or-less well-connected aquifer levels where the series is complete (Figure 1B). These are, from the bottom to the top [34]:

* The Sublithographic limestones (25 m thick) are lithographic micritic limestones with rare bioclastic beds and rare marly interbeds of decimetric thickness.

* The Dommartin limestones (63 m thick) consist of about 25 m of micritic limestones with bioclastic beds and marl interbeds at the base, an interval of about 25 m of micritic limestones with very rare bioclasts and bioclastic packstones at the top of this interval, and about 13 m of micritic limestones with very little bioclasticity at the top of the Dommartin limestones.

* The Decayed limestones (36 m—thick), which correspond to micritic limestones with little bioclasticity.

These levels are separated by aquitards such as the "Pierre Chaline" formation, a predominantly marly level, and the "Oolithe de Bure" formation, formed by a massive bed of beige limestone with oolites and a sparitic matrix. Depending on the degree of fracturing, these aquitards can be relatively permeable. The Kimmeridgian marl with exogyres that underlie the "Calcaires du Barrois" series forms an impermeable layer. The upper and lower limestone layers contained in these marls are also aquifers to a lesser extent. The aquifers in this series supply numerous springs with low flow rates [34].

The water table in the Sublithographic limestones is unconfined in the upstream part of the catchment but becomes captive north of the town of Bure under the "Pierre Chaline" formation. The water table in the Dommartin limestones is semi-captive due to the limited extent of the Oolithe de Bure. Over the extent of the study area, the water table in the Decayed limestones is unconfined.

Karstification is well developed within the "Calcaires du Barrois" series. In addition, typical karst features can be observed at the surface (sinkholes, wells). Speleological explorations to the north-west of the area [35] have highlighted the extensive development of karstification in the "Calcaires du Barrois" series. In the study area, numerous losses and resurgences of the Orge and Ormançon have been observed where these rivers cross the "Pierre Chaline" marls [36].

# 3. Data and Methods

The proposed methodology is based on spatial analysis of piezometric and fracture records, coupled with imaging and monitoring of water content using the SNMR method.

## 3.1. Hydrogeological Investigations: Piezometry Analyses

The extensive instrumentation at the site provides a considerable amount of borehole data. Our first objective was to extract from these records statistical or correlative indicators of the hydrodynamic behavior of the structures, and then to use these indicators to establish groups of similar hydrodynamic functioning.

Seven piezometric time series indicators were defined in this project:

1. The variance, which is measured via time series variation coefficient, an indicator of the dispersion of the distribution of piezometric levels.

2. The maximum amplitude (difference between the "Calcaires du Barrois" series base and the water level extremum piezometric levels). The base of the "Calcaires du Barrois" series is used as a reference to compare the height of piezometric variations measured in the boreholes.

3. The slopes of piezometric rises and recession. The positive and negative slopes of the piezometric curves provide information on the rate of piezometric rise and the rate at which wells are drying up. These indicators are determined by studying the derivative of piezometric variations as a function of time. These derivatives are calculated at all points of a hydrograph using the finite difference method. The indicators provide
information on the hydrodynamic parameters of the system. They therefore provide
an insight into the structure of the system, such as the relative fracture density in the
vicinity of the studied piezometers.

4. System inertia or memory effect. This parameter quantifies the time taken for an event to affect the chronicle or, in other words, the time *n* taken for an event occurring at time *t* to not affect the event recorded at time *t* + *n*. The inertia of the system is expressed by the autocorrelation function [37]. This function is a measure of the degree of dependence on events in the same chronicle by assessing their repetitive nature. A significant memory effect reflects a certain degree of flow regulation by the system and, therefore, the importance of a capacitive character [38,39]. The memory effect is an indicator of the structure of the system and fracturing.

5. Response time and cross-correlation coefficient. These indicators examine the strength of the relationship between the input signal, rainfall, and what can be considered as an approximation of the output signal, piezometry. It is assumed that there is a causal relationship between these two variables. The cross-correlation function corresponds to the transfer function between the two signals and provides information about the structure of the system.

These indicators were grouped by hydrodynamic behavior using principal component analysis and hierarchical ascending classification (HAC).

This approach was applied to 19 piezometers spread over the Descenderie area and over two time periods (Appendix B). Four of the 19 piezometers have a continuous series of over 25 years (piezometers 16, 17, 18, and 19) and the other piezometers have a continuous series since 2016. Datasets from recordings on the 19 piezometers between November 2016 and July 2021 were selected for use in this project because they provide the longest records to harmonize the processing, and they have been homogenized at hourly time steps. Before statistical processing, an analysis of the raw data was performed to quantify the missing data for all records (water samples or well tests). Over the 8 years of recording, the missing data represent 6.5% of the total data collected from the 19 piezometers.

Data on piezometric variations (mNGF) have been continuously recorded at quarter-hourly or hourly intervals by pressure sensors (Diver, Ven Essen Instrument) (Appendix B). Ranges are 10 m for pressure with an accuracy of 0.5 cm H₂O, and between -20 °C and 80 °C with an accuracy of +/-0.1 for temperature.

3.2. *Geological Investigations: Fracturing Analysis*

The fracturing and diffuse karstification degrees of the "Calcaires du Barrois" series have been described using borehole logging coupled with core analysis when possible. One of the objectives was to estimate the effect of lithology on the density and vertical extent of fracturing and karstification along the marl–limestone sequence where the flows seem to occur.

We used the following logs:

1. Optical Televiewer borehole wall imaging (sensor OPTV), which provides a complete, continuous scan of the borehole and, thus, provides a developed, oriented image of the borehole walls that can be used for subsequent image processing allowing measure dips and identify and orient planes of stratification, fracturing, and schistosity, (shown in polar plots). The characteristics of the Optical Televiewer are 360° camera 2 m long and 50 mm in diameter with 1 mm resolution. The quality of the images obtained depends mainly on the acquisition conditions in the borehole: the turbidity of the water, the diameter of the borehole, and the stability of the ground.

2. Natural radioactivity logging or “Gamma Ray” (GR) measures the natural gamma radioactivity of rocks. It measures the content of naturally radioactive elements in the rocks and thus clarifies and correlates the lithological cross-section obtained by sampling the ground during drilling. In the study area, it can be used to distinguish between marl intervals (high values) and limestone intervals (low values). Correlations between boreholes are based on this measurement.
We used 54 boreholes that had been drilled in the study area; they traversed both the Sublithographic limestones and the “Pierre Chaline” formation, depending on their location (Figure 1C). Forty-eight (48) of the 54 boreholes were cored and investigated using two types of logging: optical and gamma-ray. Thirty were logged optically and by gamma-ray; one was optically logged and 17 were investigated by gamma-ray logging. Of the 6 boreholes drilled destructively, only the gamma-ray logs could be used. Drilling and logging were conducted in 2015 and 2016. In most cases, wall imaging in boreholes is not relevant below the water level due to issues with turbidity.

### 3.3. *Geotechnics Investigations: Hydraulic Tests*

Two types of hydraulic tests were carried out: (i) Lefranc tests and (ii) pumping tests, carried out after the groundwater table had stabilized.

#### 3.3.1. Lefranc Tests

The Lefranc water test is an in situ test used to determine the value of the local permeability coefficient (KL), which may differ significantly from the bulk permeability coefficient. The Lefranc water test applies to all fine or gravelly soils below the water table, where the assumed permeability coefficient is greater than approximately $10^{-6}$ m/s. It can be conducted in exploratory boreholes as work progresses.

Lefranc tests are conducted either under constant load (water is injected into the borehole at a constant rate until the level in the borehole stabilizes) or under variable load (a volume of water is injected into the borehole and the drop in level is monitored as a function of time). These tests are carried out in 1 m increments as the borehole progresses.

One month after drilling completion, two constant flow injection tests were conducted throughout the borehole ({28}, {29}; Figure 1C). For {28}, the tests were carried out by injecting water at a constant rate of 0.02 l/s for 2 h and then, after a recovery time of 21 h, an injection of 0.03 L/s for 7 h. The protocol used for {29} was more or less the same in terms of injection time and recovery time. However, the injected flow rates are doubled, i.e., 0.04 L/s and then 0.06 L/s. Pressure measurements in the borehole with a DIVER pressure probe were continued for 450 h after the end of the tests. The different flow rates used for these tests were chosen by the companies conducting the geotechnical investigations.

Results of the constant load tests were based on the expression of [40], which expresses the permeability coefficient as a function of flow rate, a factor that depends only on the borehole shape and its position with respect to the aquifer boundaries, the water column height, and the borehole diameter. Only measurements recorded after stabilization were used to determine the permeability coefficient.

Twenty-six hydraulic tests (Lefranc tests) were conducted [41] by a company specialized in geotechnical studies during the destructive drilling of 4 boreholes (Figure 1C). A total of 16 slug tests were performed by injecting water at a constant rate of 1 L/s for 1 h into three boreholes {25}, {26} and {27}. Tests that did not produce significant measurements (very low infiltration) were continued for one hour with a variable load. This consisted of 10 tests, including 5 in borehole {24}, 1 in borehole {25}, 1 in borehole {27}, and 3 in borehole {26}.

#### 3.3.2. Pumping Tests

The constant flow test, which is the most commonly used pumping test, consists of lowering the piezometric surface of the water table by pumping and measuring changes in the level of this surface, as a function of time, as well as the flow rate pumped. This allows the permeability coefficient of the tested layer and the storage factor to be determined.

The pumping was carried out in borehole {20} (source borehole) and the changes in the piezometric surface over time were monitored using 3 piezometers (observation boreholes) located 3 m {21}, 8 m {22}, and 20 m {23} from the borehole (Figure 1C).

The pump test was carried out in December 2015 according to the NF P94-130 standard [42]. The lowering of the groundwater level was measured using pressure sensors in the 3 control piezometers. The pumping test was stopped after 1 h to prevent the pump
from draining. The piezometric level was measured in the source borehole being pumped and in the 3 observation piezometers during a recovery period of 21.5 h after pumping was stopped. The interpretations of the pumping test are derived from the analytical solution developed by Moench [43] for the source borehole and from the modified Theis [44] solution for one of the observation piezometers. The analytical solution developed by Moench [43] takes into account the incomplete borehole effect and the parietal effect in a homogeneous and anisotropic-free aquifer with a drip. The analytical solution developed by Theis [44] was originally developed for use in confined aquifers, but can also be used for unconfined aquifers by correcting the drawdown data [45].

The permeability coefficient k obtained from the pumping test is the horizontal permeability coefficient of the layer. It is a global or "large" coefficient, representative of the average hydraulic behavior of the volume of the layer under test. The storage factor or S coefficient is dimensionless. It gives the change in the volume of water contained in a vertical prism with a horizontal section of unit area for a unit change in the hydraulic load.

### 3.4. *Hydrogeophysics Investigations: Surface Nuclear Magnetic Resonance (SNMR)*

We used SNMR to analyze the spatial and temporal variation of water content to estimate the location of water volumes passing through these fractured environments [30,32]. The measuring device was NUMISPoly equipment from IRIS Instruments, with a coincident transmitting/receiving eight-square loop configuration composed of two squares with 40 m sides.

SNMR signal parameters (initial amplitude $E_0$ and observed transverse relaxation time $T2^*$) are derived from the SNMR signal envelope assuming a mono exponential decay [33]. The SNMR water content (mobilizable water content) was derived from $E_0$. The effective porosity in which the water is stored has been derived from the apparent transverse relaxation time $T2^*$ [32]. The hydraulic parameters of the aquifers are then estimated qualitatively: The water content gives an image of the drainage porosity of the reservoir, and the $T2^*$ constant provides information on the characteristics of water-filled pores [46]. Measuring $T2^*$ is easier and faster than measuring $T1$ because the $T1$ assessment requires the use of a two-pulse protocol with a variable delay between pulses [47]. In the absence of magnetic heterogeneities, the $T2^*$ estimate is considered reliable [47], so we used $T2^*$ in our survey.

Several procedures have been used to improve the quality of SNMR data in this study [48]. To compensate for the fact that electromagnetic (EM) noise has a large influence on the measurement, an EM noise measurement campaign was conducted prior to drilling to select the areas with the least EM noise. An eight-square loop shape was favored to further reduce the EM noise effect on the soundings. The excitation pulse was repeated and the recording sequences were averaged to reduce the influence of random noise.

SNMR signal filtering is indispensable to improve the signal-to-noise ratio and to allow the estimation of SNMR parameters. The filtering procedure was chosen based on tests carried out using the raw noise data on a synthetic signal and performed using NumPro 12 × 2 from the SAMOVAR suite. A notch filter of the type global Narrow was used to attenuate the disturbance field due to the 50 Hz harmonics of the power lines present near the sounding [49]. The maximum noise level for the use of a stack was set at 500 nV, and the maximum average noise level was set to 50 nV. Damage to the signal caused by filtering was calculated using the RMSE (Root Mean Square Error). The RMSE compares the square root of differences between predicted and observed values.

SNMR data quality was estimated based on the average signal-to-noise ratio ($S/N$) [50] derived from stacked and filtered measurements [51]. The higher the $S/N$, the higher the reliability of the SNMR survey. We consider SNMR soundings to be of acceptable quality when the $S/N$ is ≥ 2.

A total of 19 SNMR soundings were implemented in 3 acquisition zones, under different hydraulic regimes. Acquisition zones covered the entire sector with varying thicknesses of Sublithographic limestone (3 m, 12 m, and 25 m). In each zone, a single low-water campaign (summer 2021) and a high-water campaign (spring 2022) were conducted
to characterize the seasonal water level variations. To better assess SNMR water content dynamics, additional measurements were performed at site 38 prior and after rainfall events during both low water (7 soundings over 38 days from 23 June 2021 to 30 July 2021) and high water (9 soundings over 8 days from 5 April 2022 to 12 April 2022).

Only soundings with an average number of recording sequences (stacks) greater than 200 and a signal-to-noise ratio greater than 2.5 are presented in this paper (Table 1).

<table><caption>**Table 1.** Table of SNMR surveys and their main characteristics (name, site, date, average number of stacks, signal-to-noise (S/N) ratio, maximum depth with good resolution, rainfall in the 7 days preceding the sounding (Rainfall D-7), and the hydraulic regime.</caption><thead><tr><th>Sounding Name</th><th>Site</th><th>Date</th><th>Average Number of Stacks</th><th>Ratio S/N</th><th>Maximum Depth Resolution (m)</th><th>Rainfall D-7 (mm)</th><th>Hydraulic Regime High-Water (HW) or Low-Water (LW)</th></tr></thead><tbody><tr><td>{37}</td><td>{37}</td><td>1 March 2022</td><td>299</td><td>2.51</td><td>24.41</td><td>7.8</td><td>HW</td></tr><tr><td>{38}_1</td><td>{38}</td><td>23 June 2021</td><td>274</td><td>2.67</td><td>25.77</td><td>10.4</td><td>LW</td></tr><tr><td>{38}_2</td><td>{38}</td><td>28 June 2021</td><td>288</td><td>2.7</td><td>25.77</td><td>21.0</td><td>LW</td></tr><tr><td>{38}_3</td><td>{38}</td><td>11 April 2022</td><td>367</td><td>9.18</td><td>21.6</td><td>48.4</td><td>HW</td></tr><tr><td>{38}_4</td><td>{38}</td><td>12 April 2022</td><td>243</td><td>6.63</td><td>21.6</td><td>48.4</td><td>HW</td></tr><tr><td>{39}</td><td>{39}</td><td>3 March 2022</td><td>289</td><td>7.78</td><td>24.41</td><td>4.4</td><td>HW</td></tr></tbody></table>

A smooth inversion was performed with SVD (singular value decomposition) type regularization. We chose a relatively low regularization (15 for water content and 20 for T2*) based on a preconception of the geology of the soundings. A sharp contrast in hydrodynamic properties is indeed expected between the upper part, the limestone part and the marl, and the lower part of the sounding. The uncertainty in the water content is reflected in the 95% confidence interval of the SAMOVAR plots.

## 4. Results

### 4.1. Piezometric Data Processing (2016–2021)

Table 2 shows the values of the piezometric indicators determined for each borehole, together with the percentage of missing data for that borehole. These results will be used as the basis for the multi-criteria analysis that will characterize the comparative phase of this study.

<table><caption>**Table 2.** Table of values of piezometric indicators for each structure.</caption><thead><tr><th>Name</th><th>1-Variance</th><th>2-Amplitude (m)</th><th>3-Flood Rise Gradient</th><th>3-Recovery Gradient</th><th>4-Inertia (Days)</th><th>5-Response Time (Hours)</th><th>5-Rainfall-Piezometry Correlation</th><th>Missing Data (%)</th></tr></thead><tbody><tr><td>{1}</td><td>1.71</td><td>8.7</td><td>1.72 × 10<sup>-6</sup></td><td>-1.25 × 10<sup>-7</sup></td><td>69</td><td>25</td><td>0.14</td><td>1.59</td></tr><tr><td>{2}</td><td>1.14</td><td>8.88</td><td>1.13 × 10<sup>-6</sup></td><td>-1.24 × 10<sup>-7</sup></td><td>61</td><td>16</td><td>0.24</td><td>3.23</td></tr><tr><td>{3}</td><td>1.52</td><td>7.02</td><td>6.61 × 10<sup>-6</sup></td><td>-1.26 × 10<sup>-7</sup></td><td>62</td><td>37</td><td>0.18</td><td>7.05</td></tr><tr><td>{4}</td><td>1.7</td><td>15.93</td><td>7.88 × 10<sup>-7</sup></td><td>-1.91 × 10<sup>-7</sup></td><td>115</td><td>30</td><td>0.03</td><td>1.17</td></tr><tr><td>{5}</td><td>0.65</td><td>8.82</td><td>3.31 × 10<sup>-7</sup></td><td>-9.18 × 10<sup>-8</sup></td><td>40</td><td>18</td><td>0.23</td><td>0.02</td></tr><tr><td>{6}</td><td>0.82</td><td>6.71</td><td>5.57 × 10<sup>-7</sup></td><td>-1.23 × 10<sup>-7</sup></td><td>56</td><td>199</td><td>0.09</td><td>0.35</td></tr><tr><td>{7}</td><td>2.08</td><td>18.87</td><td>1.84 × 10<sup>-6</sup></td><td>-2.98 × 10<sup>-7</sup></td><td>74</td><td>16</td><td>0.17</td><td>3.77</td></tr><tr><td>{8}</td><td>0.45</td><td>8.11</td><td>9.26 × 10<sup>-8</sup></td><td>-3.98 × 10<sup>-9</sup></td><td>2</td><td>15</td><td>0.24</td><td>0</td></tr><tr><td>{9}</td><td>0.92</td><td>10.13</td><td>3.41 × 10<sup>-7</sup></td><td>-1.10 × 10<sup>-7</sup></td><td>54</td><td>18</td><td>0.2</td><td>1.59</td></tr><tr><td>{10}</td><td>0.88</td><td>12.89</td><td>1.49 × 10<sup>-7</sup></td><td>-3.35 × 10<sup>-8</sup></td><td>30</td><td>17</td><td>0.21</td><td>0</td></tr></tbody></table>
<table><caption>Table 2. *Cont.*</caption><thead><tr><th>Name</th><th>1-Variance</th><th>2-Amplitude (m)</th><th>3-Flood Rise Gradient</th><th>3-Recovery Gradient</th><th>4-Inertia (Days)</th><th>5-Response Time (Hours)</th><th>5-Rainfall-Piezometry Correlation</th><th>Missing Data (%)</th></tr></thead><tbody><tr><td>{11}</td><td>0.23</td><td>3.71</td><td>9.99 × 10<sup>-8</sup></td><td>-6.81 × 10<sup>-9</sup></td><td>8</td><td>4</td><td>0.36</td><td>0</td></tr><tr><td>{12}</td><td>0.4</td><td>5.11</td><td>7.88 × 10<sup>-8</sup></td><td>-1.32 × 10<sup>-9</sup></td><td>7</td><td>13</td><td>0.21</td><td>0.02</td></tr><tr><td>{13}</td><td>1.4</td><td>10.86</td><td>2.31 × 10<sup>-6</sup></td><td>-1.04 × 10<sup>-7</sup></td><td>74</td><td>16</td><td>0.18</td><td>0.28</td></tr><tr><td>{14}</td><td>0.64</td><td>4.83</td><td>4.93 × 10<sup>-9</sup></td><td>-5.02 × 10<sup>-10</sup></td><td>42</td><td>20</td><td>0.23</td><td>0.06</td></tr><tr><td>{15}</td><td>0.55</td><td>2.44</td><td>3.95 × 10<sup>-7</sup></td><td>-3.12 × 10<sup>-8</sup></td><td>68</td><td>176</td><td>0.1</td><td>0</td></tr><tr><td>{16}</td><td>0.82</td><td>10.91</td><td>6.49 × 10<sup>-7</sup></td><td>-1.26 × 10<sup>-8</sup></td><td>6</td><td>13</td><td>0.29</td><td>0.76</td></tr><tr><td>{17}</td><td>0.99</td><td>5.28</td><td>1.43 × 10<sup>-6</sup></td><td>-1.31 × 10<sup>-7</sup></td><td>65</td><td>100</td><td>0.16</td><td>2.28</td></tr><tr><td>{18}</td><td>3.2</td><td>11.66</td><td>3.79 × 10<sup>-6</sup></td><td>-3.63 × 10<sup>-7</sup></td><td>71</td><td>68</td><td>0.11</td><td>70.61</td></tr><tr><td>{19}</td><td>0.99</td><td>5.16</td><td>4.87 × 10<sup>-6</sup></td><td>-5.66 × 10<sup>-8</sup></td><td>60</td><td>39</td><td>0.19</td><td>29.38</td></tr></tbody></table>

### 4.1.1. Monocriteria Statistical Analysis

#### Strict Signal Indicators

For the 19 piezometers studied during the period 2016–2021, the coefficient of variation ranges from 0.23 to 3.2. This indicator is difficult to interpret in isolation from the other indicators. However, it is relevant to describe the variability of the signal. The maximum amplitude between the lowest and highest levels in each series is highly variable. It ranges from 2.4 to 18.9 m at piezometers {15} and {7}, which are about 750 m apart. The average maximum amplitude is 8.73 m with a dispersion of values of 4.2. These indicators reflect the heterogeneity of piezometric behavior in the structures.

#### Hydrodynamic Parameter Indicators

Piezometer {14} has the lowest rates of piezometric rise and recovery with values of $4.93 \times 10^{-9}$ and $-5.02 \times 10^{-10}$. The highest velocities are observed for piezometer {3}, a piezometric rise with a value of $6.6 \times 10^{-6}$, and structure {18}, for recovery with a value of $-3.63 \times 10^{-7}$. The slope values give a relative indication of the rate at which changes occur once a response begins after a rainfall event. The variation in these values supports the hypothesis of a heterogeneous environment with quite different behavior from one structure to another. The variation gives an idea of the hydrodynamic characteristics of the environment around the structures.

#### System Structure and Fracture Indicators

The inertia values derived from interpretation of the autocorrelation functions (Table 2) again highlight a divergence in the memory effect of the system in the vicinity of the investigated structures. The values range from 2 to 115 days, 400 m apart. In general, the interpretation of the autocorrelation curves shows that the memory effect can be quite different between one structure and another, despite the short distance between them. This is the case for holes {16} and {6}, 80 m apart, with memory effects of 6 and 69 days, respectively. These results are considered to be statistically robust.

The phase shift values between input and output signals, derived from the interpretation of the rain–piezometry cross-correlation curves, reflect variable response times ranging from 4 h {11} to 8.5 days {6} (Table 2).

The associated correlation coefficients are also variable, ranging from 0.03 {4} to 0.36 {11}. In absolute terms, these “rainfall–piezometric” correlation values appear to be quite low, reflecting the inertia of the system.
### 4.1.2. Multi-Criteria Statistical Analysis

Multi-criteria analyses have the advantage of being able to characterize how the indicators work together and to determine which indicators are likely to discriminate between the piezometric signals from different boreholes (Principal Component Analyses (PCA)). Hierarchical Ascending Classifications classify the piezometers based on the coordinates of the points in the PCA.

The PCA was performed on all seven indicators for the 19 piezometers located in the study area. On the F1F2 axis, the percentage of information explained was 75.48%.

The distribution of the indicators in the variable space shows that all the indicators are significant on the F1F2 axes (due to their proximity to the circle of radius 1), except the slope of the flood rise (Figure 2C).

The amplitude of the variations correlates with the variance and the slopes of the piezometric rises in the series and is anti-correlated with the recovery slopes (Figure 2C). The inertia is anti-correlated with the rainfall-piezometric parameter. Response time shows no correlation with other indicators.

According to the distribution of structures in the space of individuals (Figure 2B,D), points {6} and {15} are strongly expressed on the F2 axis and seem to be discriminated by their response time.

Points {4}, {7}, and {18}, which are strongly expressed on the F1 axis, can be explained by their dispersion and variation amplitude values, which vary inversely to recovery slopes, which are also strongly expressed in the variation of these structures.

Points {8}, {11}, {12}, {14}, and {16} are expressed on the F1 axis and are strongly explained by their "rainfall-piezometry" correlation value.

On the other hand, the remaining points {1}, {2}, {3}, {5}, {9}, {10}, {13}, {17}, and {19} form a median group that is moderately expressed on both axes (Figure 2B,D).

A HAC was then used to discriminate between groups of piezometers based on the coordinates of the points in the PCA. Five groups stood out on the dissimilarity dendrogram. The HAC was performed for a significance threshold of 0.08, symbolized by a solid line on the dendrogram. At the end of this analysis, five groups stand out (Figure 2A).

The first group (A1, pink) discriminates between points {4}, {7} and {18}. This can be explained by the F1 axis, which is mainly responsible for the "rainfall-piezometry" correlation, amplitude, and dispersion. The representative piezometer in this group is {18}. The group is representative of highly transmissive behavior.

A second group (A2, red), which essentially shows "rainfall-piezometry" correlation values expressed primarily on the F1 axis, consists of points {8}, {11}, {12}, {14}, and {16}. The representative piezometer in this group is {16}. This group is related to the fracturing of the environment and therefore has a transmissive behavior.

A third group (A3, green), representing the response time indicator, is strongly expressed on the F2 axis and consists of points {6} and {15}. The representative piezometer in this group is {15}, which is typical of a capacitive environment.

A fourth group (A4, purple), the median, differs from the first by the best flood gradients, an indicator that is not very representative on the F1F2 axes and consists of points {3} and {19}. The most representative piezometer in this group is {19}. This group presents both transmissive and capacitive behavior.

A fifth and final group (A5, black), median, which has no axis and is therefore a factor of influence, consists of points {1}, {2}, {5}, {9}, {10}, {13}, and {17}. The representative piezometer in this group is {17}. This group has an above transmissive and capacitive behavior.
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Result of multi-criteria statistical analysis. (A)—Mapping of piezometer positions according to the groups resulting from the hierarchical classification. The colors are associated with the groups determined by the hierarchical classification. (B)—Hierarchical classification of piezometers according to the selected criteria described in methodology. (C)—PCA variable graph. (D)—PCA individuals' graph.</figcaption></figure>
## 4.2. Processing Geological and Geotechnical Data

### 4.2.1. Fracturing Data

Gamma Ray logs are used as a reference for all geological and fracture studies in the study area and are used to correlate Sublithographic limestone levels from one borehole to another. Areal fracturing has led to the development of a karst network that is visible on optical logs, as shown in Figure 3C,E,F).

<figure><img src="image_6.png" /><figcaption>**Figure 3.** (A)—Weathered zone. (B)—Fracturing in a piezometer of group A1 (poorly developed fracturing). (C)—Fracturing in a borehole near a piezometer of group A2 (16), with highly developed fracturing. (D)—Fracturing in a borehole near a Group A3 piezometer (29), with vertical and horizontal fracturing visible. (E)—Dissolution marks in a borehole distant from the piezometers and therefore not able to be connected to a functioning group. (F)—Highly developed and partly karstified fracturing (presence of conduit) in a borehole that cannot be connected to a piezometric functioning group.</figcaption></figure>

Below and sometimes at some distance from the altered surface zone (Figure 3A), centimeter-scale horizontal drains are most commonly found in the limestone intervals of the most marly levels, following the bedding planes (Figure 3B,D). These drains, which show clear traces of dissolution, are often found in the limestone interbeds of N4 and N2 (Figure 1B) or in the upper part of the Exogyres marls (Figure 3B,E,F).

Fracture orientation detected by wall imagery from 25 boreholes (3 of which are to the north of the study area) is summarized in Figure 4. A primary direction centered on N150°E is visible, with a scatter of measurements around this value. In detail (Figure 4), perpendicular directions dominate in some of the boreholes.
<figure><img src="image_7.png" /><figcaption>**Figure 4.** Map of the distribution of oriented fracturing based on optical logs, with differentiation between open fracturing (in red) and closed fracturing (in blue). The synthesis of all rose diagrams is represented by the large rose.</figcaption></figure>
The fractures identified by the methods used are local and difficult to match from one borehole to another.

### 4.2.2. Hydraulics Tests

#### Lefranc Tests: Permeability on the Scale of the Borehole

Results of the tests conducted on boreholes {24}, {25}, {26}, and {27} are shown in Figure 5.

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Permeability tests (Lefranc) at various normalized depths in relation to the natural terrain (NT). Legend colors correspond to those shown in Figure 1B.</figcaption></figure>

Stratigraphically, the tests were conducted in the "Pierre Chaline" formation to a depth of 4.5 m, in the upper part of the Sublithographic limestone (N5) at depths from 4.5 m to 13.5 m, in the middle part of the Sublithographic limestone (N4) at depths from 17 to 18.5 m, and in the lower part of the Sublithographic limestone (N3 to N1) at depths from 22.5 m to 30 m. Kimmeridgian marls were also encountered in a single borehole {28}.

Permeability results from constant load tests range from $10^{-4}$ to $10^{-6}$ m/s, with permeability decreasing with depth. Permeability results from variable load tests range from $10^{-8}$ to $10^{-10}$ m/s, regardless of the formation and test depth.

The formations are permeable to semi-permeable ($10^{-5}$ to $10^{-6}$ m/s) down to a depth of 17 m (N5). Then, in the middle to lower part of the Sublithographic limestone (N4 to N1), they reach values of $10^{-8}$ to $10^{-10}$ m/s. A more permeable zone at $10^{-5}$ m/s in the lower part of the limestone (N3) was encountered in {25}.

The permeability differences measured reflect the general anisotropy of the various lithologies of the Sublithographic limestones. However, they also reflect anisotropy on a larger scale. For example, borehole {24} crosses the upper and middle parts of the Sublithographic limestones (N5 and N4), which have low permeabilities in the order of $10^{-8}$ to $10^{-10}$ m/s, whereas it is located 800 m from the other boreholes, which have higher permeability ranges (Figure 5). The borehole {24} is also very close to Group A2 and A5 piezometers, which have more-or-less transmissive behavior.

The results for the N5 formation remain consistent among nearby boreholes {25}, {26} and {27}. For the N3 formation, the variation in values at the same depth is significant, with differences of $10^{-4}$ m/s between two nearby boreholes ({25} and {26}). The N3 formation contains many more marly layers than the N5 formation. This explains the much greater permeability variation within this formation and also between the two formations.
# Pumping Tests: Permeability at Aquifer Scale

Analyses of the pumping test and the disturbance in the only responsive observation piezometer give permeability values between $10^{-8}$ and $10^{-7}$ m/s. The difference between these values can be explained by the fact that different volumes were studied due to the heterogeneity of the medium being excited by the tests.

The storage coefficient for both pumped and interference boreholes obtained from these tests is approximately $10^{-3}$ m/s.

## 4.3. SNMR Data

### 4.3.1. SNMR Signal

The lowest initial amplitude (E0) over all soundings is between 5 and 10 nV. Initial amplitudes associated with the highest pulse moments range between 1 and 5 nV. These very low initial amplitude values do, however, allow some interpretation, since the signal-to-noise ratio of the soundings (Table 1) is always greater than 2.5.

The SNMR signals are provided in Figure 6. Uncertainty bounds are derived by considering minimum and maximum values from all SNMR signal parameter sets (signal amplitude e and T2* relaxation time) that yield a 10% increase in the root mean squared error (RMSE) fitting error as compared to the optimal RMSE, as proposed by [51].

<figure><img src="image_9.png" /><figcaption>**Figure 6.** SNMR signal parameters as a function of the energizing pulse, for all soundings (Σ Pj to Pj – 7 = sum of rainfall over the 7 days preceding the sounding). (A)—Initial signal amplitude. (B)—Apparent transverse relaxation time T2*.</figcaption></figure>

Two groups of responses stand out. The first group, consisting of soundings {38}_1 and {38}_2, has initial signal amplitudes approximately 3 nV lower at low pulse moments than the group consisting of soundings {38}_3, {38}_4, {37} and {39}, which suggests lower water content in the upper part of the sounding. High uncertainties are associated with T2* estimates. Sounding {37} stands out with relatively lower T2*, which may be associated with small porosity dimensions. These original data (Figure 6) are unaffected by any inversion artefacts. It is therefore interesting to compare these trends with the inverted results (Figure 7).
<figure><img src="image_10.png" /><figcaption>**Figure 7.** SNMR inversion for all sites. (A)—SNMR water content. (B)—Apparent transverse relaxation time T2*.</figcaption></figure>

### 4.3.2. Inversion Results

The inversion results are displayed in Figure 7 with a 95% confidence interval on SNMR water content derived from the SVD inversion [52]. This approach provides an estimate of the uncertainty bounds of the water content in each layer considered separately. It accounts for the fact the buried water layers that have equal water volumes can, to some extent, yield equivalent SNMR responses. The results show very low to low water contents, which do not exceed 0.9% (1.5% considering the maximum inversion uncertainties), and that have an average value of 0.42%, with highest water content being reached in the first 0-5m (Figure 7A). These results indicate and quantify the low water storage capacity of the Sublithographic limestone matrix. In addition, the trends seen in Figure 6 are confirmed with these inverted data, since the same trends are present (two response groups and a very limited quantity of water content).

The T2* estimates are consistent with a filtered, non-inverted SNMR signal, with significantly lower values associated with Sounding {37} (see Section 4.3.1). High values associated with Sounding {38}_2 are not significant as they are associated with high signal uncertainties.

Soundings performed under different hydrologic conditions exhibit contrasting behaviors in the upper 0–5 m. Soundings {38}_1 and {38}_2 were performed in low water levels with significant antecedent precipitations. They have lowest water content (0.3% in average). Soundings {37} and {39} were performed in high water levels with little antecedent
precipitations. The have intermediate water content (0.55% in average). Soundings {38}_3 and {38}_4 were performed in high water levels with important antecedent precipitations. The have the highest water content (0.8% in average). SNMR inversions thus evidence both a primarily seasonal time-scale variation on the limestone water content, and a significant but secondary short (daily) time-scale variation.

The uncertainties provided in Figure 7 show the minimum and maximum water content that may be associated with each water layer, taken separately. Maximum water content for a specific water layer is associated with minimum or at least reduced water content in other layers. Such representation of uncertainty is not adequate to assess the minimum and maximum water content over the whole sounding. As a complement to Figure 7, we thus provide uncertainty bounds for sounding {38} for equivalent models that yield the minimum and the maximum water volumes over the investigated depth [53] (Figure 8). The average water content is between 0.57–0.61% and 0.51–0.63% for soundings {38}_4 and {38}_3, respectively, while it is only between 0.19–0.25% and 0.33–0.35% for soundings {38}_2 and {38}_1, respectively.

<figure><img src="image_11.png" /><figcaption>**Figure 8.** SNMR inversion for soundings performed at location {38}. Temporal water content variability between high ({38}_3 and {38}_4) and low ({38}_1 and {38}_2) water levels can be significant.</figcaption></figure>

**5. Discussion**

Large borehole networks are often installed with geotechnical objectives in mind, and the data collected as part of these studies are often not linked to hydrogeological studies as such. Although these boreholes can be reused as piezometers, they were not originally designed as such. In many cases, these boreholes were dug throughout the structure and are not cased. This greatly complicates the interpretation of water level variations in response to internal or external forcing. Our study site benefits from such a network of discrete boreholes set up for geotechnical purposes. In the following discussion, we examine the strengths and weaknesses of this network, as well as possible ways of overcoming them.
5.1. *Piezometric Response and Fracture Network Connectivity*

Statistical analysis of these chronicles requires a large dataset and, depending on the
indicators used, makes it possible to classify the different piezometric variations according
to their response to rainfall. This descriptive analysis made it possible to characterize five
behaviors, grouping 19 piezometers and identifying transmissive (groups A1 and A2),
an inertial (group A3), and a group with mixed behaviors (groups A4 and A5) (Figure 9).

<figure><img src="image_12.png" /><figcaption>**Figure 9.** Classification of hydrodynamic behavior of water levels in piezometers in response to rainfall events.</figcaption></figure>

Considered alone, piezometric response does allow some relevant insights into the
hydrodynamics behavior. This is the case for Groups A4 and A5, where the analysis of
the piezometric response provides information on the more-or-less inertial functioning
of the boreholes, and suggests the presence of fractures in the process of development or
a poorly connected fracture network. It is also the case for Group A1, for which neither
geological drillings allowing fracture analysis, nor geophysical logs, are available near the
piezometers. According to the descriptive analysis, Group A1 behaves like a permeable
medium with a rapid response time to rainfall events but a relatively long memory effect
with a slow decrease in water level. We associated such behavior with a moderate degree
of fracturing (Figure 3B), with the rapid water level rise in the piezometer being simply
due to the filling of the piezometer and not significative of the hydrodynamic response
of the water table in the piezometer area. The SNMR results confirm this interpretation,
as the significantly lower T2* associated with sounding {37} hints at a significative matrix
contribution to the effective porosity.

Detailed fracture analysis may evidence discrete features that are not relevant for
understanding hydrodynamic responses. For example, group A3 is associated with matrix
functioning, despite the presence of open fractures near one piezometer in this group
(Figure 3D). The orientation of these open fractures is perpendicular to the main flow,
which may explain the inertial behavior of this piezometer.

In other cases, detailed fracture analysis using geological boreholes and SNMR can
avoid erroneous water level dynamics interpretation. In the example provided in Figure 10,
optical imaging allows to identify a karst conduit that is responsible for the regulation of the
groundwater level. SNMR sounding (39), located close to piezometer (16) (Figure 10), did
not record evidence of significative water content below this fracture. We thus formulate
the hypothesis that the water level below this base level does not contribute to the flow
and may not be considered as representative from the water table.
<figure><img src="image_13.png" /><figcaption>**Figure 10.** Link between conduit and piezometric level: the karst conduit maintains a base level in the piezometer. Only water levels above this base level contribute to the flow.</figcaption></figure>

## 5.2. Hydraulic Tests and Local Matrix Characterization

In karstified environments, permeabilities can vary by several orders of magnitude depending on the scale of measurement (laboratory, borehole, and catchment) [27]. For each of these permeability ranges, different measurement methods are used. Hydraulic tests carried out in this study are mostly sensitive to the near borehole media.

Permeability measurements carried out indicate higher permeability at the surface of the Sublithographic limestones and lower permeability at depth in these limestones, in the marls, and in the marlier limestone layers. These values highlight a semi-permeable behavior for limestones to be impermeable when marls predominate, whereas the piezometric response tends to highlight the reactivity of the system.

In our study, the duration of the Lefranc tests is too short and the pumping test flow rates are too low to attract the more distant water contained within the carbonate matrix or connected fractures, so they mainly provide local matrix characterization.

## 5.3. SNMR and Water Content Quantification

Information on water mobilization in response to forcing obtained from piezometric responses, fracturing analysis, and hydraulic tests remains incomplete as they do not quantify the water saturation of the medium. SNMR measurements provide access to the water content and its variation over a wider access radius than those investigated by borehole tests or piezometric variations.

In our study, very low (3 and 8 nV) SNMR signals associated with high signal-to-noise ratio allow us to quantify very low (<0.6%) water content with significative variability.

We relate this variability to both a primarily seasonal time-scale variation of the limestone water content and a significant but secondary short (daily) time-scale variation. This result reinforces the interpretation of the piezometric behavior of matrix-like and mixed matrix/transmissive groups, with mostly seasonal matrix dynamics and daily-to-weekly fissures and fracture dynamics.
## 6. Conclusions

This study utilized a variety of methods to better understand the functioning of karst groundwater systems. With karst outcrops covering 10-15% of the world's continental surface [21], understanding these hydrosystems is crucial for addressing the challenges associated with the exploitation of groundwater resources. It is estimated that around 25% of the world's population relies on groundwater extracted from karst areas for drinking water [52]. Additionally, studying these systems is essential during geotechnical or civil engineering investigations [4].

Karst hydrosystems are complex due to their triple porosity (matrix, fractures, and conduits) and inherent heterogeneity. Both the fracture properties and the hydrodynamic properties of the carbonate matrix determine the functioning of karst aquifers. Piezometric information is often the basis for understanding flow organization. Many studies rely on discrete piezometric data to characterize these complex environments where flows are distributed according to permeability constraints.

Our study benefits from a large network of boreholes initially implemented for geotechnical purposes. Despite the density of this observation network, the collection of water level variations is not easily linked to the hydrogeological properties of the aquifer. While piezometric responses alone provide some insights into hydrodynamic behavior, detailed fracture analysis can reveal discrete features that are not always relevant for understanding these responses. However, in some cases, combining detailed fracture analysis from geological boreholes with SNMR is necessary to avoid erroneous interpretations of water level dynamics.

High-density and high-quality data in karst environments are rarely available due to the complexity and heterogeneity of these systems. Comprehensive datasets are crucial for accurate analysis and modeling, yet they are challenging to obtain. By leveraging the extensive and high-quality data from this unique site, we can gain valuable insights into the hydrogeologic complexities of fractured-rock karst systems, providing a strong basis for effectively understanding and managing these environments.

The combination of geological, hydrogeological, and geophysical measurements in this project highlights the contribution of water masses with two types of aquifer responses: a rapid response at hourly-to-daily timescales to rainfall events in well-connected fractures, and a slower response at seasonal timescales in the matrix. The characterization of the dynamics and the quantification of the water content in the matrix and fissured/fractured compartments were enabled by the combined interpretation of piezometric response analyses and SNMR data.

Thus, coupling various separately collected measurements provides a key to understanding the aquifer as a whole. Future work could include longer pumping tests or modelling fluctuations of base groundwater levels to better assess matrix flows.

**Author Contributions:** Conceptualization, M.B., C.B. and N.M.; methodology, M.B., L.-M.G., V.B., N.M., R.V., K.C., C.B., S.D., S.G. and G.C.; software, M.B., N.M. and G.C.; validation, N.M., R.V. and K.C.; formal analysis, M.B.; investigation, M.B., L.-M.G. and V.B.; resources, S.G.; data curation, M.B., L.-M.G. and V.B.; writing—original draft preparation, M.B., C.B. and N.M.; writing—review and editing, C.B., N.M., K.C., R.V., S.D. and S.G.; Supervision, C.B., N.M., K.C. and S.G.; project administration, C.B., K.C. and S.G.; funding acquisition, S.G., C.B., K.C., N.M. and S.D. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by ANDRA PhD grant number 20080158.

**Data Availability Statement:** Data are contained within the article.

**Acknowledgments:** We would like to thank the reviewers for their work and constructive comments and suggestions. We would also like to express our sincere thanks to Patricia Bobeck for her proof-reading and for her contribution to the English version of the manuscript.

**Conflicts of Interest:** The authors declare no conflict of interest.
## Appendix A

<table><caption>Table A1. Links between number, name, and information.</caption><thead><tr><th>Number</th><th>Name</th><th>Information</th><th>Number</th><th>Name</th><th>Information</th></tr></thead><tbody><tr><td>1</td><td>CIG1028</td><td>Piezometry</td><td>31</td><td>CIG1031</td><td>Borehole imaging</td></tr><tr><td>2</td><td>CIG1029</td><td>Piezometry</td><td>32</td><td>CIG1080</td><td>Borehole imaging</td></tr><tr><td>3</td><td>CIG1030</td><td>Piezometry</td><td>33</td><td>CIG1088</td><td>Borehole imaging</td></tr><tr><td>4</td><td>CIG1031</td><td>Piezometry</td><td>34</td><td>CIG1120</td><td>Borehole imaging</td></tr><tr><td>5</td><td>CIG1227</td><td>Piezometry</td><td>35</td><td>CIG1122</td><td>Borehole imaging</td></tr><tr><td>6</td><td>CIG1228</td><td>Piezometry</td><td>36</td><td>CIG1254</td><td>Borehole imaging</td></tr><tr><td>7</td><td>CIG1229</td><td>Piezometry</td><td>37</td><td>EST3BIS</td><td>SNMR</td></tr><tr><td>8</td><td>CIG1230</td><td>Piezometry</td><td>38</td><td>EST14</td><td>SNMR</td></tr><tr><td>9</td><td>CIG1237</td><td>Piezometry</td><td>39</td><td>EST18</td><td>SNMR</td></tr><tr><td>10</td><td>CIG1245</td><td>Piezometry</td><td>40</td><td>CIG1005</td><td>Fracturing</td></tr><tr><td>11</td><td>CIG1246</td><td>Piezometry</td><td>41</td><td>CIG1009</td><td>Fracturing</td></tr><tr><td>12</td><td>CIG1247</td><td>Piezometry</td><td>42</td><td>CIG1011</td><td>Fracturing</td></tr><tr><td>13</td><td>CIG1248</td><td>Piezometry</td><td>43</td><td>CIG1077</td><td>Fracturing</td></tr><tr><td>14</td><td>CIG1249</td><td>Piezometry</td><td>44</td><td>CIG1078</td><td>Fracturing</td></tr><tr><td>15</td><td>CIG1250</td><td>Piezometry</td><td>45</td><td>CIG1079</td><td>Fracturing</td></tr><tr><td>16</td><td>EST1012</td><td>Piezometry</td><td>46</td><td>CIG1080</td><td>Fracturing</td></tr><tr><td>17</td><td>EST1021</td><td>Piezometry</td><td>47</td><td>CIG1086</td><td>Fracturing</td></tr><tr><td>18</td><td>EST1037</td><td>Piezometry</td><td>48</td><td>CIG1088</td><td>Fracturing</td></tr><tr><td>19</td><td>EST5071</td><td>Piezometry</td><td>49</td><td>CIG1090</td><td>Fracturing</td></tr><tr><td>20</td><td>CIG1023</td><td>Pumping tests</td><td>50</td><td>CIG1108</td><td>Fracturing</td></tr><tr><td>21</td><td>CIG1024</td><td>Pumping tests</td><td>51</td><td>CIG1114</td><td>Fracturing</td></tr><tr><td>22</td><td>CIG1025</td><td>Pumping tests</td><td>52</td><td>CIG1116</td><td>Fracturing</td></tr><tr><td>23</td><td>CIG1026</td><td>Pumping tests</td><td>53</td><td>CIG1117</td><td>Fracturing</td></tr><tr><td>24</td><td>CIG1137</td><td>Lefranc tests</td><td>54</td><td>CIG1118</td><td>Fracturing</td></tr><tr><td>25</td><td>CIG1138</td><td>Lefranc tests</td><td>55</td><td>CIG1120</td><td>Fracturing</td></tr><tr><td>26</td><td>CIG1139</td><td>Lefranc tests</td><td>56</td><td>CIG1121</td><td>Fracturing</td></tr><tr><td>27</td><td>CIG1140</td><td>Lefranc tests</td><td>57</td><td>CIG1252</td><td>Fracturing</td></tr><tr><td>28</td><td>CIG1237</td><td>Lefranc tests</td><td>58</td><td>CIG1253</td><td>Fracturing</td></tr><tr><td>29</td><td>CIG1250</td><td>Lefranc tests</td><td>59</td><td>CIG1254</td><td>Fracturing</td></tr><tr><td>30</td><td>CIG1011</td><td>Borehole imaging</td><td></td><td></td><td></td></tr></tbody></table>
## Appendix B

<table><caption>**Table A2.** Characteristics of the piezometers used in this study. Type and data-acquisition start dates.</caption><thead><tr><th>Name</th><th>Depth (m)</th><th>Borehole Screen</th><th>Geological Level</th><th>Data</th><th>Data Acquisition Start On</th></tr></thead><tbody><tr><td>{16}</td><td>25.5</td><td>Ø 112/122.5 mm + gravel packing</td><td>Sublithographic limestone (N1)</td><td>Pressure Temperature</td><td>13 April 1996<br>30 June 1999</td></tr><tr><td>{17}</td><td>30</td><td>Ø 113/125 mm opening 1 mm + gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature Sampling</td><td>4 April 2000<br>6 July 2000<br>6 July 2000</td></tr><tr><td>{18}</td><td>20.2</td><td>Ø 113/125 mm opening 1 mm + gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature Sampling</td><td>5 April 2000</td></tr><tr><td>{19}</td><td>39.06</td><td>Ø 119.8/125 mm</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>20 December 1995<br>10 July 1999</td></tr><tr><td>{1}</td><td>20</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure</td><td>25 January 2016</td></tr><tr><td>{2}</td><td>14</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure</td><td>17 February 2016</td></tr><tr><td>{3}</td><td>28</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure</td><td>25 February 2016</td></tr><tr><td>{4}</td><td>30</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure</td><td>27 January 2016</td></tr><tr><td>{5}</td><td>30.3</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>31 August 2016</td></tr><tr><td>{6}</td><td>30.3</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>31 August 2016</td></tr><tr><td>{7}</td><td>30.3</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>22 August 2016</td></tr><tr><td>{8}</td><td>32</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>31 August 2016</td></tr><tr><td>{9}</td><td>30.3</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>31 August 2016</td></tr><tr><td>{10}</td><td>20.3</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>18 August 2016</td></tr><tr><td>{11}</td><td>20.3</td><td>gravel packing</td><td>Sublithographic limestone (N1)</td><td>Pressure Temperature</td><td>16 June 2016</td></tr><tr><td>{12}</td><td>15.4</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>29 September 2016</td></tr><tr><td>{13}</td><td>30.35</td><td>gravel packing</td><td>Kimmeridgian marls</td><td>Pressure Temperature</td><td>22 August 2016</td></tr><tr><td>{14}</td><td>15.3</td><td>gravel packing</td><td>Sublithographic limestone (N5)</td><td>Pressure Temperature</td><td>24 November 2016</td></tr><tr><td>{15}</td><td>30.3</td><td>gravel packing</td><td>Sublithographic limestone (N4)</td><td>Pressure Temperature</td><td>3 August 2016</td></tr></tbody></table>

## References

1. Stevanović, Z. Global Distribution and Use of Water from Karst Aquifers. Geol. Soc. Lond. Spec. Publ. **2018**, *466*, 217–236. [CrossRef]

2. Ravbar, N.; Šebela, S. The Effectiveness of Protection Policies and Legislative Framework with Special Regard to Karst Landscapes: Insights from Slovenia. Environ. Sci. Policy **2015**, *51*, 106–116. [CrossRef]

3. Hartmann, A.; Goldscheider, N.; Wagener, T.; Lange, J.; Weiler, M. Karst Water Resources in a Changing World: Review of Hydrological Modeling Approaches: Karst Water Resources Prediction. Rev. Geophys. **2014**, *52*, 218–242. [CrossRef]

4. Chen, C.-T.; Hu, J.-C.; Lu, C.-Y.; Lee, J.-C.; Chan, Y.-C. Thirty-Year Land Elevation Change from Subsidence to Uplift Following the Termination of Groundwater Pumping and Its Geological Implications in the Metropolitan Taipei Basin, Northern Taiwan. Eng. Geol. **2007**, *95*, 30–47. [CrossRef]
5. Mangin, A. Pour une meilleure connaissance des systèmes hydrologiques à partir des analyses corrélatoire et spectrale (The use of autocorrelation and spectral analyses to obtain a better understanding of hydrological systems). *J. Hydrol.* **1984**, *67*, 25-43.
[CrossRef]

6. Nistor, M.M.; Rahardjo, H.; Satyanaga, A.; Hao, K.Z.; Xiaosheng, Q.; Sham, A.W.L. Investigation of Groundwater Table Distribution Using Borehole Piezometer Data Interpolation: Case Study of Singapore. *Eng. Geol.* **2020**, *271*, 105590. [CrossRef]

7. De Luca, D.A.; Lasagna, M.; Debernardi, L. Hydrogeology of the Western Po Plain (Piedmont, NW Italy). *J. Maps* **2020**, *16*, 265-273. [CrossRef]

8. Fiorillo, F. Tank-Reservoir Drainage as a Simulation of the Recession Limb of Karst Spring Hydrographs. *Hydrogeol. J* **2011**, *19*, 1009-1019. [CrossRef]

9. Bennett, G.; Van Camp, M.; Shemsanga, C.; Kervyn, M.; Walraevens, K. Assessment of Spatial and Temporal Variability of Groundwater Level in the Aquifer System on the Flanks of Mount Meru, Northern Tanzania. *J. Hydrol. Reg. Stud.* **2022**, *44*, 101212. [CrossRef]

10. Giese, M.; Reimann, T.; Bailly-Comte, V.; Maréchal, J.-C.; Sauter, M.; Geyer, T. Turbulent and Laminar Flow in Karst Conduits Under Unsteady Flow Conditions: Interpretation of Pumping Tests by Discrete Conduit-Continuum Modeling. *Water Resour. Res.* **2018**, *54*, 1918-1933. [CrossRef]

11. Praamsma, T.; Novakowski, K.; Kyser, K.; Hall, K. Using Stable Isotopes and Hydraulic Head Data to Investigate Groundwater Recharge and Discharge in a Fractured Rock Aquifer. *J. Hydrol.* **2009**, *366*, 35-45. [CrossRef]

12. Maréchal, J.-C.; Ladouche, B.; Dörfliger, N.; Lachassagne, P. Interpretation of Pumping Tests in a Mixed Flow Karst System: Pumping Tests in Karst System. *Water Resour. Res.* **2008**, *44*, W05401. [CrossRef]

13. Banks, E.W.; Simmons, C.T.; Love, A.J.; Shand, P. Assessing Spatial and Temporal Connectivity between Surface Water and Groundwater in a Regional Catchment: Implications for Regional Scale Water Quantity and Quality. *J. Hydrol.* **2011**, *404*, 30-49. [CrossRef]

14. Valois, R.; Rau, G.C.; Vouillamoz, J.-M.; Derode, B. Estimating Hydraulic Properties of the Shallow Subsurface Using the Groundwater Response to Earth and Atmospheric Tides: A Comparison With Pumping Tests. *Water Resour. Res.* **2022**, *58*, e2021WR031666. [CrossRef]

15. Bon, A.F.; Aoudou Doua, S.; Banakeng, L.A.; Narke, C.; Chouto, S.; Mbouombouo Ndam, A. Contribution of a Geostatistical Model of Electrical Conductivity in the Assessment of the Water Pollution Index of the Quaternary Aquifer of the Lake Chad Basin (Kousseri-Cameroon). *Arab. J. Geosci.* **2020**, *13*, 170. [CrossRef]

16. Dewandel, B.; Aunay, B.; Maréchal, J.C.; Roques, C.; Bour, O.; Mougin, B.; Aquilina, L. Analytical Solutions for Analysing Pumping Tests in a Sub-Vertical and Anisotropic Fault Zone Draining Shallow Aquifers. *J. Hydrol.* **2014**, *509*, 115-131. [CrossRef]

17. Castellazzi, P.; Arroyo-Domínguez, N.; Martel, R.; Calderhead, A.I.; Normand, J.C.L.; Gárfias, J.; Rivera, A. Land Subsidence in Major Cities of Central Mexico: Interpreting InSAR-Derived Land Subsidence Mapping with Hydrogeological Data. *Int. J. Appl. Earth Obs. Geoinf.* **2016**, *47*, 102-111. [CrossRef]

18. Haimson, B.C.; Doe, T.W. State of Stress, Permeability, and Fractures in the Precambrian Granite of Northern Illinois. *J. Geophys. Res.* **1983**, *88*, 7355-7371. [CrossRef]

19. Dausse, A.; Leonardi, V.; Jourde, H. Hydraulic Characterization and Identification of Flow-Bearing Structures Based on Multi-Scale Investigations Applied to the Lez Karst Aquifer. *J. Hydrol. Reg. Stud.* **2019**, *26*, 100627. [CrossRef]

20. Goldscheider, N.; Drew, D. *Methods in Karst Hydrogeology*; International Contributions to Hydrogeology; Taylor & Francis: London, UK, 2007; ISBN 978-0-415-42873-6.

21. Kiraly, L. Karstification and Groundwater Flow. *Speleogenesis Evol. Karst Aquifers* **2003**, *1*, 1-26.

22. Ford, D.C.; Williams, P.W. *Karst Hydrogeology and Geomorphology*; Chapman and Hall: London, UK, 2007. [CrossRef]

23. Bakalowicz, M. Karst Groundwater: A Challenge for New Resources. *Hydrogeol. J.* **2005**, *13*, 148-160. [CrossRef]

24. Bakalowicz, M. Le milieu karstique: Études et perspectives, identification et caractérisation de la ressource. In Proceedings of the CFH—Colloque Hydrogéologie et Karst au Travers des Travaux de Michel Lepiller 17 Mai 2008, Orléans, France, 16-17 May 2008.

25. Atkinson, T.C. Diffuse Flow and Conduit Flow in Limestone Terrain in the Mendip Hills, Somerset (Great Britain). *J. Hydrol.* **1977**, *35*, 93-110. [CrossRef]

26. Shuster, E.T.; White, W.B. Seasonal Fluctuations in the Chemistry of Lime-Stone Springs: A Possible Means for Characterizing Carbonate Aquifers. *J. Hydrol.* **1971**, *14*, 93-128. [CrossRef]

27. Clauser, C. Permeability of Crystalline Rocks. *Eos* **1992**, *73*, 233-238. [CrossRef]

28. Jeannin, P.-Y. Structure et Comportement Hydraulique des Aquifères Karstiques (Structure and Hydraulic Behavior of Karst Aquifers). Ph.D. Thesis, Université de Neuchâtel, Neuchâtel, Switzerland, 1996.

29. Delbart, C.; Valdés, D.; Barbecot, F.; Tognelli, A.; Couchoux, L. Spatial Organization of the Impulse Response in a Karst Aquifer. *J. Hydrol.* **2016**, *537*, 18-26. [CrossRef]

30. Chalikakis, K.; Plagnes, V.; Guerin, R.; Valois, R.; Bosch, F.P. Contribution of Geophysical Methods to Karst-System Exploration: An Overview. *Hydrogeol. J.* **2011**, *19*, 1169-1180. [CrossRef]

31. Chalikakis, K.; Nielsen, M.R.; Legchenko, A. MRS Applicability for a Study of Glacial Sedimentary Aquifers in Central Jutland, Denmark. *J. Appl. Geophys.* **2008**, *12*, 176-187. [CrossRef]

32. Legchenko, A.; Valla, P. A Review of the Basic Principles for Proton Magnetic Resonance Sounding Measurements. *J. Appl. Geophys.* **2002**, *50*, 3-19. [CrossRef]
33. Mazzilli, N.; Chalikakis, K.; Carrière, S.D. Surface Nuclear Magnetic Resonance Monitoring Reveals Karst Unsaturated Zone Recharge Dynamics during a Rain Event. *Water* **2020**, *12*, 3183. [CrossRef]

34. Scholz, E. *Synthèse Hydrogéologique des Calcaires du Barrois*; ANDRA: Bure, France, 2011; p. 97.

35. Jaillet, S. Un Karst de Bas Plateau: Le Barrois, Structure, Fonctionnement, Evolution (A Low Plateau Karst: The Barrois, Structure, Functioning, Evolution). Ph.D. Thesis, Université de Bordeaux, Bordeaux, France, 2000.

36. Decloux, J.-P. *Hydrologie et Hydrogeologie des Aquifères de Surface sur les Bassins de la Saulx et de L'ornain—Rapport de Tâche 1 (RT1)*; ANDRA: Bure, France, 2018; p. 47.

37. Mangin, A. Contribution à L'étude Hydrodynamique des Aquifères Karstiques. Ph.D. Thesis, Université de Dijon, Dijon, France, 1975.

38. Cholet, C. Fonctionnement Hydrogéologique et Processus de Transport dans les Aquifères Karstiques du Massif du Jura. Ph.D. Thesis, Université Bourgogne Franche-Comté, Besançon, France, 2017.

39. Delbart, C. Variabilité Spatio-Temporelle du Fonctionnement d'un Aquifère Karstique du Dogger: Suivis Hydrodynamiques et Géochimiques Multifréquences; Traitement du Signal des Réponses Physiques et Géochimiques. Ph.D. Thesis, Université Paris-Sud—Paris XI, Paris, France, 2013.

40. Hvorslev, J. *Time Lag and Soil Permeability in Ground*; Waterways Experiment Station: Vicksburg, MS, USA, 1951; Volume 50.

41. ISO 22282-2; Reconnaissance et Essais Géotechniques—Essais Géohydrauliques—Partie 2: Essai de Perméabilité à L'eau Dans un Forage en Tube Ouvert. ISO: Geneva, Switzerland, 2013.

42. *NF P94-130 Standard*; Soils: Investigation and Testing—Pumping Test. AFNOR: Paris, France, 2002.

43. Moench, A.F. Flow to a Well of Finite Diameter in a Homogeneous, Anisotropic Water Table Aquifer. *Water Resour. Res.* **1997**, *33*, 1397–1407. [CrossRef]

44. Theis, C.V. The Relation between the Lowering of the Piezometric Surface and the Rate and Duration of Discharge of a Well Using Ground. *Eos Trans. AGU* **1935**, *16*, 519–524. [CrossRef]

45. Kruseman, G.P.; De Ridder, N.A. *Analysis and Evaluation of Pumping Test Data*; International Institute for Land Reclamation and Improvement: Wageningen, The Netherlands, 1991; Volume 47.

46. Vouillamoz, J.-M. La Caractérisation des Aquifères par une Méthode Non Invasive: Les Sondages par Resonance Magnetique Protonique. Ph.D. Thesis, Université de Paris XI, Paris, France, 2003.

47. Legchenko, A.; Baltassat, J.-M.; Beauce, A.; Bernard, J. Nuclear Magnetic Resonance as a Geophysical Tool for Hydrogeologists. *J. Appl. Geophys.* **2002**, *50*, 21–46. [CrossRef]

48. Legchenko, A.; Valla, P. Processing of Surface Proton Magnetic Resonance Signals Using Non-Linear Fitting. *J. Appl. Geophys.* **1998**, *39*, 77–83. [CrossRef]

49. Legchenko, A.; Valla, P. Removal of Power-Line Harmonics from Proton Magnetic Resonance Measurements. *J. Appl. Geophys.* **2003**, *53*, 103–120. [CrossRef]

50. Legchenko, A. MRS Measurements and Inversion in Presence of EM Noise. *Bol. Geol. Min.* **2007**, *118*, 489–508.

51. Mazzilli, N.; Boucher, M.; Chalikakis, K.; Legchenko, A.; Jourde, H.; Champollion, C. Contribution of Magnetic Resonance Soundings for Characterizing Water Storage in the Unsaturated Zone of Karst Aquifers. *Geophysics* **2016**, *81*, WB49–WB61. [CrossRef]

52. Legchenko, A. *Magnetic Resonance Imaging for Groundwater*; Wiley-ISTE: Hoboken, NJ, USA, 2013; ISBN 978-1-84821-568-9.

53. Legchenko, A.; Comte, J.-C.; Ofterdinger, U.; Vouillamoz, J.-M.; Lawson, F.M.A.; Walsh, J. Joint Use of Singular Value Decomposition and Monte-Carlo Simulation for Estimating Uncertainty in Surface NMR Inversion. *J. Appl. Geophys.* **2017**, *144*, 28–36. [CrossRef]

**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.