Crustal Structure of the Seismogenic Volume of the 2010–2014 Pollino (Italy) Seismic Sequence From 3D P- and S-Wave Tomographic Images

A tomographic analysis of Mt. Pollino area (Italy) has been performed using earthquakes recorded in the area during an intense seismic sequence that occurred between 2010 and 2014. 870 local earthquakes with magnitude ranging from 1.8 to 5.0 were selected considering the number of recording stations, the signal quality, and the hypocenter distribution. P- and S-wave arrival times were manually picked and used to compute 3D velocity models through tomographic seismic inversion. The resulting 3D distributions of V_P and V_S are characterized by high resolution in the central part of the investigated area and from surface to about 10 km below sea level. The aim of the work is to obtain high-quality tomographic images to correlate with the main lithological units that characterize the study area. The results will be important to enhance the seismic hazard assessment of this complex tectonic region. These images show the ductile Apennine platform (V_P = 5.3 km/s) overlaying the brittle Apulian platform (V_P = 6.0 km/s) at depth of around 5 km. The central sector of the area shows a clear fold and thrust interface. Along this structure, most of the seismicity occurred, including the strongest event of the sequence (M_W 5.0). High V_P (>6.8 km/s) and high V_P/V_S (>1.9) patterns, intersecting the southern edge of this western seismogenic volume, have been interpreted as water saturated rocks, in agreement with similar geological context in the Apennines. These fluids could have played a role in nucleation and development of the seismic sequence. A recent study revealed the occurrence of clusters of earthquakes with similar waveforms along the same seismogenic volume. The hypocenters of these cluster events have been compared with the events re-located in this work. Jointly, they depict a 10 km × 4 km fault plane, NW-SE oriented, deepening towards SW with a dip angle of 40–45°. Instead, the volume of seismicity responsible for the M_L 4.3 earthquake developed as a mainshock-aftershock sequence, occurring entirely within the average-to-low V_P/V_S Apennine platform. Our results agree with other independent geophysical analyses carried out in this area, and they could significantly improve the actual knowledge of the main lithologic units of this complex tectonic area.

Introduction

Mt. Pollino area (Italy) is located at the transition between the southern edge of the Apennines chain and the Calabrian subduction forearc (Figure 1A). It is one of the most structurally complex tectonic regions in the intra Apennine seismogenic belt (Michetti et al., 1997; Cinti et al., 2002). The area is considered a seismic gap in the Apennine chain because no earthquakes of magnitude M > 6 occurred in historical time, in striking contrast with the nearby north and south sectors where the occurrence of many strong earthquakes is well documented (Figure 1B). During the last centuries, the seismicity of Pollino area has been of moderate magnitude (M_L < 6.0, CPTI15 database; Rovida et al., 2016), but paleo-seismological studies have recognized active faults south-east of the study area that generated at least two earthquakes of magnitude larger than M_L 6.5 in the last 10,000 years (Cinti et al., 2002; Cheloni et al., 2017). The geology of the area is characterized by the overlay of a low strain-rate extensional deformation occurred during Early Pleistocene (D’Agostino et al., 2011) over pre-existing strike-slip Late Pliocene–Early Pleistocene structures (Ghisetti and Vezzani, 1982; Schiattarella, 1998). The current NE-SW extensional stress regime (Totaro et al., 2016; Napolitano et al., 2021) is driven by the simultaneous retrograde S-Eward motion of the subduction zone and by the opening of the Tyrrhenian back-arc basin (D’Agostino et al., 2011). The result of such lithospheric evolution is a complex structure characterized by W- and E-dipping sub-parallel normal faults in the upper crust (above 15 km depth) and by strike-slip faults below 15 km depth (Totaro et al., 2015; Brozzetti et al., 2017; Ferranti et al., 2017; Napolitano et al., 2021).

FIGURE 1

FIGURE 1. (A) Map of the study area, located between the Southern Apennines and the Calabrian subduction fore arc (CFA). Seismic events analyzed are shown with red circles. Yellow stars represent the two largest events of the 2010–2014 Pollino sequence. Beach balls represent their focal mechanism. Permanent and temporary seismic stations operated by INGV (blue triangle) and UNICAL (cyan circle) and shared between INGV and UNICAL (blue triangle in cyan circle) and by GFZ (green diamond) are shown. (B) The historical earthquakes of M_L>6 (CPTI15, Rovida et al., 2016) are shown with blue squares. The target area, clearly a seismic gap, is enclosed in a magenta square and the orange circle is the location of the swarm area. (C) The local stress field, modified after Napolitano et al. (2021).

The recent seismicity of the study area has been characterized by low and moderate magnitude earthquakes. A seismic sequence followed the M_W 5.6 earthquake that occurred on September 9, 1998, in the Mercure basin (Guerra et al., 2005; Brozzetti et al., 2009). Between 2010 and 2014 a seismic sequence of more than 10,000, low-to-moderate earthquakes occurred in the study area, on the western side of the Pollino Massif, with a time-varying rate of seismicity (Passarelli et al., 2015; Totaro et al., 2015). The sequence developed within two seismogenic volumes and showed a combination of both swarm-like and classic mainshock-aftershock sequence (Passarelli et al., 2015; Napolitano et al., 2021). Two main events occurred in 2012: a M_L 4.3 in the eastern volume on 28 May, and a M_L 5.0 (the largest magnitude earthquake of the sequence) in the western volume on 25 October (yellow stars in Figure 1A). These events were characterized by normal fault mechanism (black beach balls in Figure 1A) in perfect agreement with the local extensional stress field of the area (Figure 1C, modified after Napolitano et al., 2021). In addition, a slow earthquake, located at the western seismogenic volume of the sequence, occurred from mid-2012 to mid-2013, releasing aseismically an amount of strain energy equivalent to M_W 5.5 (Cheloni et al., 2017).

Many studies have been carried out to assess the seismic hazard of the Pollino area seismic gap. Earthquakes from the 2010–2014 sequence have been analyzed by many authors to provide new insights about the structure of the seismogenic faults (Totaro et al., 2015; Brozzetti et al., 2017) and to study site effects (Napolitano et al., 2018). Part of them have been used by Napolitano et al. (2021) to relocate with high precision clusters of events with similar waveforms to compute focal mechanisms of small earthquakes, to obtain a very local stress field, and to investigate the swarm evolution in space and time. Geological-structural, seismic, and geophysical data have been combined by Ferranti et al. (2017) to characterize the seismotectonic frame of the study area from surface to deep crust. The presence of fluid-filled connected faults has been inferred from 2D imaging of scattering and absorption contributions to the total attenuation of high frequency coda waves (Napolitano et al., 2020; Sketsiou et al., 2020). The coda-normalization method has been applied to perform a 3D attenuation image of fluid storage and tectonic interactions across the Pollino fault network (Sketsiou et al., 2021).

A detailed study of the elastic properties of the crust is useful to assess the seismic hazard of complex tectonic areas. When many well recorded earthquakes are available, this goal can be achieved by passive seismic tomography, based on the use of P- and S-wave travel times. The technique is particularly effective in areas characterized by appropriate distributions of hypocenters and seismic stations. Until now, the Pollino area has been only marginally included in large scale tomographic studies, whose main target area was the Southern Apennines—Calabrian arc region (Barberi et al., 2004; Orecchio et al., 2011). These papers analyzed earthquakes occurred before the 2010–2014 Pollino seismic sequence. Barberi et al. (2004) used earthquakes recorded from 1978 to 2001 and 84 artificial sources, while Orecchio et al. (2011) used data recorded between 1981 and 2008. In both papers, the authors performed the inversion using a horizontal and vertical grid spacing of 40 and 10 km, respectively. More recently, Totaro et al. (2014) used earthquakes occurred from 1981 to October 2012, partially including events of the 2010–2014 seismic sequence, to perform a seismic tomography at regional scale. Also in this case, the authors focused on the likely transition from the thinner Tyrrhenian crust to the thicker crust beneath Calabria. Finally, Pastori et al. (2021) provided a new local 1D velocity model with the goal of performing refined 1D and double-difference locations using the events from the last Pollino sequence.

Our work aims to provide a highly detailed 3D tomographic image of the shallow crust in the Pollino area, taking advantage of the 2010–2014 seismic sequence. We mean to enhance our knowledge about the elastic properties of the seismogenic volume responsible for the long-lasting sequence and to investigate the possible role of fluids in triggering the sequence, as suggested by many authors (Passarelli et al., 2015; Napolitano et al., 2020; Sketsiou et al., 2020; Sketsiou et al., 2021). To reach this goal, we preliminarily checked the quality of the available waveforms recorded at the permanent and temporary stations installed in the area between 2010 and 2014. For each of these waveforms, we manually picked P- and S-wave arrivals and checked their quality. A preliminary location of these events has been obtained using one of the large scale models available for Southern Italy. Using this reliable dataset, through a statistical procedure, we obtained a 1D model, used as starting model of the 3D tomography which simultaneously inverts the arrival time of body waves for both velocity model and earthquake location parameters. Synthetic tests were carried out to validate the results. Then, the results were compared with other independent geological and geophysical information available for this area. A reliable interpretation of the properties of the shallow crust beneath the Pollino area might be an essential starting point for future more detailed analysis about seismogenesis and evolution of the 2010–2014 sequence.

Data

We selected 870 earthquakes occurred in the Pollino area, within a volume of 100 × 120 × 25 km³ (Figure 1A). These events, characterized by magnitude between 1.8 and 5.0, were recorded between January 30, 2010, and November 25, 2014. We manually picked P- and S-wave arrival times of each waveform to guarantee the robustness and reliability of the analyses to be carried out. We collected 9981 P- and 6862 S-wave arrival times recorded by 29 seismic stations (Figure 1A) operated in the area by three institutes: Università della Calabria (Unical), Istituto Nazionale di Geofisica e Vulcanologia (INGV), and Deutsches GeoForschungsZentrum of Potsdam (GFZ) (Passarelli et al., 2012; Margheriti et al., 2013). Weights from 0 to 4 were assigned to each picking based on its accuracy (Table 1).

TABLE 1

TABLE 1. Weights associated to P- and S-wave pickings. Increasing weight is related to decreasing accuracy in manual picks of the direct wave arrival times.

We applied the modified Wadati diagram method to evaluate the picking accuracy of the P- and S-wave arrival times (Chatelain, 1978). Furthermore, the Wadati diagram provides an estimate of the average V_P/V_S value for the study area. The average V_P/V_S, estimated considering all points within one standard deviation from the best fit line, results to be $\frac{V_P}{V_s}=1.786\pm 0.004$ (slope of the red line in Figure 2A). The V_P/V_S obtained from this analysis is very similar to those of previous studies characterized by different scales and datasets (Barberi et al., 2004; Piana Agostinetti and Amato, 2009; Ferranti et al., 2017; Pastori et al., 2021).

FIGURE 2

FIGURE 2. (A) Modified Wadati diagram of the arrival time data set of this study. Black points represent retained data to perform the fit (within one standard deviation); the fitted line is shown in red. The slope corresponds to V_P/V_S = 1.786 ± 0.004. (B) Histogram representing the quality factor (q_f, Michele et al., 2019), i.e., the four classes of the earthquake location (A-, B-, C-, and D-classes) used in this work. Only A-class location and B-class location are obtained. (C) The seven location parameters used to obtain q_f: Err (km) is the error on horizontal (red) and vertical (green) hypocenter locations, azimuthal gap in degrees (Gap), root mean square error in seconds (RMS), number of phases (Nph), distance between the probability density function (pdf) expected value and its maximum likelihood (LocDist), and radius of the volume described by the scatter points radius of pdf (PDFRad).

As a first step, we performed a preliminary hypocenter location of the initial dataset using the probabilistic approach implemented in the software NonLinLoc (Lomax et al., 2000). This preliminary location takes advantage of the 1D P-wave velocity model by Barberi et al. (2004) and of the V_P/V_S estimated from the modified Wadati diagram (Figure 2A). We established a lower threshold of 6 direct wave arrival times per event including at least 2 S-wave picks and a root-mean-square (RMS) smaller than 0.6 s to retain location results. We combined seven location parameters (Figure 2C) in order to evaluate the earthquake location quality factor following the formula proposed by Michele et al. (2019) (q_f, Figure 2B). We obtained 779 A-class and 75 B-class locations, which confirm the high quality of the dataset and the reliability of signal pre-processing. The quality assessment allowed us to define a high-quality dataset comprising 854 well located earthquakes, for which 9643 P- and 6649 S-wave arrival times are available, useful for the following analysis.

Tomographic Method and Inversion Strategy

In this section, we describe the procedure that was applied to obtain the 3D tomographic inversion and explain how we tuned the selected parameters. We first describe the selection of the most reliable 1D model used as starting model for the tomographic inversion, and then we focus on the inversion strategy.

A reliable final 3D velocity model must necessarily start from a trustworthy 1D model (Kissling et al., 1994). We chose the P-wave velocity model by Barberi et al. (2004) as the starting reference 1D velocity model for the preliminary hypocenter location and data pre-processing. Then, we applied a statistical procedure that allows to find the most reliable local reference 1D velocity model to be used as an input in the 3D tomography of the study area (Vanorio et al., 2005; De Matteis et al., 2010). Following this procedure, we generated 120 1D P-wave models characterized by velocity values randomly chosen in the range ±500 m/s around the value of the reference model for each layer (Figure 3A). For each of these 120 models, we performed a preliminary relocation of the dataset using NonLinLoc (Lomax et al., 2000) and then we performed a 3D tomographic inversion using the same velocity model (Latorre et al., 2004). We then compared the RMS time residuals reduction between the first and 12th iteration (red and light blue points in Figure 3B). Even though the starting RMS values span from 0.17 s up to 0.44 s, all the final solutions reach RMS time residuals smaller than 0.1 s. We selected the 47 velocity models whose starting RMS time residuals were smaller than 0.25 s. The velocity at each layer was averaged and a new local reference 1D P-wave velocity model was obtained (solid red line in Figure 3A). This new velocity model turns out to be faster (0.23 km/s on average) than the starting model by Barberi et al. (2004) for almost any layers and agrees with the most recent 1D model obtained by Pastori et al. (2021). The S-wave velocity model has been obtained dividing the P-wave velocity by V_P/V_S = 1.786, estimated from the modified Wadati diagram (Figure 2A).

FIGURE 3

FIGURE 3. (A) 120 randomly generated initial 1D P-wave velocity models (grey lines). Blue line represents the starting P-wave velocity model by Barberi et al. (2004), black lines are the upper and lower velocity bounds to perform the random selection. Red line is the new local P-wave velocity model achieved in this work and used as starting velocity model for the 3D velocity tomography. (B) Starting (light and dark red) and final (cyan) RMS achieved from the 3D inversion of the 120 randomly generated velocity models. Darker red points below 0.25 s represent the 47 averaged models to achieve the solid red line in panel (A). The two stars mark the starting and final RMS of the final inversion, using the new velocity model (red line in panel A). (C,D) P-wave (green) and S-wave (violet) residuals of the starting (lighter colors) and final (darker colors) inversion. Triangles on the x-axes of the corresponding color show the mean value of each distribution.

We used a linearized, iterative travel time tomography which simultaneously inverts the arrival time of body waves for both velocity model and earthquake locations parameters (Latorre et al., 2004). This method has been used in many applications in tectonic regions and volcanic and geothermal areas using data from passive (e.g., Vanorio et al., 2005; Amoroso et al., 2014; Amoroso et al., 2018) and active surveys (De Landro et al., 2017; De Landro et al., 2020). The hypocenter and station distributions allow us to set the mesh spacing of the inversion grid to 5 × 5 × 1 km³. A finer grid of 0.5 × 0.5 × 0.5 km³ was set to compute the theoretical travel times at each station through a finite-difference solution of the eikonal equation (Podvin and Lecomte, 1991). Travel times were recalculated by numerical integration of the slowness along the ray path. The model parameters (velocities, hypocenter coordinates, and origin time) were inverted using the least squares root (LSQR) method of Paige and Saunders (1982).

The damping parameter was set to 0.3 through the empirical approach of the trade-off curve between the model parameters (P- and S-wave model, spatial hypocenter location, and origin time) and the data variances L-curve (Supplementary Figure S1). The velocity distribution in a continuous medium is described by a trilinear interpolating function based on a grid of regularly spaced nodes.

The derivative weight sum (DWS), that provides a measure of the density of the rays within each grid cell, and the checkerboard test were performed to mark the well-resolved areas of the map. Supplementary Figures S2, S3 in the Supplementary Material show the checkerboard tests performed on both V_P (maximum amplitude 400 m/s) and V_S (maximum amplitude 300 m/s) models. Areas well resolved by the checkerboard test strongly match with the DWS contour line. Hence, the tomographic figures are presented with the DWS outline to show the well resolved parts of the map.

Results

Travel-time residual distributions for P- and S-wave velocity models are represented in Figure 3C and Figure 3D, respectively. The mean value of both final distributions (3D model) is closer to zero than the initial distributions (1D model), particularly that of the P-wave. The RMS is reduced by 59% after 12 iterations, from the initial value of 0.175 s to final 0.072 s (yellow stars in Figure 3B). These results demonstrate that the 3D velocity model corresponds to a significant improvement in the earthquake location quality.

Figure 4 shows horizontal slices at 2, 3, 4, 5, 6, 7, 8, and 10 km depth of the V_P (left panels) and V_S (right panels) models. The corresponding maps of the derived V_P/V_S ratio, ranging from 1.48 to 2.08, whose well resolved area is contoured by the DWS corresponding to the V_S model, are shown in Supplementary Figure S4. The sectors resolved by the DWS are contoured by a dark gray line. The V_P values span from 3.6 to 6.9 km/s, while V_S values are between 2.1 and 3.9 km/s. Each horizontal slice shows dark grey circles that represent seismic events above and below 500 m from the selected depth (labeled in the upper right corner of the maps). Due to the event and station configurations, our results are well resolved mainly in the central part of the study area. Starting from depth of 2 km, the well resolved area becomes broader down to 7 km depth and then it shrinks at larger depths. Thereby, we do not interpret anomalies located at depth greater than 10 km.

FIGURE 4

FIGURE 4. Maps of V_P(A) and V_S(B) computed at 8 depths: 2, 3, 4, 5, 6, 7, 8, and 10 km. Color scales are shown at the bottom. Areas crossed by no rays are shown in white. Area resolved by DWS are contoured by the dark grey line. Labels from a to e indicate the velocity anomalies identified and interpreted in the text.

Between 2 and 5 km depth, a strongly homogeneous low V_P (4.5–4.6 km/s) anomaly appears (see label a in Figure 4). This low velocity anomaly becomes narrower with increasing depth until it disappears at depth of 6 km in the SW sector of the resolved area (15.9°E, 39.85°N). At the same depths, V_S maps show a low velocity anomaly (see label b in Figure 4) in the SW corner of the resolved area (15.9°E, 39.87°N). A low V_P (5.2–5.4 km/s) and low V_S (2.9–3.1 km/s) anomaly appears in the central sector of the map. It deepens from NW to SE between 5 and 8 km depth, following the direction of the western cloud of seismicity (label c in Figure 4). In addition, the eastern cloud of seismicity (see label d in Figure 4) is marked by low V_P (5.2–5.4 km/s) and medium-to-low V_S (3.0–3.3 km/s).

Between 6 and 8 km depth, the central part and eastern side of the map are characterized by a high V_P anomaly (6.7–6.9 km/s, see label e in Figure 4). At 6 km depth, this anomaly intercepts the western seismogenic volume of the sequence, characterizing it almost completely. Average values of V_S have been found in the same location (3.3–3.5 km/s). The western side of the western cloud of seismicity (see label f in Figure 4) is marked by a high V_P (6.5–6.6 km/s) and the highest values of V_S (3.7–3.9 km/s). These values of V_S also appear at 8 km depth, but no V_P anomalies are evident at the same depth.

Figures 5, 6 show the E-W and N-S V_P (left panels) and V_P/V_S (right panels) cross sections, respectively. These cross sections are 50 km wide and include depths between 2 km a.s.l and 13 km b.s.l. They are also projected on the map located on top of Figure 5 with red lines. The range of values of V_P/V_S spans from 1.45 to 2.1. All the anomalies in terms of the achieved P- and S-wave velocity models are retrieved and shown at depth also in terms of V_P/V_S.

FIGURE 5

FIGURE 5. Vertical cross sections (AA’,BB’,CC’,DD’) along the EW direction of V_P and V_P/Vs models. The cross sections are shown on map in the uppermost panel of this figure, containing also seismic stations (cyan diamonds) and the events re-located during the 3D tomographic inversion (grey circles). The cross sections are 50 km long and 13 km deep. The contour of the DWS is shown using the dark grey line. The transition between IAP and AP (at average of 5.3 km/s) is shown by a black thick line. Labels from a to e indicate the velocity anomalies identified and interpreted in the text.

FIGURE 6

FIGURE 6. Vertical cross sections (EE’,FF’,GG’,HH’) along the SN direction of V_P and V_P/Vs models. The cross sections are shown on map in the uppermost panel of Figure 5. The cross sections are 50 km long and 13 km deep. The contour of the DWS is shown using the dark grey line. The transition between IAP and AP (at average of 5.3 km/s) is shown by a black thick line. Labels from a to e indicate the velocity anomalies identified and interpreted in the text.

The low V_P anomaly found between 2 and 5 km depth results also in a low V_P/V_S anomaly (1.55–1.65, see label a in Figure 5 AA’, BB’; Figure 6). The shape of this low V_P and low V_P/V_S anomaly in the E-W direction adheres to the slope of the southern edge of the western fault plane that was the most active during the 2010–2014 seismic sequence. From both directions, it is quite clear that as the low V_P, low V_P/V_S body emerges at around 5 km depth, the seismicity of the western cluster ends. This seismicity basically occurred in the deeper layer characterized by high V_p (V_P = 6.0 km/s on average), as shown by the horizontal slices. We also find out that this high V_P leads to high V_P/V_S values (V_P/V_S > 1.9, see label e in Figure 5 BB’, CC’ and Figure 6 FF’, GG’ HH’).

On the other hand, the resulting V_P/V_S for the medium surrounding the eastern cluster of seismicity, responsible for the M_L 4.3 earthquake, shows average values (1.65–1.75, see label d in Figure 5, BB’; Figure 6, HH’).

Figure 6 EE’ and FF’ show that the low V_P anomaly, deepening from 5 to 8 km depth, is also characterized by low V_P/V_S (1.65–1.75, see label c). The slope of this low V_P anomaly follows part of the western earthquake group until about half of the cloud of hypocenters (Figure 6, FF’).

Finally, the western side of the western cluster of seismicity is marked by high V_P/V_S (>1.85, see label f in Figure 5 BB’ and Figure 6, EE’).

Discussion

The tomographic images obtained from our analysis show a strong V_P increase at an average depth of 5 km, as marked by a black line in Figures 5, 6 (left panels). We interpret this velocity contrast as the lithological transition from the shallower ductile mesozoic sequence of the Apennine platform (AP) to the deeper limestones of the inner Apulian platform (IAP) (Improta et al., 2000; Improta et al., 2017). This interpretation is confirmed by the structural sections inferred for the Mercure–Pollino region (Ferranti et al., 2017). The shallow AP formations include both the areas labeled a and d (Figures 5, 6). Usual V_P values of the Apennine formations found in other similar geological setting of the Southern Apennines are around 5.3–5.5 km/s (Improta et al., 2000; Improta et al., 2017). These values have been found on the eastern side of the study area (label d in Figures 4–6) and decrease moving toward west (label a in Figures 4–6), as well as V_P/V_S. This might suggest a different rheologic behavior moving within the same AP. The substantial absence of seismicity within the volume characterized by anomaly a suggests that this ductile and highly fractured volume acts as a barrier, preventing the seismicity generated by the western volume to propagate upward. Combination of low V_P and low-to-medium V_P/V_S indicates absence of significant fluid saturation within the AP crustal layer. It is noticeable, however, that a further reduction of V_P/V_S comes out around 0–2 km depth (Figure 6 FF’), but further analysis is needed to fully explain the nature of this anomaly.

The central sector of the study area shows the shallower AP deepening between 5 and 8 km depth in the NNW-SSE direction within the deeper IAP (label c in Figure 5 CC’, Figure 6 FF’). We interpret the peculiar shape of this interface as a portion of the AP–IAP fold and thrust interface. From the geological point of view, this result is not surprising, relying on the fact that Southern Italy is floored by the fold and thrust belt in the western and central part (Ferranti et al., 2017). More specifically, regarding the Pollino area, it has been found that during the late contraction stage (Late Pliocene to Early Pleistocene), the thick Apulia crust was involved in the thrust system and then reversely reactivated as normal faults in the actual extensional regime, responsible for the seismic hazard of the area (Ferranti et al., 2017; Filice and Seeber, 2019). For the first time, a three-dimensional tomographic image shows with impressive detail what until now had been hypothesized by geological sections. This low V_P-low V_P/V_S extension of the AP within the IAP (label c) agrees also with the low attenuation volume observed by Sketsiou et al. (2021) between the western and eastern cluster of seismicity of the 2010–2014 Pollino sequence.

High V_P and high V_P/V_S anomalies characterize the IAP west (label f) and east (label e) of the western seismogenic volume of the 2010–2014 Pollino seismic swarm (Figures 4–6). The western f anomaly is confined between 5 and 8 km depth and seems not to reach the western volume of seismicity. Instead, the eastern e anomaly, mostly confined between 7 and 10 km depth, intersects the southern edge of the western seismogenic volume at 6–7 km depth. This means that fluids might be likely involved in the nucleation and development of the western swarm-like sequence. Indeed, similar values of V_P/V_S have been associated with high pore fluid pressure and fluid migration along active normal fault systems in the central and Southern Apennines (Chiarabba et al., 2009; Amoroso et al., 2014; Baccheschi et al., 2019; Chiarabba et al., 2020). These values of V_P and V_P/V_S provide a Poisson ratio in the range 0.3–0.4 which are marker of water saturated rocks in the investigated crustal volume (Dvorkin et al., 1999; Dutta, 2002). This resolution boundary, due to the shallow seismicity in the swarm area, does not allow us to understand whether the high V_P/V_S extends to larger depths. We speculate about a possible deeper extension to the SE of the swarm area, relying on the high attenuation volume on the SE side of the map at 12.5 km depth found by Sketsiou et al. (2021), perfectly matching with the eastern high V_P/V_S found in this work. In the same paper, another high attenuation body has been found west of the study area that can be related to the f anomaly found in this work. The SE sector of the Pollino area has also been identified by Napolitano et al. (2020) as high scattering and high absorption area in a broad frequency band. These patterns depicted in their 2D maps have been interpreted as a complex fluid-filled connected fault network, where the largest and deepest structures are located SE of the target area (Cinti et al., 2002; Cheloni et al., 2017). Sketsiou et al. (2020) using a novel technique confirm the results in terms of absorption anomalies location.

The two seismogenic volumes involved in the 2010–2014 Pollino sequence belong to two different geological and geophysical settings. The eastern seismogenic volume is entirely included in the low-to-average V_P/V_S AP, while the western volume is included within the high V_P and high V_P/V_S brittle IAP, along the IAP-AP fold and thrust interface (Figure 5, BB’; Figure 6, FF’, HH’). Different geological settings featured the seismicity of these two volumes. The former sequence developed as a classic mainshock-aftershock sequence, following the M_L 4.3 earthquake on May 28, 2012, while the latter developed as combination of swarm-like and mainshock-aftershock sequence (Passarelli et al., 2015). Regarding the western seismogenic volume, Napolitano et al. (2021) applied the master-slave relative location technique (Got et al., 1994) to clusters of earthquakes with local magnitude ranging between 0.6 and 2.7. Their results allowed to reconstruct with high detail the most active fault of the western seismogenic volume. In this work, we computed the absolute location of the cluster reference events using the 3D velocity model obtained from our tomographic inversion and added the hypocenters of the 423 earthquakes analyzed by Napolitano et al. (2021). Figure 7 shows the hypocenters of relative located cluster events (red circles) and the hypocenters of events used in this paper (grey circles). These hypocenters lay within the same seismogenic volume and depict a 10 km × 4 km wide fault system, located between 4 and 8 km depth, orientated NW-SE, and deepening towards SW with a dip of around 40–45°, as found by Napolitano et al. (2021). We argue that the structural barrier represented by the AP did not allow the rupture to further widen, preventing it from generating events larger than M5.0. The likely presence of pressurized fluids stopped by a structural barrier is also suggested by Pastori et al. (2021).

FIGURE 7

FIGURE 7. Hypocenter locations from three different views (panels A–C) of the western seismogenic volume achieved in this work (grey circles) and cluster events (red circles) by Napolitano et al. (2020), achieved relocating the reference event using the new 3D velocity model.

Conclusion

The detailed knowledge of the crustal structure of an area considered as seismic gap, but affected in the recent past by a seismic sequence, is fundamental to evaluate and mitigate the seismic risk. The widest seismic gap in Italy, the Pollino area, recently experienced a long-lasting seismic sequence, characterized by both mainshock-aftershock and swarm behavior. Despite the fact that many studies have been recently carried out about this area, a detail of the crustal structure involved in the sequence and of the surrounding area has not been provided yet. The use of manually picked P- and S-wave arrival times of local earthquakes allowed to obtain three-dimensional images of the crustal volume involved in the 2010–2014 Pollino seismic sequence and near surrounding crust with unprecedented detail. The images, well resolved within the first 10 km depth, point out the main lithological units and lateral variations of velocity within this seismic gap area. Although the crustal structure of the Pollino area is mainly composed of two large overlapping geological units, namely, the Apennine platform and the underlying inner Apulian platform, our tomographic results highlight the complexity achieved by them in this area. Between 5 and 8 km depth in the central sector of the study area, the Apennine platform deepens in the NNW-SSE direction within the inner Apulian platform, forming what we interpreted as a fold and thrust, a common feature in the Southern Apennines framework. Our images clearly show that the two seismogenic volumes of the last Pollino sequence occurred in two different geological settings. The eastern cluster of seismicity, characterized by an aftershock sequence following the M_L 4.3 earthquake, entirely developed within the low V_P and low-to-average V_P/V_S Apennine platform. The larger western cluster occurred at the boundary of the Apennine platform–inner Apulian platform fold and thrust. The southern edge of this cluster is characterized by high V_P and high V_P/V_S, whose values provide a Poisson’s ratio indicative of water saturated rocks.

In addition, cluster of similar waveforms, analyzed in a previous paper by Napolitano et al. (2021), have been re-located with the new 3D model and added to the hypocenters located in this work. These repeaters occurred along or very near the same western fault plane well depicted in this work: a 10 km × 4 km wide structure, located between 4 and 8 km depth, orientated NW-SE, and deepening towards SW with a dip of around 40–45°. The structural barrier represented by Apennine platform stopped the propagation of the seismicity toward north and upward.

Data Availability Statement

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

Author Contributions

FN, AG, and ML performed the data processing and the manual picking of the arrival times. FN and OA conceived this work and worked at the tomographic inversion. FN, SG, OA, ML, and PC were involved in model analysis and interpretation of results. FN, OA, SG, and ML wrote the paper. PC provided funding for this work and coordinated the working group. All co-authors were involved in the review of the manuscript.

Funding

This work has been supported by PRIN-MATISSE (20177EPPN2) project funded by the Italian Ministry of Education and Research.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/feart.2021.735340/full#supplementary-material

References

Amoroso, O., Ascione, A., Mazzoli, S., Virieux, J., and Zollo, A. (2014). Seismic Imaging of a Fluid Storage in the Actively Extending Apennine Mountain belt, Southern Italy. Geophys. Res. Lett. 41 (11), 3802–3809. doi:10.1002/2014GL060070

CrossRef Full Text | Google Scholar

Amoroso, O., Festa, G., Bruno, P. P., D’Auria, L., De Landro, G., Di Fiore, V., et al. (2018). Integrated Tomographic Methods for Seismic Imaging and Monitoring of Volcanic Caldera Structures and Geothermal Areas. J. Appl. Geophys. 156, 16–30. doi:10.1016/j.jappgeo.2017.11.012

CrossRef Full Text | Google Scholar

Baccheschi, P., De Gori, P., Villani, F., Trippetta, F., and Chiarabba, C. (2019). The Preparatory Phase of the Mw 6.1 2009 L’Aquila (Italy) normal Faulting Earthquake Traced by Foreshock Time-Lapse Tomography. Geology 48 (1), 49–55. doi:10.1130/G46618.1

CrossRef Full Text | Google Scholar

Barberi, G., Cosentino, M. T., Gervasi, A., Guerra, I., Neri, G., and Orecchio, B. (2004). Crustal Seismic Tomography in the Calabrian Arc Region, South Italy. Phys. Earth Planet. Interiors 147 (4), 297–314. doi:10.1016/j.pepi.2004.04.005

CrossRef Full Text | Google Scholar

Brozzetti, F., Cirillo, D., de Nardis, R., Cardinali, M., Lavecchia, G., Orecchio, B., et al. (2017). Newly Identified Active Faults in the Pollino Seismic gap, Southern Italy, and Their Seismotectonic Significance. J. Struct. Geology. 94, 13–31. doi:10.1016/j.jsg.2016.10.005

CrossRef Full Text | Google Scholar

Brozzetti, F., Lavecchia, G., Mancini, G., Milana, G., and Cardinali, M. (2009). Analysis of the 9 September 1998 Mw 5.6 Mercure Earthquake Sequence (Southern Apennines, Italy): a Multidisciplinary Approach. Tectonophysics 476, 210–225. doi:10.1016/j.tecto.2008.12.007

CrossRef Full Text | Google Scholar

Chatelain, J. L. (1978). Etude fine de la sismicité en zone de collision continentale au moyen d’un réseau de stations portables: la région Hindu-Kush Pamir. Thèse de doctorat de 3ème cycle Université scientifique et médicale de Grenoble. Berlin: Springer, 219.

Google Scholar

Cheloni, D., D’Agostino, N., Selvaggi, G., Avallone, A., Fornaro, G., Giuliani, R., et al. (2017). Aseismic Transient during the 2010-2014 Seismic Swarm: Evidence for Longer Recurrence of M ≥ 6.5 Earthquakes in the Pollino gap (Southern Italy)? Sci. Rep. 7 (576), 576. doi:10.1038/s41598-017-00649-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Chiarabba, C., Buttinelli, M., Cattaneo, M., and De Gori, P. (2020). Large Earthquakes Driven by Fluid Overpressure: The Apennines normal Faulting System Case. Tectonics 39, e2019TC006014. doi:10.1029/2019TC006014

CrossRef Full Text | Google Scholar

Chiarabba, C., De Gori, P., and Boschi, E. (2009). Pore-pressure Migration along a normal-fault System Resolved by Time-Repeated Seismic Tomography. Geology 37, 67–70. doi:10.1130/G25220A.1

CrossRef Full Text | Google Scholar

Cinti, F. R., Moro, M., Pantosti, D., Cucci, L., and D’Addezio, G. (2002). New Constraints on the Seismic History of the Castrovillari Fault in the Pollino gap (Calabria, Southern Italy). J. Seismol. 6, 199–217. doi:10.1023/A:1015693127008

CrossRef Full Text | Google Scholar

D’Agostino, N., D’Anastasio, E., Gervasi, A., Guerra, I., Nedimović, M. R., Seeber, L., et al. (2011). Forearc Extension and Slow Rollback of the Calabrian Arc from GPS Measurements. Geophys. Res. Lett. 38 (17), a–n. doi:10.1029/2011GL048270

CrossRef Full Text | Google Scholar

De Landro, G., Amoroso, O., Russo, G., and Zollo, A. (2020). 4D Travel-Time Tomography as a Tool for Tracking Fluid-Driven Medium Changes in Offshore Oil-Gas Exploitation Areas. Energies 13 (22), 5878. doi:10.3390/en13225878

CrossRef Full Text | Google Scholar

De Landro, G., Serlenga, V., Russo, G., Amoroso, O., Festa, G., Bruno, P. P., et al. (2017). 3D Ultra-high Resolution Seismic Imaging of Shallow Solfatara Crater in Campi Flegrei (Italy): New Insights on Deep Hydrothermal Fluid Circulation Processes. Sci. Rep. 7, 3412. doi:10.1038/s41598-017-03604-0

PubMed Abstract | CrossRef Full Text | Google Scholar

De Matteis, R., Romeo, A., Pasquale, G., Iannaccone, G., and Zollo, A. (2010). 3D Tomographic Imaging of the Southern Apennines (Italy): A Statistical Approach to Estimate the Model Uncertainty and Resolution. Stud. Geophys. Geod 54 (3), 367–387. doi:10.1007/s11200-010-0022-x

CrossRef Full Text | Google Scholar

Dutta, N. C. (2002). Geopressure Prediction Using Seismic Data: Current Status and the Road Ahead. Geophysics 67 (6), 2012–2041. doi:10.1190/1.1527101

CrossRef Full Text | Google Scholar

Dvorkin, J., Prasad, M., Sakai, A., and Lavoie, D. (1999). Elasticity of marine Sediments: Rock Physics Modeling. Geophys. Res. Lett. 26 (12), 1781–1784. doi:10.1029/1999GL900332

CrossRef Full Text | Google Scholar

Ferranti, L., Milano, G., and Pierro, M. (2017). Insights on the Seismotectonics of the Western Part of Northern Calabria (Southern Italy) by Integrated Geological and Geophysical Data: Coexistence of Shallow Extensional and Deep Strike-Slip Kinematics. Tectonophysics 721, 372–386. doi:10.1016/j.tecto.2017.09.020

CrossRef Full Text | Google Scholar

Filice, F., and Seeber, L. (2019). The Culmination of an Oblique Time‐Transgressive Arc Continent Collision: The Pollino Massif between Calabria and the Southern Apennines, Italy. Tectonics 38, 3261–3280. doi:10.1029/2017TC004932

CrossRef Full Text | Google Scholar

Ghisetti, F., and Vezzani, L. (1982). Different Styles of Deformation in the Calabrian Arc (Southern Italy): Implications for a Seismotectonic Zoning. Tectonophysics 85 (3-4), 149–165. doi:10.1016/0040-1951(82)90101-9

CrossRef Full Text | Google Scholar

Got, J.-L., Fréchet, J., and Klein, F. W. (1994). Deep Fault Plane Geometry Inferred from Multiplet Relative Relocation beneath the South Flank of Kilauea. J. Geophys. Res. 99 (B8), 15375–15386. doi:10.1029/94JB00577

CrossRef Full Text | Google Scholar

Guerra, I., Harabaglia, P., Gervasi, A., and Rosa, A. B. (2005). The 1998-1999 Pollino (Southern Apennines, Italy) Seismic Crisis: Tomography of a Sequence. Ann. Geophys. 48 (6), 995–1007. doi:10.4401/ag-3249

CrossRef Full Text | Google Scholar

Improta, L., Bagh, S., De Gori, P., Valoroso, L., Pastori, M., Piccinini, D., et al. (2017). Reservoir Structure and Wastewater-Induced Seismicity at the Val d’Agri Oilfield (Italy) Shown by Three-DimensionalVpandVp/VsLocal Earthquake Tomography. J. Geophys. Res. Solid Earth 122, 9050–9082. doi:10.1002/2017jb014725

CrossRef Full Text | Google Scholar

Improta, L., Iannaccone, G., Capuano, P., Zollo, A., and Scandone, P. (2000). Inferences on Upper Crustal Structure of the Southern Apennines (Italy) from Seismic Refraction Investigations and Subsurface Data. Tectonophysics 317, 275–298. doi:10.1016/S0040-1951(99)00267-X

CrossRef Full Text | Google Scholar

Kissling, E., Ellsworth, W. L., Eberhart-Phillips, D., and Kradolfer, U. (1994). Initial Reference Models in Local Earthquake Tomography. J. Geophys. Res. 99 (B10), 19635–19646. doi:10.1029/93JB03138

CrossRef Full Text | Google Scholar

Latorre, D., Virieux, J., Monfret, T., Monteiller, V., Vanorio, T., Got, J.-L., et al. (2004). A New Seismic Tomography of Aigion Area (Gulf of Corinth, Greece) from the 1991 Data Set. Geophys. J. Int. 159 (3), 1013–1031. doi:10.1111/j.1365-246x.2004.02412.x

CrossRef Full Text | Google Scholar

Lomax, A., Virieux, J., Volant, P., and Berge-Thierry, C. (2000). “Probabilistic Earthquake Location in 3D and Layered Models,” in Advances in Seismic Event Location (Dordrecht: Springer), 101–134. doi:10.1007/978-94-015-9536-0_5

CrossRef Full Text | Google Scholar

Margheriti, L., Amato, A., Braun, T., Cecere, G., D’Ambrosio, C., De Gori, P., et al. (2013). Emergenza nell’area del Pollino: le attività della Rete Sismica Mobile. Rome: Rapporti tecnici INGV.

Google Scholar

Michele, M., Latorre, D., and Emolo, A. (2019). An Empirical Formula to Classify the Quality of Earthquake Locations. Bull. Seismol. Soc. Am. 109 (6), 2755–2761. doi:10.1785/0120190144

CrossRef Full Text | Google Scholar

Michetti, A. M., Ferreli, L., Serva, L., and Vittori, E. (1997). Geological Evidence for strong Historical Earthquakes in an “Aseismic” Region: the Pollino Case (Southern Italy). J. Geodyn. 24 (1/4), 61–86. doi:10.1016/s0264-3707(97)00018-5

CrossRef Full Text | Google Scholar

Napolitano, F., De Siena, L., Gervasi, A., Guerra, I., Scarpa, R., and La Rocca, M. (2020). Scattering and Absorption Imaging of a Highly Fractured Fluid-Filled Seismogenetic Volume in a Region of Slow Deformation. Geosci. Front. 11 (3), 989–998. doi:10.1016/j.gsf.2019.09.014

CrossRef Full Text | Google Scholar

Napolitano, F., Galluzzo, D., Gervasi, A., Scarpa, R., and La Rocca, M. (2020). Fault Imaging at Mt Pollino (Italy) from Relative Location of Microearthquakes. Geophys. J. Int. 224 (1), 637–648. doi:10.1093/gji/ggaa407

CrossRef Full Text | Google Scholar

Napolitano, F., Gervasi, A., La Rocca, M., Guerra, I., and Scarpa, R. (2018). Site Effects in the Pollino Region from the HVSR and Polarization of Seismic Noise and Earthquakes. Bull. Seismol. Soc. Am. 108 (1), 309–321. doi:10.1785/0120170197

CrossRef Full Text | Google Scholar

Orecchio, B., Presti, D., Totaro, C., Guerra, I., and Neri, G. (2011). Imaging the Velocity Structure of the Calabrian Arc Region (Southern Italy) through the Integration of Different Seismological Data. Bollettino di Geofisica Teorica Ed. Applicata 52 (4), 625–638. doi:10.4430/bgta0023

CrossRef Full Text | Google Scholar

Paige, C. C., and Saunders, M. A. (1982). LSQR: An Algorithm for Sparse Linear Equations and Sparse Least Squares. ACM Trans. Math. Softw. 8 (1), 43–71. doi:10.1145/355984.355989

CrossRef Full Text | Google Scholar

Passarelli, L., Hainzl, S., Cesca, S., Maccaferri, F., Mucciarelli, M., Roessler, D., et al. (2015). Aseismic Transient Driving the Swarm-like Seismic Sequence in the Pollino Range, Southern Italy. South. Italy. Geophys. J. Int. 201 (3), 1553–1567. doi:10.1093/gji/ggv111

CrossRef Full Text | Google Scholar

Passarelli, L., Roessler, D., Aladino, G., Maccaferri, F., and Moretti, M. (2012). “Pollino Seismic experiment (2012–2014),” in Deutsches GeoForschungsZentrum GFZ (Other/seismic network). doi:10.14470/9N904956

CrossRef Full Text | Google Scholar

Pastori, M., Margheriti, L., De Gori, P., Govoni, A., Lucente, F. P., Moretti, M., et al. (2021). The 2011-2014 Pollino Seismic Swarm: Complex Fault Systems Imaged by 1D Refined Location and Shear Wave Splitting Analysis at the Apennines-Calabrian Arc Boundary. Front. Earth Sci. 9, 90. doi:10.3389/feart.2021.618293

CrossRef Full Text | Google Scholar

Piana Agostinetti, N., and Amato, A. (2009). Moho Depth andVp/Vsratio in Peninsular Italy from Teleseismic Receiver Functions. J. Geophys. Res. 114 (B6). doi:10.1029/2008JB005899

CrossRef Full Text | Google Scholar

Podvin, P., and Lecomte, I. (1991). Finite Difference Computation of Traveltimes in Very Contrasted Velocity Models: a Massively Parallel Approach and its Associated Tools. Geophys. J. Int. 105 (1), 271–284. doi:10.1111/j.1365-246x.1991.tb03461.x

CrossRef Full Text | Google Scholar

Rovida, A., Locati, M., Camassi, R., Lolli, B., and Gasperini, P. (2016). CPTI15, the 2015 Version of the Parametric Catalogue of Italian Earthquakes. Istituto Nazionale di Geofisica e Vulcanologia 7, 112. doi:10.6092/INGV.IT-CPTI15

CrossRef Full Text | Google Scholar

Schiattarella, M. (1998). Quaternary Tectonics of the Pollino ridge, Calabria-Lucania Boundary, Southern Italy. Geol. Soc. Lond. Spec. Publications 135 (1), 341–354. doi:10.1144/GSL.SP.1998.135.01.22

CrossRef Full Text | Google Scholar

Sketsiou, P., De Siena, L., Gabrielli, S., and Napolitano, F. (2021). 3-D Attenuation Image of Fluid Storage and Tectonic Interactions across the Pollino Fault Network. Geophys. J. Int. 226 (1), 536–547. doi:10.1093/gji/ggab109

CrossRef Full Text | Google Scholar

Sketsiou, P., Napolitano, F., Zenonos, A., and De Siena, L. (2020). New Insights into Seismic Absorption Imaging. Phys. Earth Planet. Interiors 298, 106337. doi:10.1016/j.pepi.2019.106337

CrossRef Full Text | Google Scholar

Totaro, C., Koulakov, I., Orecchio, B., and Presti, D. (2014). Detailed Crustal Structure in the Area of the Southern Apennines-Calabrian Arc Border from Local Earthquake Tomography. J. Geodynamics 82, 87–97. doi:10.1016/j.jog.2014.07.004

CrossRef Full Text | Google Scholar

Totaro, C., Orecchio, B., Presti, D., Scolaro, S., and Neri, G. (2016). Seismogenic Stress Field Estimation in the Calabrian Arc Region (South Italy) from a Bayesian Approach. Geophys. Res. Lett. 43 (17), 8960–8969. doi:10.1002/2016GL070107

CrossRef Full Text | Google Scholar

Totaro, C., Seeber, L., Waldhauser, F., Steckler, M., Gervasi, A., Guerra, I., et al. (2015). An Intense Earthquake Swarm in the Southernmost Apennines: Fault Architecture from High-Resolution Hypocenters and Focal Mechanisms. Bull. Seismological Soc. America 105 (6), 3121–3128. doi:10.1785/0120150074

CrossRef Full Text | Google Scholar

Vanorio, T., Virieux, J., Capuano, P., and Russo, G. (2005). Three-dimensional Seismic Tomography fromPwave andSwave Microearthquake Travel Times and Rock Physics Characterization of the Campi Flegrei Caldera. J. Geophys. Res. 110 (B3), 1. doi:10.1029/2004JB003102

CrossRef Full Text | Google Scholar