About the shallow resistivity structure of Vesuvius volcano

Magnetotelluric (MT) soundings performed in the past in the volcanic area of Mt. Vesuvius by two independent research groups showed in the same places MT apparent resistivity curves with very similar shape, but statically shifted by one order of magnitude, at least. To try to resolve this ambiguity new controlled source audio-magnetotelluric (CSAMT) measurements have been carried out in the same MT sites. The interpretation of the CSAMT dataset, combined with that of two shallow dipole-dipole geoelectrical resistivity tomographies previously carried out in the area have allowed a reliable electrical structure to be recovered down to a few km of depth, which will next be used for a best constrained re-interpretation of the deep MT soundings. Mailing address: Dr. Antonio Troiano, Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Vesuviano, Napoli, Italy; e-mail: troiano@ov.ingv.it Vol51,1,2008_DelNegro 16-02-2009 21:28 Pagina 181


Introduction
Mt. Vesuvius is considered a highly hazardous volcano and is among the most surveyed active volcanoes in the world for the great concern due to the high level of urbanization existing all around its slopes, mainly those facing closest to the south-eastern suburbs of Naples.
As is well known, geophysics is considered a powerful tool both for modelling the structure and controlling the dynamics of active volcanoes.In particular, the application of electromagnetic (EM) geophysical methods to vol-canology is a topic of great interest, given the strong dependence of the electrical resistivity on the shallow and deep physical characteristics of a volcanic apparatus.
The EM method used to explore deep crustal structures is magnetotellurics (MT), whose investigation depth ranges from a few hundred meters to tens of kilometres.MT was applied to Vesuvius by Di Maio et al. (1998) and subsequently by Manzella et al. (2004).The MT apparent resistivity curves drawn by the two research teams, while showing essentially the same pattern and an intermediate minimum localized at about the same frequency, appear vertically shifted by one up to two orders of magnitude even in some common stations, without significant differences in phases.
Indeed, the problem of finding the correct level of the MT curves has been, and still is a crucial step for a reliable resistivity modelling of the explored structures.As is well known, a change in level of a curve leads to a dramatic change in the depth of allocation of the resistivity sequences, just as documented by the inter-preted sections provided by the two research groups across a nearly coincident N-S profile through the Vesuvius edifice.In fact, Di Maio et al. (1998) interpreted the intermediate apparent resistivity minimum they found in the MT curves, as due to a deep conductor around 12 km of mean depth b.s.l.(below sea level), rising up to about 8 km of mean depth bsl beneath the Vesuvius cone.Manzella et al. (2004), instead, interpreted the similarly shaped, but downshifted minimum of their MT curves, as due to a conductive zone at around 1 km of mean depth bgl (below ground level).As a consequence, the volcanological attribute of this intermediate conductive zone turned out to be quite different, depending on the depth at which it was located, i.e. possible evidence of an ascending intracrustal molten layer (Di Maio et al., 1998), or a shallow hydrothermal reservoir (Manzella et al., 2004).
Di Maio et al. (1998) maintain that the level of their MT curves is confirmed by the comparison with a nearby N-S Dipole-Dipole (DD) Electrical Resistivity Tomography (ERT) reported in their same paper.Manzella et al. (2004), after making recourse to TDEM measurements to control the static level of their MT curves, maintain that their joint TDEM-MT interpreted resistivity model is in good agreement with the same N-S DD ERT reported in Di Maio et al. (1998).It is thus evident that the two groups made a quite different use of the N-S DD ERT image reported in Di Maio et al. (1998), which actually is nothing but a pseudosection.
To try to solve this problem, a shallow 2D resistivity imaging of the Vesuvius area, where the MT stations were located is provided in this paper, by combining the interpretation of the already existing DD ERT data and a new Controlled Source Audio-Magnetotelluric (CSAMT) dataset.In particular, the 2D inversion of the already mentioned N-S DD ERT profile (Di Maio et al., 1998) and an E-W DD ERT profile, also previously performed in the area (Patella and Mauriello, 1999) is illustrated and discussed, followed by that of the new CSAMT profile performed in the Vesuvius central western area.
The choice of the CSAMT method is motivated by the fact that it uses an artificial EM field to estimate the impedances in the 10-10 5 Hz frequency band, where the natural MT signals are typically weak, thus allowing shallow targets to be most reliably detected.
The combined CSAMT and DD ERT interpretation down to a few km of depth is here proposed for the level of the MT curves to be constrained in an optimal way, and hence the deeper electrical structures to be more reliably imaged.The MT re-interpretation will be the topic of a future paper.

Previous geoelectrical surveys revisited
Besides the shallow DD ERT reported in Di Maio et al. (1998) and Patella and Mauriello (1999), a near-surface DD ERT was also performed in the Vesuvius area by Supper and Seiberl (2000).We consider here only the nearsurface curved ERT profile, about 7.9 km long (black line in fig.1), located very close to the CSAMT survey.Following Supper and Seiberl's (2000) interpretation down to the investigation   depth of 150 m bgl, the top northwestern branch of the profile is generally dominated by a bottom layer with high resistivity ranging between 750 and 2000 Ωm, overlain by an outcropping layer showing a still higher resistivity in the range 6000-20000 Ωm and a thickness ranging between a few m and about 60 m.Along the southeastern branch of the profile, the general appear-ance changes completely.Only in some segments is there a high resistivity top layer as in the former branch.Generally, the top portion of the section is represented by a lower resistivity layer (200-2000 Ωm) extending down to a depth of about 60 m, followed by an extremely high resistivity layer (12000-60000 Ωm).From this nearsurface resistivity overview, it clearly emerges that, in this zone of the volcano and within the top 150 m of depth, the relatively lower resistivity zones have a limited extension and an average resistivity of about 1000 Ωm, whereas the highest resistivity zones have average resistivity of even some ten thousands Ωm.
To infer the resistivity behavior below 150 m of depth bgl down to a few km of depth bgl, at most, we now show the 2D inversion models which have been elaborated across the shallow N-S (Di Maio et al., 1998) and E-W (Mauriello and Patella, 1999) DD ERT, located as in fig. 1.Both ERT profiles were surveyed using a dipole spread of 500 m, and the spacing between the centre of the emitting and receiving dipoles was enlarged at the constant step of 500 m, starting from 1 km up to 4.5 km, at most.The two DD ERT have been inverted using the method and computer code proposed by Oldenburg andLi (1994, 1999).We have chosen a uniform half-space as background reference model and made subsequent inversions with different values for its resistivity, ranging between 10 Ωm and 10000 Ωm.This procedure allows the depth of investigation (termed DOI index) to be evaluated until the model is controlled by the experimental data.In practice, the DOI index, which ranges between 0 and 1, quantifies the change between the models obtained by the different background resistivity.DOI values not exceeding 0.1, for reference resistivities differing by two orders of magnitude, indicate that the area is well constrained by the data.We have found that reference resistivities below 100 Ωm do not let the inversion converge.The resistivity models finally obtained along the two DD ERT profiles are shown in fig.2, while fig. 3 shows the comparison between experimental and modelled responses.
Assigning a 25% average error on the data, a rms values of 3.9 and 3.6 has been recovered respectively for the N-S and WE profiles.
The two sections display a nearly continu- ous conductive horizon with a resistivity ranging from 30 up to 300 Ωm, underlying a more resistive and discontinuous cover, undersampled by the dipole length.The thickness of the conductor, probably related to the aquifer, increases coming towards the NW sector, reaching depths of about 1500 m b.s.l.A bedrock with a resistivity not less than about 2000 Ωm and up 20000 Ωm is recovered in the deepest portion of the model.In the N-S section, beneath Vesuvius, a more conductive zone interrupts this horizon, likely ascribable to local hydrothermal systems (Caliro et al., 2005).

The CSAMT survey
CSAMT data were acquired in single site mode at eight stations distributed along a nearly straight NW-SE transect, about 5 km long and passing about 2 km south of the Vesuvius crater, as shown in fig. 1, in the same central area previously explored by Di Maio et al.
(1998) and Manzella et al. (2004).A Stratagem (Geometrics/EMI) system has been used.The controlled source emits in a band from 1 kHz to 100 kHz and for the lower frequencies the natural source is utilised.The antennas were placed near the receiver, at a distance such that both the s/n were high and no near source effects were present.The DD ERT dataset was used to validate the level of the apparent resistivity curves.
The rotationally invariant skew parameter (Swift, 1967) was derived from the MT impedance tensor in order to measure the deviation from 1D or 2D structures.As shown in fig.4, at almost all the periods the skew values are less than 0.2, with the exception of sites 3 and 4, thus globally indicating a 1D/2D nature of the data.
The occurrence of near-surface anomalies also poses a problem for the determination of the true regional strike direction.We used Groom and Bailey's approach (Groom and Bailey, 1989) to recover the 2D regional strike and the telluric distortion parameters.In particular, we have analyzed the three decomposition parameters, namely shear, twist and regional strike angle, by the rose diagrams in fig. 5. Shear and twist nearly vanish at almost all sites (no local distortion), with the exception of the sites with high skew values.The regional strike angle is approximately 60°N or, equivalently, 30°N, due to the 90°ambiguity in the strike direction, at all the sites.The ambiguity about the strike direction has been removed using the ERT 2D inversion (fig.2), which across the E-W profile shows a more homogeneous behavior than across the N-S profile.The CSAMT apparent resistivity and phase pseudo-sections relative to the TM and TE modes are shown in fig.8.
The CSAMT soundings were inverted to find the best model mb justifying the whole dataset along the profile.We employed the regularized algorithm by Rodi and Mackie (2001), which is based on the conjugate gradients method for the iterative minimization of the following objective function The idea behind the minimization of such a functional, due to Tikhonov and Arsenin (1977), who named the general procedure regularization, was introduced to overcome mathematical difficulties in the theory of ill-posed problem.Regularization has enormous practical benefits in the interpretation of experimental data.Specifically, Φm has a smoothing effect on the inversion results, in the sense that, minimizing a functional of this kind, one finds a minimum structure model (low Φm) with a reasonably low data misfit Φd.Other models will be able to satisfy the observation, due to the non-uniqueness of the problem, but this approach avoids unnec- essary complex solutions.The best model mb depends strongly on the structural information of the initial model m 0 , which may lead to biased inversion results if it is not consistent with the data.A homogeneous half space with a properly chosen electrical resistivity ρ 0 can be used as initial model m 0 (Tauber et al., 2003).
The inversion results depend significantly also on the trade-off parameter τ.
In order to perform a correct evaluation of the best model mb, we applied a minimum error floor of 10 per cent on log(ρ a ) and 5 per cent on the phases using different starting models for the inversion, and found that the least rms error corresponds with ρ 0 =10000 Ωm.Moreover this high ρ 0 value agrees with the very resistive bedrock indicated by DD ERT section shown in fig. 2 and by the near-surface DD ERT performed in the Vesuvius area by Supper and Seiberl (2000).Furthermore, we analyzed the behavior of the rms error versus the regularization term Φm, at different τ values, as plotted in fig.6.The best compromise between a low data misfit and a smooth model is obtained for τ=1.Therefore, we have at last accepted as best model mb the one with ρ 0 =10000 Ωm and τ=1, shown in fig. 7. The final rms achieved is of about 3.3, while the response curves of the best fit model are compared with the data in fig.8.
As previously said, due to the non-uniqueness of the CSAMT data inversion problem, other 2D models may result consistent with the given dataset.To explore this aspect, we consider now the behavior of the so-called sensitivity matrix S in the accepted best model mb, defined as (Scwalenberg et al., 2002) (3.4) which describes the influence of an infinitesimal resistivity variation in any block of the 2D best model, on the data.Although the sensitivity matrix calculation represents the most timeconsuming phase of the inversion, its analysis is indispensable to understand which of the model parameters are less resolved by the data, in the sense that different values of these model parameters may generate a final model fitting the data as adequately as the chosen best model.In practice, each element of the sensitivity matrix S is normalized by data error, grid size and the maximum value (Scwalenberg et al., 2002).The results are shown in fig.9.
To fix a threshold for the sensitivity, an analysis of the resolved part of the model was performed.In fig. 10 we have plotted the rms obtained imposing different values of resistivity below a chosen isoline of the sensitivity.Isolines above 10 −5 show some variation of the relative rms, when the resistivity of part of the model is changed, while below 10 −5 the rms remains quite constant, showing that the relative part of the model is not well resolved.Therefore, zones with sensitivity below this threshold are likely to be excluded from the interpretation.
Sensitivity generally decreases rapidly with depth, as shown in fig.9, where, excluding some zones with sensitivity values below the said threshold, the structures down to 2 km of depth appear well imaged.
Accounting for this last conclusion, we have redrawn in fig.11 the final best model derived  from the inversion of the CSAMT, where the lowest sensitivity zones have been excluded from model parameterization.
The uppermost part of the section appears characterized by a heterogeneous layer, consisting of a conductive region, with resistivities not exceeding 200 Ωm, locally as low as a few tens Ωm, associated with hydrothermal or a mix of groundwater and hydrothermal systems (Caliro et al., 2005), broken by zones of resistivities up to some thousands Ωm, related to cold, compact volcanic rocks.A nearly continuous bedrock with a high resistivity extends at depths below 1 km b.s.l.The general resistivity pattern across the CSAMT profile appears reasonably similar to that emerging from the shallow DD ERT profiles roughly in the same area (see fig. 2).

Conclusions
We have approached the inversion of DD ERT and CSAMT datasets collected in the Mt.Vesuvius area to provide a basis to resolve the ambiguity related to static shift effects in some MT apparent resistivity curves, previously obtained in the same area by two different research teams (Di Maio et al., 1998;Manzella et al., 2004).MT apparent resistivity curves vertically shifted by one up to two orders of magnitude were indeed displayed even in some common stations, while showing the same pattern, an intermediate minimum localized at about the same frequency and no significant differences in phases.
The choice of the ERT and CSAMT methods was motivated as follows.ERT is not afflicted by static shift effects, and, in conjunction with CSAMT, which in the 10-10 5 Hz band is characterized by a high signal-to-noise ratio, provides the optimum way to fix unambiguously the level of the MT apparent resistivity curves.
The shallow ERT and CSAMT interpreted sections have shown in the whole volcanic central area the existence of very similar behavior of the resistivity, which spans over five orders of magnitude, from less than 10 Ωm up to 10 5 Ωm.A first heterogeneous layer appeared in both sections, confined within the first km of depth bsl and consisting of a conductive zone, with resistivities ranging between a few tens Ωm and some hundreds Ωm.It was associated with a sequence of hydrothermal or mixed groundwater-hydrothermal systems, separated by more resistive parts, likely related to cold, compact volcanic rocks.This layer overlies a very highly resistive bedrock, showing resistiv- show some variation of the relative rms, while below 10 -5 the rms remains quite constant, showing that the relative part of the model was not well resolved.ities ranging from a few thousands Ωm up to some ten thousands Ωm.
The high resolution level of this shallow part of the volcanic structure is definitely assumed to provide the required strong constraints in view of a re-interpretation of the MT soundings available in this area, for the exploration of the deeper lithosphere structure, which will be the topic of a future paper.

Fig. 1 .
Fig. 1.Sketch map of the Vesuvius volcanic area.UTM coordinates are shown on both axes in figure.Red ellipses: CSAMT stations.Green lines: shallow DD ERT profiles.Black curved line: near-surface ERT profile.

Fig. 3 .
Fig. 3. Pseudosections of the N-S and E-W DD ERT profiles and responses of the best fit N-S and E-W models shown in fig. 2. Data are fitted with an rms of 3.9 and 3.6 respectively.

Fig. 2 .
Fig. 2. 2D resistivity models across the N-S and E-W shallow DD ERT profiles (green straight lines in fig.1).Black arrows indicate the profiles intersection.

Fig. 4 .Fig. 5 .
Fig. 4. Swift's skew as a function of period and CSAMT site number.Skew values below 0.2 indicate a 1D/2D nature of the dataset.

Fig. 6 .
Fig. 6.The rms error as a function of model roughness relative to different values of the trade-off parameter τ.

Fig. 7 .
Fig. 7. 2D resistivity ending models obtained by the joint inversion of the CSAMT TM and TE modes, using a homogeneous half-space with resistivity ρ0=10000 Ωm and τ=1 for the regularization parameter.

Fig. 8 .
Fig. 8. Pseudosection for the CSAMT profile and response of the best fit model shown in fig. 7. Data are fitted with an rms of 3.3.

Fig. 9 .
Fig. 9. Plot of the normalized sensitivity matrix related to the model of fig. 6. Blocks in the sections where the common log of sensitivity takes values below −5 are not interpretable.

Fig. 10 .
Fig. 10.Analysis of the resolved part of the model of fig. 6.The rms are obtained imposing constant values of resistivity below a chosen isoline of the sensitivity in the model of fig. 6. Isolines above 10 -5show some variation of the relative rms, while below 10 -5 the rms remains quite constant, showing that the relative part of the model was not well resolved.

Fig. 11 .
Fig. 11.2D joint inversion of CSAMT TM and TE modes across the profile in fig. 1, relative to 10000 Ωm for the resistivity ρ0 of the starting model, τ=1 for the trade-off parameter.A horizon of 2 km bgl is assumed.