1 Introduction

1.1 Mud volcanoes

Mud volcanoes (MVs) are geological structures that are the result of eruptions or extrusions of deep-lying fluids, such as mud, mineralized water and hydrocarbons. Eruptions involve the effusive or explosive ejection of deep-lying fluids with a significant range of intensity, magnitude, explosivity, eruption rate; extrusions involve the effusive, quiet rather than explosive ejection of deep-lying fluids from a vent or fissure where gas, if present, is released quietly, rather than explosively.

Surface morphology is usually characterized by a cone with a diameter ranging from a few meters to several hundred meters. Morphology is principally a function of duct size, interstitial fluid pressure, emission frequency and type and viscosity of the material emitted. The density of the material emitted dictates the maximum height the mud volcano (MV) can reach. Frequently, minor cones—‘griphons’—and ponds can be found on the flanks of the volcanic building.

Emissions may be continuous or intermittent. In some cases, this activity may be interrupted by violent and instantaneous explosion (“paroxysm”), with an eruptive column composed mainly of clay material, water, and gases. The erupted material, may also consist of mud breccias composed of a mud matrix with chaotically distributed angular to rounded rock clasts from a few millimetres to metres in diameter, and may reach a large distance from the emission point (Gattuso et al. 2021; Mazzini et al. 2021a).

MVs are found both on- and offshore. Most are located in areas that are subject to compressional tectonics and associated with sediment accretion at convergent margins (Aslan et al. 2001; Brown and Westbrook 1988; Di Giulio 1992; Dimitrov 2002; Etiope et al. 2002; Graue 2000; Henry et al. 1990; Higgins and Saunders 1974; Kopf et al. 1998; Le Pichon et al. 1990; Martinelli and Judd 2004; Mascle et al. 2001; Milkov 2005; Reed et al. 1990; Kopf 2002). MVs arise in the western Alps of Italy (Martinelli 1999; Martinelli and Judd 2004), in the Apennines (Canora et al. 2012; Rainone et al. 2015) and in Sicily (Etiope et al. 2007, 2002).

Mazzini and Etiope (2017) provided an updated review of the knowledge and implications of MVs. The complex mechanisms and processes involved in the formation of MVs (Nanni and Zuppi 1986; Kopf 2002) are still not completely understood. Overpressure and under-consolidation are the main features that cause the upward migration of water and gas-rich mud often along fault planes and open fracture lines (Ansted 1866; Brown, 1990; Dimitrov, 2002; Etiope and Klusman 2010; Fertl 1976; Higgins and Saunders 1974; Hovland et al. 1997; Kopf 2002; Milkov 2000; Nath et al. 2008; Niemann and Boetius 2010; Yassir 1989). The main engine leading to the formation of MVs is generally a combination of gravitative instability of sediments and fluid overpressures (e.g. Kopf 2002; Revil 2002). Additional mechanisms may include thermal effect in pore fluids, dehydration reactions and disequilibrium compaction (Revil 2002). As reported by Mazzini and Etiope (2017), the MV formation should consist of the following processes: gravitative instability (shale buoyancy); fluid overpressure; hydrofracturing; fluid flow along permeable open fractures or faults. The final stage of MV growth is its manifestation on the surface. This can happen gradually with progressive and slow release of mud and fluids, or with violent, paroxysmal events (eruption).

Another peculiarity of mud volcanism is the transport of breccia to the surface. It is well known that some clasts erupted at MV sites originate from several kilometres deep (Mazzini and Etiope 2017). Mazzini et al. (2021a), in their study of the Lokbatan MV in Azerbaijan, which is possibly the most active MV on Earth, found that powerful eruptions are also able to break apart mud breccia cones and transport blocks with the mud breccia flows. They also proposed that the Lokbatan MV effectively seals subsurface conduits after eruptions, allowing pressure to recharge more quickly than at other MVs, leading to more frequent, more energetic, and larger eruptions. The authors suggested that the correlation between the ability of an active system to self-seal and the energy of eruptions, could be valid for all MVs. For the reasons stated above, MVs are not only of significant geological and geomorphological interest; they are also non-negligible sources of hazard and risk: the most hazardous phenomena at MVs are the paroxysm events with their indirect phenomena associated (Gattuso et al. 2021). They are instantaneous, nowday unpredictable and caused different injuries/died in the world (Baloglanov et al. 2018; Davies et al. 2008; Manga and Bonini 2012; Mazzini et al. 2007). Therefore, understanding their subsurface structure and circuits is crucial to defining hazard and risk conditions in any adjacent areas (Gattuso et al. 2021).

Several scientific studies have focused on the application of geophysical surveys in MV areas. In many of these studies, electrical resistivity surveys has been the method of choice used to capture images of the near-surface structure of MVs (Keller and Rapolla 1974; Ping-Yu et al. 2011; Zeyen et al. 2011). Electrical resistivity assessments (pseudo-3D) of MVs in New Zealand and Venezuela are discussed in, respectively, Zeyen et al. (2011) and Blake et al. (2021). Zeyen et al. (2011) found a complex structure of the MV at shallow depths (< 30 m). They found that the deeper areas of rising mud are not well correlated to vent distribution: there is no correspondence between the lowest resistivities associated with the presence of mineralized water at a depth of 8–20 m and the vent concentration areas, nor the position of the main vent. Moreover, the authors found no connection between the upper and the lower mud reservoirs below the main vent. Nonetheless, the existence of a narrow mud conduit, up to 3 m wide, at this location, cannot be ruled out as the authors carried out a pseudo-3D rather than an effective 3D survey (unlike in this study), which guarantees adequate accuracy at shallow depths.

In certain instances, electrical resistivity surveys have been used alongside seismic surveys (Accaino et al. 2007; Grassi et al. 2022; Rainone et al. 2015). Electromagnetic applications (VLF-EM) carried out in Taiwan were discussed by Lin and Jeng (2010) while Grassi et al. (2022) and Imposa et al. (2016) employed passive seismic surveys to map the underground structure of an MV.

1.2 The investigated area, its previous studies and aims of the study

The study area is located in the province of Teramo, in the region of Abruzzo (Central Italy). It is approximately 1 km west of the town of Pineto and 1.4 km from the Adriatic Sea, on the north bank of the ‘Calvano’ stream (Fig. 1). The Pineto or ‘Cenerone’ MV, in this paper referred to as ‘Cenerone-Pineto’ MV, is mentioned in the Italian MV database compiled by Martinelli and Judd (2004).

Fig. 1
figure 1

a Geographical setting of the ‘Cenerone-Pineto’ mud volcano. b Outline of the main features of the ‘Cenerone-Pineto’ mud volcano in plan view (from Google EarthTM)

Other researchers studied this MV. Some researchers focused principally on the hydrochemical aspects (Desiderio and Rusi 2004; Desiderio et al. 2007, 2010); others on the environmental aspects (Colli and De Ascentis 2003; Buccolini et al. 2003; Scalella and Di Francesco 2004). Etiope et al. (2003, 2007) and Etiope & Klusman (2010) discussed geochemical data relating to MVs in Italy, including the ‘Cenerone-Pineto’ MV, particularly in relation to gas emissions. The authors found gas emissions related to N2, CO2 and CH4. Helium concentration from the gas emissions (Etiope et al. 2007; Scalella and Di Francesco 2004) suggests the presence of deep open fractures.

Rainone et al. (2015) carried out a multidisciplinary study on the ‘Cenerone-Pineto’ MV which includes geology, hydrogeology and, particularly, hydrochemistry and geophysics: through two-dimensional (2D) reflection seismic, two-dimensional (2D) and three-dimensional (3D) ERT surveys, the authors defined the underground structure of this MV.

This paper presents new results from synthetic modelling of 2D and 3D Electrical Resistivity Tomography (ERT) field data collected by Rainone et al. (2015) at the ‘Cenerone-Pineto’ MV (Central Italy). This study arises from the need to replicate and validate the observed structural features from the ERT and seismic reflection data gathered by Rainone et al. (2015). Synthetic modelling results aim to provide key indicators for the objective and quantitative interpretation of the measured resistivity pattern, which is not always straightforward (Van Schoor 2002), and to verify the underground structure of this MV proposed by Rainone et al. (2015).

1.3 The resistivity signature of mud fluids

Electrical Resistivity Tomography (ERT) provides an imaging of the subsurface electrical resistivity pattern and allows subsurface structures to be identified. Its theory (e.g., Dahlin and Loke 1998; Daily and Owen 1991; Loke et al. 2003) and application (e.g., Cassiani et al. 2009; Guérin et al. 2004; Kuras et al. 2009; Mazzini et al. 2021b; Ritz et al. 1999; Torrese 2022; Torrese et al. 2022) are well-documented in geophysical literature. ERT surveys have been successfully used for environmental purposes due to their capacity to identify resistivity changes arising from water-content variation (Torrese 2020), the migration of contaminants (Kean et al. 1987; Pilla and Torrese 2022; Torrese and Pilla 2021; Van et al. 1991), lithological variation or to MV occurrence (Accaino et al. 2007; Chang et al. 2010; Grassi et al. 2022; Imposa et al. 2016; Rainone et al. 2015). ERT surveys are appropriate for the study of the ‘Cenerone-Pineto’ MV because of the large gap (more than 4.5 times) that occurs between the electrical resistivity of the typical groundwater and the waters that emerge from the MV.

Changes in electrical resistivity correlate with variations in fluid conductivity, water saturation, solid material and porosity. These may be used to map contaminants, geological structure, open fractures, stratigraphic units, groundwater and cavities. Electrical resistivity methods are techniques well-suited to the identification and mapping of low resistivity structures such as mineralized water saturated deposits, low-middle resistivity structures such as freshwater saturated fine-grained deposits (e.g., clayey deposits), as well as high resistivity structures like freshwater saturated coarse-grained deposits (e.g., sandy deposits).

The principal geophysical targets of this study are mud fluids, i.e., mud-rich mineralized waters with electrical conductivity (EC) constantly above 8 mS/cm, which are low-resistivity units. Therefore, mud fluid migration pathways are correlated with structures with minimum resistivity values.

The resistivity signature of different hydrogeological bodies is dependent on two factors, first, the size of the target (e.g., a mud-fluid conduit) in relation to its depth; second, the contrast between the resistivity of the target and that of the host rock. The amplitude of resistivity anomalies is an inverse function of the distance between the measurement points and the target: the deeper the target, the lower the method’s reliability in identifying it.

The host rock around the target is frequently disturbed, i.e., fractured and jointed, unconsolidated, weathered or altered. The physical and chemical condition of the surrounding rock at times generates a larger, anomalous bulk volume compared to the target itself. Accordingly, the actual geophysical size of a target depends on the geologic setting, and it is usually larger than its true size (Chalikakis et al. 2011; Torrese 2020; Torrese et al. 2021). This, however, typically means that the target is easier to detect.

Further, sedimentary covers that show a high-resistivity gap compared to the host rock around the target may substantially alter the underlying target’s geophysical response. Moreover, very low resistivity values at the ground surface may, locally, reduce the penetration depth of the electrical flow pattern and disturb the target’s apparent geometry. In fact, given that the current flow pattern depends on the layering, a ground with lower resistivity in the upper layer reduces the penetration depth because the current flow tends to concentrate into the upper layer with lower resistivity.

The investigation depth, the vertical and horizontal resolutions of electrical resistivity surveys, are linked to: i) the electrode spacing, ii) the configuration array, iii) the quadrupole sequence, iv) the signal-to-noise ratio (S/N) ratio, v) the contrast between the target’s resistivity and the surrounding rock and/or background resistivity.

For the reasons stated, of the existing geophysical methods, ERT is particularly well-suited to identifying mud fluid migration pathways. ERT detection of these targets can be demanding and replicating and validating these structures through synthetic modelling may be necessary.

2 The ‘CENERONE-PINETO’ MV and its hydrogeological and hydrochemical setting

The dome of the ‘Cenerone-Pineto’ MV measures 15 m × 10 m. It is approximately 2 m high with a crater measuring 2.5 m in diameter (Figs. 1b and 2). Fluids and solids (mud, cold brine, gas) are ejected through effusive activity from the crater (Fig. 2). Recent studies (Colli and De Ascentis, 2003) have mentioned five MVs in the province of Teramo characterized by effusive activity. There is no evidence of paroxysm events occurred in the past. Over the last decades, the emissions have undergone an intense reduction, as regards the number of MVs present, as well as their activity. The ‘Cenerone-Pineto’ MV is the only one that continues to have effusive activity, of the known MVs in the area. Mud emissions are estimated to be a few liters per day.

Fig. 2
figure 2

Pictures and description of the ‘Cenerone-Pineto’ mud volcano. a Overall view. b Mud pond

Desiderio and Rusi (2004) have dated the clay-bearing sediment that emerges from the MV to the Pleistocene. These sediments come from the Mutignano formation (Bigi et al. 1997; Casnedi et al. 1977; Casnedi 1983; Crescenti et al 1980). They are pelitic deposits which consist of over-compacted clays with thin silty-sandy layers, 6°-10° ENE dipping. Below the MV area there are 2000–2500 m thick foredeep pelitic deposits (ENI 1972).

The area is characterized by a compressional tectonic setting with folds and blind thrusts buried under Middle Pliocene deposits (Fig. 3). The main tectonic element in the area is the Campomare NNW-SSE buried thrust (Casnedi 1983; Crescenti et al. 1980; Finetti et al. 2005; Scisciani and Montefalcone 2005), which is located approximately below the MV area and represents a structural high at a depth of approximately 550 m (Fig. 3).

Fig. 3
figure 3

Geological setting of the ‘Cenerone-Pineto’ mud volcano area, schematic and simplified cross section which shows the MV location referred to the main tectonic element in the area, the Campomare NNW-SSE buried thrust

In the MV area, the pelitic deposits are found at a depth of approximately 10 m below clayey-silty eluvial-colluvial deposits that crop out. Ditches and streams in the area follow fault systems oriented NW–SE and NE-SW (Bigi et al. 1997).

The piezometric head is 13.6 m above sea level (asl) in correspondence to the MV, that is approximately 2 m above ground level, whereas the main fresh groundwater level is about 11 m asl as observation boreholes surrounding the MV revealed.

The waters that emerge from the MV’s pond have electrical conductivity (EC) constantly above 8 mS/cm and oxidation–reduction (REDOX) potential always negative. Groundwater from surrounding observation boreholes shows completely different values with EC and REDOX potential ranging between 0.94 and 2.29 mS/cm and between + 89 and + 149 mV, respectively. The chemistry of the waters that emerge from the MV’s pond is sodium chloride with ionic formula Na, K, Mg, Cl, HCO3, SO4, with a low sulphate content and pH always slightly alkaline. These chemical features are similar to those of groundwaters found in the alluvial and foredeep deposits in the Apennines that emerge at surface in correspondence with structural highs and thrusts or with tectonic lineaments that crop out or are buried (Nanni and Vivalda 1999; Desiderio and Rusi 2004; Desiderio et al. 2005, 2010). These structural discontinuities represent preferential flow paths for basal saline waters and facilitate the flow towards the surface (Torrese and Pilla 2021) aided by N2, CO2 and CH4 gases.

The isotopic data for 18O and 2He from the waters that emerge from the MV’s pond is less negative than that of the groundwater found in the valleys’ alluvial deposits and the carbonate structures of the Apennines. The waters from the MV are not at all connected to recent recharge waters as δ18O and δ2H values from the MV’s groundwater are more positive than the water that surrounds it (Rainone et al. 2015).

All these hydrochemical observations suggest that there is no connection between mud fluids that rise to the surface through the clayey lithological sequence and fresh groundwater, according to the “mud dome” type A model (Kopf 2002).

3 The subsurface structure of the mV based on seismic and ERT field data

Based on 2D-3D ERT and shear wave reflection seismic (Fig. 4), Rainone et al. (2015) suggested:

-the presence of a shallow mud fluid reservoir (MFR), not located beneath the crater, but laterally offset (Fig. 5); its low resistivity suggests that it is congruent with the occurrence of a high permeability layer saturated with mud-rich mineralized waters;

-a connection between the MFR and the MV’s crater through a narrow feeder channel (FC) (Fig. 5); the subvertical migration of mud fluids seems to occur along an open fracture (Fig. 5); indeed, the shear wave reflection seismic indicates that the shallow root of the MV is affected by upward bending of layers (Fig. 10 in Rainone et al. (2015)) and shows a main central open fracture that coincides with the MV’s crater.

-evidence of an open fracture situated just below the reservoir (fracture below the reservoir, FBR) where the shear wave velocity model suggests the occurrence of a low-velocity, fractured zone (Fig. 10 in Rainone et al. (2015)); this open fracture would feed the MFR with mud fluids (Fig. 5).

Details on the interpretative conceptual model (Fig. 5) based on the geophysical results along with hydrochemical results are provided by Rainone et al. (2015).

Fig. 4
figure 4

Survey map (from Google Earth™). Location of 2D and 3D ERT surveys together with some electrode settings; spread and surface profile of the stack cross-section of the reflection seismic survey (RFL) together with the first and last station settings

Fig. 5
figure 5

a Conceptual model of the ‘Cenerone-Pineto’ mud volcano based on the joint interpretation of shear wave reflection seismic and, 2D and 3D ERT field data models according to Rainone et al. (2015). The conceptual model shows the occurrence of (i) the mud fluid reservoir (MFR) located approximately 17–33 m below ground level with a 50 m offset towards ENE of the crater, (ii) the open fracture located below the mud fluid reservoir (FBR) within the over-compacted clay deposits that would feed the MFR with mud fluids, (iii) the fracture located just below the crater (FBC) that would not feed the crater with mud fluids and (iv) the narrow feeder channel (FC) that would connect the MFR to the crater. The purple arrow indicates mud fluid migration pathways

4 Need for synthetic modelling

Synthetic modelling is a numerical simulation aimed at evaluating the capacity of the implemented model to detect and delineate predefined targets. Synthetic modelling is the inversion of synthetic (simulated) datasets generated by forward modelling. Unlike inverse modelling, which aims to estimate the model parameters of a given dataset, forward modelling is aimed at calculating the response of a given model. Forward modelling produces datasets (synthetic datasets) that can be compared to the actual observation (field datasets). Comparison between synthetic and field datasets can be made by comparing synthetic and field subsurface models generated by inversion of synthetic and field datasets, respectively.

This study arises from the need to verify whether the geophysical targets detected by 2D and 3D ERT field data are true underground features. This depends upon the response of the employed experimental setup to the degree of heterogeneity in the subsoil and the eventuality that the resistivity pattern obtained could be affected by measure and inversion artifacts.

The targets the study is focused on, identified by the joint interpretation of shear wave reflection seismic and 2D-3D ERT, were defined by Rainone et al. (2015) in their conceptual model shown in Fig. 5. Geophysical targets to be verified by synthetic modelling play a key role in defining mud fluid migration pathways of the ‘Cenerone-Pineto’ MV. The predefined targets the synthetic modelling has been focused on are (Fig. 5):

  • the occurrence of the mud fluid reservoir (MFR), located approximately 17–33 m below ground level and with a 50 m offset towards ENE of the crater;

  • to rule out of the current occurrence an open fracture located just below the crater (FBC) that would feed the MFR with mud fluids;

  • the occurrence of the open fracture located below the mud fluid reservoir (FBR) within the over-compacted clay deposits that would feed the MFR with mud fluids;

  • the occurrence and width of the narrow feeder channel (FC) that would connect the MFR to the crater.

Two-dimensional synthetic modelling is aimed at replicating, simulating and verifying the entire section of subsoil already investigated by 2D ERT field data, essential in defining the mud fluid migration pathways of the MV down to a depth of, approximately, 40 m. Three-dimensional modelling is focused on a localized and restricted volume of subsoil that crosses the MV, specifically the crater, the FC and the uppermost part of the MFR: this portion of subsoil, already investigated by 3D ERT field data, hosts shallow, target structures that play a key role in accurately defining the shallow mud fluid migration pathways of the MV down to a depth of approximately 20 m.

Detection of small sized structures such as the FC, or deep structures like the FBR, according to the implemented conceptual model (Fig. 5), is often challenging. Furthermore, shallow lateral heterogeneity (e.g., sandy-clayey units and freshwater filled clays in Fig. 6a) as well as uneven topography such as the crater rim, may mask the true resistivity pattern of the subsoil model and generate artifacts and depth or lateral displacement of the geophysical targets (Torrese et al. 2021; Torrese 2022). For all these reasons, synthetic modelling plays a key role in validating the implemented conceptual model inherent mud fluid migration pathways of the ‘Cenerone-Pineto’ MV.

Fig. 6
figure 6

Inverse model from 2D ERT field data showing a a full-scale resistivity range and b a low resistivity range extraction related to the occurrence of mud fluids and clayey-silty deposits filled with mineralized water; a narrow feeder channel (FC) seems to connect the MFR to the crater

5 Material and methods

5.1 Overview

Two-dimensional (2D) and three-dimensional (3D) synthetic modelling was performed to exactly replicate the experimental setup involved in the field dataset collection of 2D ERT and 3D ERT surveys (Fig. 4), respectively, and simulate the measured resistivity patterns (Figs. 6, 10, 11a). The steps involved in the modelling were:

  • designing of 2D and 3D conceptual models which simulate the measured resistivity patterns;

  • conversion of the 2D and 3D conceptual models into 2D (Figs. 7a, 8a, 9a) and 3D (Figs. 12a, 13a1–e1) synthetic models;

Fig. 7
figure 7

a 2D ERT synthetic model which simulates a 5 m wide narrow feeder channel (FC) and a 25 m wide open fracture located below the mud fluid reservoir (FBR); b the related inverse model which clearly shows the FC that would connect the MFR to the crater but doesn’t clearly show the presence of the FBR that would feed the MFR with mud fluids. The width of the FC and the FBR must be considered as the width of the fracture zones through which mud fluids spread during their uprising

Fig. 8
figure 8

a 2D ERT synthetic model which simulates a 5 m wide narrow feeder channel (FC) and a 25 m wide open fracture located just below the crater (FBC); b the related inverse model which clearly shows the FC that would connect the MFR to the crater but does not show the presence of the FBC that would feed the crater with mud fluids directly from below. The width of the FC and the FBC must be considered as the width of the fracture zones through which mud fluids spread during their uprising

Fig. 9
figure 9

a 2D ERT synthetic model which simulates a 5 m wide narrow feeder channel (FC) and a 25 m wide open fracture located below the mud fluid reservoir (FBR); b the related inverse model obtained by combined (hybrid) arrays: the combination of different arrays delivered better detectability and imaging of FB and FBR (to compare with Fig. 7). The width of the FC and the FBR must be considered as the width of the fracture zones through which mud fluids spread during their uprising

  • generation of 2D and 3D synthetic datasets by forward modelling through the application of the synthetic models, according to the experimental setup involved in the field dataset collection;

  • inversion of the synthetic datasets to generate realistic, synthetic 2D (Figs. 7b, 8b, 9b) and 3D (Figs. 11b, c, 12b, 13a2–e2) inverse resistivity models to be compared with field 2D and 3D inverse resistivity models.

Subsurface models obtained from synthetic datasets were compared with subsurface models obtained from field datasets. Both synthetic and field subsurface models were generated by an iterative inversion process, composed of forward and inverse modelling.

5.2 Design of conceptual and synthetic models

The conceptual model (Fig. 5), which was implemented based on the joint interpretation of 2D and 3D ERT models and shear wave reflection seismic, according to that proposed by Rainone et al. (2015), was converted into 2D (Figs. 7a, 8a, 9a) and 3D (Figs. 12a, 13a1-e1) synthetic models aimed at replicating the measured resistivity pattern revealed by 2D ERT (Fig. 6) and 3D ERT (Figs. 10, 11a). These synthetic models are inhomogeneous models discretized in different resistivity units that represent the geometry of the MV along with the mud fluid migration pathway, as well as the geometry of the features underlying and surrounding the MV.

Fig. 10
figure 10

Inverse model from 3D ERT field data. Low resistivity range extraction (volume) is related to the occurrence of mud fluids compared to full-scale range extracted slices; a narrow feeder channel (FC) seems to connect the MFR to the crater

Fig. 11
figure 11

Inverse models from 3D ERT field and synthetic data: comparison between the slices y = 20 m extracted by a the inverse model from field data, b the inverse model from synthetic data which simulates a 3 m wide narrow feeder channel (FC) and c the inverse model from synthetic data which simulates a 23 m wide FC. Synthetic data modelling revealed that the FC appears to be detectable if its width is at least 3 m; however, it is well defined only if its width is at least 23 m

These models include the geophysical targets that play a key role in defining the mud fluid migration pathways:

  • the mud fluid reservoir (MFR);

  • the open fracture located just below the crater (FBC);

  • the open fracture located below the mud fluid reservoir (FBR).

  • the narrow feeder channel (FC).

Various synthetic models were implemented to test and verify different widths of the FBC, the FBR and the FC. This allowed an assessment of the detectability developed by the experimental setup in identifying these geophysical targets.

The conceptual model (Fig. 5) was converted into 2D (Figs. 7a, 8a, 9a) and 3D (Figs. 12a, 13a1–e1) synthetic models by assigning the measured electrical resistivity value to each unit.

5.3 Forward and inverse modelling

Same electrode arrays (Wenner-Schlumberger (WS) and dipole–dipole (DD) for 2D ERT and 3D ERT, respectively), quadrupole sequence (596 and 270 quadrupoles for 2D ERT and 3D ERT, respectively) and experimental layout (48 electrodes spaced 5 m apart and a snake grid comprised of 12 × 4 electrodes spaced 10 m apart, both along the X and Y axes for 2D ERT and 3D ERT, respectively), which were employed in the field dataset collection (Rainone et al. 2015), were used in 2D ERT and 3D ERT synthetic modelling. A test involving 2D joint inversion of synthetic data from combined arrays (Wenner-Schlumberger, dipole–dipole, pole-dipole and pole-pole) was also performed to evaluate any advantages in detecting demanding targets such as the FC.

2D ERT profile intersects the MV’s crater at the 12th electrode (Fig. 4). Electrode 18 of 3D ERT was placed on the MV’s crater (Fig. 4).

The computational domain was discretized into tetrahedral cells in both forward and inverse modelling. A grid 235 m × 5 m in size, with a maximum depth of ≈40 m, was designed for 2D ERT modelling and a grid 110 m × 30 m in size, with a maximum depth of ≈20 m, was designed for 3D ERT. Foreground region was discretized using a 2.5 m element size for 2D ERT and a 5 m element size for 3D ERT. The background region was discretized using an increasing element size towards the outermost parts of the domain, according to the sequence 1×, 1×, 2×, 4× and 8× the foreground element size.

ERTLab Solver (by Multi-Phase Technologies LLC, Geostudi Astier srl), based on tetrahedral Finite Element Modelling (FEM), was used for both forward and inverse modelling. Two-dimensional and 3D synthetic datasets were generated by forward modelling through the application of the synthetic models.

The forward modelling was performed using mixed-boundary conditions (Dirichlet–Neumann) and a tolerance (stop criterion) of 1.0E-7 for a Symmetric Successive Over-Relaxation Conjugate Gradient (SSORCG) iterative solver. Of the various algorithms (e.g., Sensitivity Conjugate Gradients (SCG), Incomplete Cholesky Conjugate Gradient (ICCG), Preconditioned Conjugate Gradient Method (CGPC)), some of which may be even more robust (e.g., Levenberg–Marquardt (Levenberg 1944; Marquardt 1963; Nelles 2001)), the SSORCG iterative solver was used as it provides fast convergence by decreasing the spectral condition number of the matrix (Bing and Greenhalgh 2001; Spitzer 1995) and is particularly efficient and accurate in the modelling of complex resistivity structures (meter-sized structures in this study) requiring a computational domain discretized in a large number of nodes (Spitzer 1995).

Data inversion was based on a least squares smoothness constrained approach (LaBrecque et al. 1996). Noise was appropriately managed using a data-weighting algorithm (Morelli and LaBrecque 1996) that allows the adaptive changes of the variance matrix after each iteration for those data points that are poorly fit by the model: after each inversion iteration, the weights contained in the variance matrix are decreased for those data which are poorly fit by the reconstructed model; the weights of data with good fits are kept at their original values. In this algorithm, a variation of the least-absolute deviations method described by Mosteller and Tukey (1977) was implemented.

Since field data always contain noise in practice, the obtained synthetic datasets were corrupted with 1% random noise to provide realistic results and were inverted to generate inverse resistivity models. This value was considered appropriate to simulate the low level of field noise in the MV’s area, also considering the depth, volume and resistivity of the geophysical target (Liu et al. 2020).

The inverse modelling was performed using a maximum number of internal inverse Preconditioned Conjugate Gradient (PCG) iterations of 15 and a tolerance (stop criterion) for inverse PCG iterations of 0.001. The amount of roughness from one iteration to the next was controlled in order to assess maximum layering along the X and Z directions for 2D ERT and along the X, Y, and Z directions for 3D ERT: a low value of reweight constant (0.1) was set with the objective of generating maximum heterogeneity.

Both 2D (Figs. 7, 8, 9) and 3D (Figs. 10, 11, 12, 13) inverse resistivity models were chosen with a criterion based on the achievement of a minimum data residual (misfit error). Misfit in terms of chi-squared error and root mean square (RMS) error of inverted resistivity data are shown in Table 1.

Table 1 Misfit in terms of chi-squared error and root mean square (RMS) error of inverted resistivity data

6 Results

6.1 2D synthetic inverse model

6.1.1 Evidence of fluid migration pathways

Two-dimensional synthetic data modelling allowed the clear recognition of the uprising of mud fluids, down to a depth of 39 m; this is characterized by areas of low resistivity, ranging between 1.9 and 5.6 Ω m.

The resistivity pattern revealed by the 2D synthetic inverse model (Fig. 7b) is highly analogous with the resistivity pattern shown by the 2D field inverse model (Fig. 6a) and, according to the interpretation proposed by Rainone et al. (2015), consistent with the presence of:

  • a low resistivity (1.88–3.47 Ω m) shallow body representing the crater’s root: here mud fluids stagnate in deposits of a clayey-silty matrix;

  • a low resistivity (3.4–3.9 Ω m) body located 50 m to the ENE of the MV at a depth between 17 and 33 m representing the mud fluid reservoir (MFR) or mud chamber; this is a high permeability layer that is saturated with mud fluids isolated from the adjacent groundwater system and representing the last phase of mud accumulation before final emission;

  • a narrow feeder channel (FC) located 11–31 m below ground level that would connect the MFR to the crater. It is worth highlighting that unlike the crater’s root and the MFR which appear very well defined by the resistivity contrast, the presence of the FC is not very clear on the 2D field inverse model with full-scale resistivity range (Fig. 6a): in fact it has a change in resistivity small compared to the host rock and appears with gradient outlines (i.e., diffuse or not well defined) due to its low width/depth ratio (in this case the higher accuracy obtained from the 3D field inverse model is decisive). The FC appears better exposed on the 2D field inverse model using a resistivity range extraction as shown in Fig. 6b in addition to the 2D synthetic inverse model. This suggests that mud fluids flow sub-horizontally along the MFR at depth and then flow upward towards the MV’s crater; the subvertical migration of mud fluids appears to be very localized as it occurs along an open fracture, as also suggested by the shear wave reflection seismic.

6.1.2 Detectability of fluid migration pathways

Two-dimensional synthetic data modelling confirmed ERT’s difficulty in identifying the open fracture located below the mud fluid reservoir (FBR) along which mud fluids would rise up to the reservoir, supplying it. As field data models, synthetic models seem to indicate its presence, even if these do not allow a clear identification of the open fracture at such depths. Synthetic data models suggest that the supposed low resistivity band related to this open fracture, and apparent on field data models, should have a minimum width of 25 m (Fig. 7) to be able to generate a perturbation in the resistivity distribution at a depth of approximately 35 m.

Otherwise, the possibility that ERT is able to identify the open fracture located just below the crater (FBC) that would supply it directly from below can be excluded. ERT is unable to detect this structure at depth even if the low resistivity band shows a width of 25 m (Fig. 8).

Two-dimensional synthetic data modelling revealed that the narrow feeder channel (FC) appears to be detectable if its width is at least 3 m. However, this channel is well defined only if its width is at least 5 m, as shown in Fig. 7.

In addition, 2D synthetic data modelling revealed that the joint inversion of data from combined arrays (Wenner-Schlumberger, dipole–dipole, pole-dipole and pole-pole) delivered better detectability and imaging of the narrow feeder channel (FC) and the open fracture located below the mud fluid reservoir (FBR), as shown in Fig. 9.

6.2 3D synthetic inverse model

6.2.1 Evidence of fluid migration pathways

Three-dimensional synthetic data modelling allowed the clear recognition of the uprising of mud fluids at shallow depths, down to a depth of 19 m, characterized by areas of very low resistivity, ranging between 1.6 and 4.2 Ω m.

The resistivity pattern revealed by the 3D synthetic inverse model (Figs. 11b,c, 12, 13) is highly analogous with the resistivity pattern shown by the 3D field inverse model (Figs. 10, 11a) and, according to the interpretation proposed by Rainone et al. (2015), consistent with the presence of:

  • a low resistivity (1.88–3.47 Ω m) shallow body representing the crater’s root: here mud fluids stagnate in deposits of a clayey-silty matrix;

  • a low resistivity (3.4–3.9 Ω m) body located 50 m to the ENE of the MV at a depth of 17 m representing the mud fluid reservoir (MFR); this is a high permeability layer that is saturated with mud fluids isolated from the adjacent groundwater system;

  • a narrow feeder channel (FC) located 11 m below ground level, which would connect the MFR to the crater; this appears well exposed on the resistivity range extraction shown in Fig. 10.

Fig. 12
figure 12

Low resistivity range extraction from a 3D ERT synthetic model which simulates a 23 m wide narrow feeder channel (FC) and b the related inverse model: synthetic data modelling revealed that the FC appears to be well defined only if its width is at least 23 m. The width of the FC must be considered as the width of the fracture zones through which mud fluids spread during their uprising

Fig. 13
figure 13

Low resistivity range extraction from a 3D ERT synthetic model which simulates a 23 m wide narrow feeder channel (FC) and b the related inverse model, compared to full scale range extracted slices; synthetic data modelling revealed that the FC appears to be well defined only if its width is at least 23 m. The width of the FC must be considered as the width of the fracture zones through which mud fluids spread during their uprising

This confirms that the shallow, subvertical migration of mud fluids appears to be very localized as it occurs along an open fracture, as also suggested by the shear wave reflection seismic. The narrow feeder channel (FC) is apparent down to the maximum depth of the 3D model, 19 m. Due to its size, the 3D model is able to show only the uppermost and westernmost part of the mud fluid reservoir (MFR).

6.2.2 Detectability of fluid migration pathways

Three-dimensional synthetic data modelling revealed that the narrow feeder channel (FC) appears to be detectable if its width is at least 3 m (Fig. 11b), as also revealed by 2D modelling. However, this channel is well defined only if its width measures at least 23 m, as shown in Figs. 11c, 12 and 13. This is apparent by the comparison between slices and low resistivity ranges extracted by the synthetic inverse models Fig. 13.

7 Discussion

This is one of the few studies focused on 2D and 3D detailed mapping and imaging of shallow fluid flow pathways in terrestrial MVs. It may provide further insights into the shallow structure of terrestrial MV given that, in literature, the internal structure of mud volcanoes is usually only schematically reconstructed and represented on a larger scale (Deville et al. 2003; Stewart and Davies 2006), often making a comparative analysis difficult.

Synthetic modelling allowed target structures already noted in ERT field data and consistent with seismic data to be recognized as true resistivity features (Rainone et al. 2015). Modelling results allowed the ruling out of factors such as uneven topography at the crater rim, electrode configuration and shallow vertical and lateral heterogeneities, having masked the true resistivity pattern of the subsoil and generated artifacts and/or depth or lateral displacement of the geophysical targets. The resistivity imaging derived from ERT field data seems unaffected by these disturbing factors and measure or inversion artifacts.

The resistivity pattern revealed by 2D and 3D synthetic inverse models are highly analogous with the resistivity pattern shown by 2D and 3D field inverse models. This suggests that synthetic models, forward and inversion procedures accurately replicated and simulated the actual resistivity pattern of the internal structure of the ‘Cenerone-Pineto’ MV and the experimental setup involved in the field dataset collection.

The structural features replicated and validated in this study, according to the conceptual model proposed by Rainone et al. (2015) (Fig. 5), play a key role in defining the internal structure of the ‘Cenerone-Pineto’ MV: (i) a mud fluid reservoir (MFR) or mud chamber (a high permeability layer saturated with mud-rich mineralized waters), (ii) a narrow feeder channel (FC) that would connect the MFR to the crater, (ii) the crater’s root.

One of the most significant findings of this study is that synthetic data models confirmed the 50 m offset between the mud fluid reservoir (MFR), at a depth between 17 and 33 m, and the crater. This means that the MV is not supplied with mud fluids directly from below as would be the case with an uprising of deep fluid along a near-vertical open fracture. These findings are consistent with the conceptual model proposed by Rainone et al. (2015). Other studies have also shown similar evidence about the internal structure of MVs.

Zeyen et al. (2011), through their pseudo-3D resistivity imaging of the Mangapakeha MV in New Zealand, found that the mud chamber at a depth of about 20 m is not well correlated to the distribution of the vents: at a depth of 8–20 m, the lowest resistivities do not correspond to the areas of vent concentration and, surprisingly, certainly not to the position of the main vent. Especially below the main vent, they found no connection between the upper and the lower mud reservoir.

Blake et al. (2021), through their pseudo-3D resistivity imaging of the Piparo MV on the island of Trinidad, found that the maximum extent of mud fluids is not located beneath the main crater; conversely, they found that the central conduit occluded and an alternative lateral pathway was formed that created the current main crater location. This evidence from Blake et al. (2021) is also in agreement with field data from this study revealing further possible emissions of mud fluids at surface before agricultural activities and other anthropogenic activities destroyed these features at the ‘Cenerone-Pineto’ MV site.

Accaino et al. (2007), through their tomographic inversion of 3D seismic data and 2D resistivity imaging of Nirano MVs in northern Italy, found shallow mud fluid reservoirs not correlated to the distribution of any MV observed on the surface.

Even the depth of the mud fluid reservoir (MFR) of the ‘Cenerone-Pineto’ MV, that is 17 m, is similar to that of other terrestrial MVs. Accaino et al. (2007) recognized at Nirano MVs in northern Italy a mud chamber located at about 25 m below the MV surface—only 8 m deeper than the MFR of the ‘Cenerone-Pineto’ MV. Zeyen et al. (2011) found at the Mangapakeha MV in New Zealand a mud chamber at a depth of approximately 20 m—only 3 m deeper than the MFR of the ‘Cenerone-Pineto’ MV.

Compared with other examples in literature (Martinelli and Dadomo 2005) generally characterized by deeper reservoirs, the MFR observed in the ‘Cenerone-Pineto’ MV represents the last phase of mud accumulation before the final emission.

Another significant finding of this study is that the narrow feeder channel (FC) observed in the ‘Cenerone-Pineto’ MV, which would connect the MFR to the crater, should have a width of at least 3–5 m, as revealed by 2D and 3D synthetic data models. The fact that the FC appears to be better exposed and defined by the 3D field data model (Fig. 10) than 3D synthetic data models (Fig. 13) suggests that the actual width of this conduit may be greater than 3–5 m. This is also supported by the comparison between 2D and the more accurate 3D synthetic data models. As expected, the FC is better exposed and defined in its nearly horizontal portion where it is overlain by more resistive (coarser) shallow deposits, granting a greater penetration depth of the electrical flow pattern in this area. In fact, given that the current flow pattern depends on the layering, a ground with higher resistivity in the upper layer increases the penetration depth because the current flow tends to concentrate into the deeper layer with lower resistivity. It is worth highlighting that Zeyen et al. (2011) speculated on the presence of a shallow, narrow, up to 3 m wide mud conduit at the Mangapakeha MV in New Zealand, a width very similar to that found for the ‘Cenerone-Pineto’ MV.

As field data models, synthetic models seem to indicate the presence of a possible deep open fracture located below the mud fluid reservoir (FBR) that would feed the MFR with mud fluids. Although synthetic models do not allow a clear identification of the open fracture at such depths, they suggest that it should have a minimum width of 25 m. A lower thickness of the fractured area would not justify the perturbation in the resistivity distribution related to this open fracture and apparent on the field data models. ERT’s difficulty in identifying this open fracture is due to its great depth and the occurrence of low resistivity mud fluids within the reservoir (MFR) located above it; this reduces the penetration depth of the electrical flow pattern in this area, preventing a clear detection of this target by ERT (Szalai et al. 2009; Torrese 2020; Torrese et al. 2021). It is worth pointing out that although the minimum width of 25 m assessed in this study for the FBR may appear too high and unrealistic a value, Zeyen et al. (2011) also found that mud fluids uprise within some 20 m-wide chimneys before reaching the surface at the Mangapakeha MV in New Zealand, which has a size comparable to the ‘Cenerone-Pineto’ MV.

Despite this evidence, it is appropriate to make some considerations regarding the width of the narrow feeder channel (FC) and the fracture located below the mud fluid reservoir (FBR). A width of at least 3–5 m for the FC and a width of at least 25 m for the FBR is way too large to be geologically realistic. Fractures with such openings are not mechanically feasible (especially in the geological context of the ‘Cenerone-Pineto’ MV which consists of clayey-silty and silty-sandy deposits) unless one is referring to a fracture zone in which both host rock and fractures are present or the damage zone of a fault. Fractures or conduits associate with MVs structures should not have an opening larger than a few decimeters to be mechanically feasible. This is the reason why the width of the FC and the FBR found by this study should be considered as the width of the fracture zones through which mud fluids spread during their uprising.

Another hypothesis that would justify the too large width of these structures is that the line along which the ERT sampled the subsurface is subparallel to the fractures. In this case ERT would detect an apparent rather than actual width. Actually, the ERT profile was carried out along an ENE–WSW direction (Fig. 4) that is roughly orthogonal to the major direction for fracture zones in the Apennines, NW–SE, the direction of the main axis of the buried anticlines in the area (as well as being the only direction that allows an adequate length of the electrode array to reach adequate investigation depths). However this direction is very similar to that of a fault system in the area which is oriented NE-SW (Bigi et al. 1997).

This hypothesis, however, is less plausible than the previous one: in fact, both the FC and the FBR were clearly detected by the seismic reflection profile (Fig. 5 of this paper, Fig. 11 of Rainone et al. 2015) which was carried out along an ENE–WSW direction (Fig. 4), i.e. a direction very similar to that of the ERT. If the line along which the seismic reflection profile sampled the subsurface had been subparallel to the fractures, the seismic profile would not have detected any fracture.

Therefore, the hypothesis that the width of the FC and the FBR must be considered as the width of the fracture zones through which mud fluids spread during their uprising, represents the currently most probable hypothesis.

Furthermore, it can be ruled out that the apparent width of these structures may be in some way related to unconsidered factors affecting detectability and resolution: the detectability provided by the experimental (electrode spacing, configuration array, quadrupole sequence) and inversion setup has been strictly verified by synthetic modelling; the misfit in terms of chi-squared error (1) for the inversion final iterations suggests that inverse models are free of artefacts due to an inversion over-fit (chi-squared error << number of measures) or excessive smoothing due to an inversion under-fit (chi-squared error >> number of measures).

In addition, with regard to the seismic reflection profile, experimental setup (acquisition spread geometry, geophone and source spacing, offset distance, folding, signal bandwith, signal-to-noise ratio (SNR)) and processing sequence provided suitable lateral detectability to identify meter sized structures at shallow depths (Torrese et al. 2019).

Synthetic models suggest that ERT is unable to detect any deep open fracture located just below the crater (FBC) that would supply it directly from below. This is due to its great depth and the occurrence of very low resistivity mud fluids within the shallow crater’s root and the narrow feeder channel (FC) located above it; as for the FBR, but to an even greater extent, this strongly reduces the penetration depth of the electrical flow pattern in this area, preventing any detection of this target by ERT. However, on the basis of both field and synthetic ERT results, the presence of the FBC cannot categorically be ruled out. It is worth highlighting that multiple upwelling pathways may have occurred over time in the ‘Cenerone-Pineto’ MV, as well as in other MVs. In fact, as found by Blake et al. (2021), the upwelling conduits can become occluded and alternative lateral pathways formed, even creating new crater locations. Furthermore, as found by Mazzini et al. (2021a), upwelling conduits can also become occluded by powerful eruptions. Although this is not the case of the ‘Cenerone-Pineto’ MV as there is no evidence of paroxysm events occurred in the past, it is worth underlining that the self-sealing of subsurface conduits after eruptions could allow the pressure to recharge more quickly, leading to more frequent, more energetic and larger eruptions (Gattuso et al. 2021; Mazzini et al. 2021a).

The complex internal structures of the ‘Cenerone-Pineto’ MV, as revealed from this study, is consistent with the structure of other terrestrial MVs, as shown by Blake et al. (2021) and Zeyen et al. (2011). Even the resistivity range measured for the mud deposits of the ‘Cenerone-Pineto’ MV, that is 2–6 Ω m, is similar to those of other terrestrial MVs. Blake et al. (2021) found < 2.5 Ω m for the Piparo MV on the island of Trinidad; Zeyen et al. (2011) found 5–10 Ω m for the Mangapakeha MV in New Zealand. These low resistivity values are corroborated with the resistivity of the mud deposits of other areas, such as the Wushanting MV area (Southwestern Taiwan) (Ping-Yu et al. 2011), Lokbatan and Bozdag-Kobijsky MVs (Eastern Azerbaijan) (Scholte et al. 2003), Sidoarjo MV (East Java, Indonesia) (Istadi et al. 2009) and Mutnovsky volcano (South Kamchatka, Russia) (Bessonova et al. 2012). It is worth highlighting that the bulk resistivity range found for the mud deposits of the ‘Cenerone-Pineto’ MV (2–6 Ω m) is given by the resistivity of the waters that emerge from the MV's pond (constantly below 1.25 Ω m, above 8 mS/cm in EC) slightly increased (clayey-silty matrix for the crater’s root) or significantly increased (silty-sandy matrix for the MFR and the FC) by the matrix resistivity.

This study recommends the use of synthetic modelling for the ERT investigation of the internal structure of terrestrial MVs. In such a complex geological setting, which includes difficult geological targets, numerical simulation aims at evaluating the capacity of the implemented model to detect and delineate predefined targets and, therefore, allows the interpretation of field data to be guided. As revealed by this study, difficult geological targets in such a complex geological setting can be represented by narrow feeder channels or deep open fractures located in zones affected by reduced electrical current penetration, preventing a clear detection of the target. In this study, 2D synthetic modelling was aimed at replicating, simulating and verifying the entire investigated area, while 3D synthetic modelling was focused on a localized and restricted volume of subsoil that hosts shallow, target structures that play a key role in defining the mechanisms and processes involved in mud volcanism.

Synthetic modelling also revealed that the experimental setup used in this study (electrodes spaced 5 m apart for 2D ERT and spaced 10 m apart for 3D ERT) provided a fairly accurate definition of the main hydrogeological units and fracture zones through which mud fluids spread during their uprising. Given the resolution provided by the experimental and inversion setup used, the geophysical models obtained, although unable to detect isolated, disconnected fractures, were able to identify hydraulically effective fracture zones (permeability network for fluids) a few meters in size, at least 3–5 m. Obviously, detectability falls as depth increases.

Given the geological/hydrogeological complexity of the internal structure of MVs with both vertical (e.g., FBR, part of the FC) and horizontal (e.g., MFR, part of the FC) structures, the use of a single array in ERT surveys, such as the Wenner-Schlumberger array, which ensure high vertical resolution and signal amplitude, or the dipole–dipole array, which provides enhanced lateral resolution, may not ensure the best imaging though they provide valuable results (e.g., Zeyen et al. 2011). Synthetic modelling (compare Fig. 7-single array with Fig. 9-combined arrays) found that the merging of different arrays differing in vertical and lateral resolution delivered better detectability and imaging and, therefore, delivered a more accurate inverse model of the MV according to Szalai et al. (2009), Torrese (2020) and Torrese and Pilla (2021).

Despite the limitations characteristic of synthetic data modelling, this study’s results validate the structural features observed by ERT field data and consistent with seismic data (Rainone et al. 2015). They are consistent with the hypothesis that mud fluids are likely to rise along deep open fractures (FBR) within the over-compacted clay deposits located to ENE of the MV. The occurrence of an over-pressurized mud fluids reservoir (MFR) or mud chamber was revealed at 17–33 m below ground level. This high permeability unit is not located directly below the crater of the MV but appears to be located 50 m to the ENE. Mud fluids would flow sub-horizontally along it at depth and before flowing upward toward the MV's crater along a narrow feeder channel (FC) located 11–31 m below ground level. The subvertical migration of mud fluids appears to be very localized as it occurs along an open fracture. Mud fluids can only reach the surface at the MV's crater due to the localized thinning of the drift deposits.

Note that in this paper the author assumes that the subvertical migration of mud fluids occurs along open fractures rather than along faults. As the two features may have very different hydraulic properties, this assumption is justified by the fact that, at this scale and in the absence of specific structural geology studies, no evidence of displacement or rock-mass movements has been observed.

Findings from this study can be extended to other terrestrial MVs located in similar geological contexts and may be helpful in survey designing and interpreting stages of future 2D, 3D and time-lapse ERT surveys aimed at detecting difficult geological targets that represent crucial elements for the definition of the mechanisms and processes involved in mud volcanism.

8 Conclusions

This paper presents new results from ERT synthetic modelling for the investigation of the shallow, internal structure of terrestrial mud volcanoes. This is one of the few studies focused on 2D and 3D detailed mapping and imaging of shallow fluid flow pathways in terrestrial mud volcanoes. Furthermore, recommendations for the use of the ERT method for the stratigraphic mapping of terrestrial mud volcanoes are provided, especially for the detection of difficult geological targets such as narrow feeder channels or deep open fractures located in zones affected by reduced electrical current penetration, preventing a clear detection of the target.

This study revealed the internal structure of the ‘Cenerone-Pineto’ mud volcano (Central Italy) which includes: a possible pressurized mud fluid reservoir or mud chamber approximately 17–33 m below ground level fed by a deeper subvertical open fracture; the apparent 50 m offset between this mud fluid reservoir and the crater; a narrow, at least 3–5 m, subvertical shallow feeder channel along which the mud fluid appears to be flowing towards the surface from the reservoir up to the crater. These results are consistent with the observed structural features already noted in ERT and seismic field data collected at the ‘Cenerone-Pineto’ mud volcano.

The main results of the study are the evidences about the presence of a mud chamber, which represents the last phase of mud accumulation before final emission, not located beneath the crater, but laterally offset, as well as the presence of a narrow, shallow feeder channel, i.e., a much more complex structure than one would expect. This means that the mud volcano is not supplied with mud fluids directly from below as would be the case with an uprising of deep fluid along a near-vertical open fracture and that shallow mud fluid reservoirs may be not correlated to the distribution of any mud volcano observed on the surface. These results provide further insights into the shallow structure of terrestrial MV given that, in literature, the internal structure of mud volcanoes is usually only schematically reconstructed and represented on a larger scale, often making a comparative analysis difficult.

Findings from this study can be extended to other terrestrial MVs located in similar geological contexts and may be helpful in survey designing and interpreting stages of future 2D, 3D and time-lapse ERT surveys aimed at detecting difficult geological targets that represent crucial elements for the definition of the mechanisms and processes involved in mud volcanism.