Simulating three-dimensional seismograms in 2.5-dimensional structures by combining two-dimensional finite difference modelling and ray tracing

Finite difference (FD) simulation of elastic wave propagation is an important tool in geophysical research. As large-scale 3-D simulations are only feasible on supercomputers or clusters, and even then the simulations are limited to long periods compared to the model size, 2-D FD simulations are widespread. Whereas in generally 3-D heterogeneous structures it is not possible to infer the correct amplitude and waveform from 2-D simulations, in 2.5-D heterogeneous structures some inferences are possible. In particular, Vidale & Helmberger developed an approach that simulates 3-D waveforms using 2-D FD experiments only. However, their method requires a special FD source implementation technique that is based on a source definition which is not any longer used in nowadays FD codes. In this paper, we derive a conversion between 2-D and 3-D Green tensors that allows us to simulate 3-D displacement seismograms using 2-D FD simulations and the actual ray path determined in the geometrical optic limit. We give the conversion for a source of a certain seismic moment that is implemented by incrementing the components of the stress tensor.
Therefore, we present a hybrid modelling procedure involving 2-D FD and kinematic ray-tracing techniques. The applicability is demonstrated by numerical experiments of elastic wave propagation for models of different complexity.     

Regional tomographic inversion of the amplitude and phase of Rayleigh waves with 2-D sensitivity kernels

In this study, we test the adequacy of 2-D sensitivity kernels for fundamental-mode Rayleigh waves based on the single-scattering (Born) approximation to account for the effects of heterogeneous structure on the wavefield in a regional surface wave study. The calculated phase and amplitude data using the 2-D sensitivity kernels are compared to phase and amplitude data obtained from seismic waveforms synthesized by the pseudo-spectral method for plane Rayleigh waves propagating through heterogeneous structure. We find that the kernels can accurately predict the perturbation of the wavefield even when the size of anomaly is larger than one wavelength. The only exception is a systematic bias in the amplitude within the anomaly itself due to a site response.
An inversion method of surface wave tomography based on the sensitivity kernels is developed and applied to synthesized data obtained from a numerical simulation modelling Rayleigh wave propagation over checkerboard structure. By comparing recovered images to input structure, we illustrate that the method can almost completely recover anomalies within an array of stations when the size of the anomalies is larger than or close to one wavelength of the surface waves. Surface wave amplitude contains important information about Earth structure and should be inverted together with phase data in surface wave tomography.     

Slowness-backazimuth weighted migration: A new array approach to a high-resolution image

We have developed a new array method combining conventional migration with a slowness-backazimuth deviation weighting scheme. All seismic traces are shifted based on the theoretical traveltime of the scattered wave from specific gridpoints in a 3-D volume. Observed slowness and backazimuth are calculated for each raypath and compared with theoretical values in order to estimate slowness and backazimuth deviations. Subsequently, stacked energy calculated by a conventional migration method is weighted by the slowness and backazimuth deviations to suppress any arrival energy whose slowness and backazimuth are inconsistent with the expected theoretical values. This new method was applied to two P- wave data sets which comprise (1) underside reflections at the 410 and 660 km mantle discontinuities and (2) D" reflections as well as their corresponding synthetic data sets. The results show that the weighting scheme dramatically increases the resolution of the migrated images and enables us to obtain well-constrained, focused images, making upper-mantle discontinuities and D" reflections more distinct by reducing their surrounding energy.     

Comparison of global synthetic seismograms calculated using the spherical 2.5-D finite-difference method with observed long-period waveforms including data from the intra-Antarctic region

We have been developing an accurate and efficient numerical scheme, which uses the finite-difference method (FDM) in spherical coordinates, for the computation of global seismic wave propagation through laterally heterogeneous realistic Earth models. In the field of global seismology, traditional axisymmetric modeling has been used widely as an efficient approach since it can solve the 3-D elastodynamic equation in spherical coordinates on a 2-D cross-section of the Earth, assuming structures to be invariant with respect to the axis through the seismic source. However, it has the severe disadvantages that asymmetric structures about the axis cannot be incorporated and the source mechanisms with arbitrary shear dislocation have not been attempted for a long time. Our scheme is based on the framework of axisymmetric modeling but has been extended to treat asymmetric structures, arbitrary moment-tensor point sources, anelastic attenuation, and the Earth center which is a singularity of wave equations in spherical coordinates. All these types of schemes which solve 3-D wavefields on a 2-D model cross-section are classified as 2.5-D modeling, so we have named our scheme the spherical 2.5-D FDM. In this study, we compare synthetic seismograms calculated using our FDM scheme with three-component observed long-period seismograms including data from stations newly installed in Antarctica in conjunction with the International Polar Year (IPY) 2007–2008. Seismic data from inland Antarctica are expected to reveal images of the Earth's deep interior with enhanced resolution because of the high signal-to-noise ratio and wide extent of this region, in addition to the rarity of sampling paths along the rotation axis of the Earth. We calculate synthetic seismograms through the preliminary reference earth model (PREM) including attenuation using a moment-tensor point source for the November 9, 2009 Fiji earthquake. Our results show quite good agreement between synthetic and observed seismograms, which indicates the accuracy of observations in the Antarctica, as well as the feasibility of the spherical 2.5-D modeling scheme.     

Quadrangle-grid velocity–stress finite difference method for poroelastic wave equations

A quadrangle-grid velocity–stress finite difference method, based on a first-order hyperbolic system that is equivalent to Biot's equations, is developed for the simulation of wave propagation in 2-D heterogeneous porous media. In this method the velocity components of the solid material and of the pore fluid relative to that of the solid, and the stress components of three solid stresses and one fluid pressure are defined at different nodes for a staggered non-rectangular grid. The scheme uses non-orthogonal grids, allowing surface topography and curved interfaces to be easily modelled in the numerical simulation of seismic responses of poroelastic reservoirs. The free-surface conditions of complex geometry are achieved by using integral equilibrium equations on the surface, and the source implementations are simple. The algorithm is an extension of the quadrangle-grid finite difference method used for elastic wave equations.     

Teleseismic waveform modelling with a one-way wave equation

It is now widely accepted that elastic properties of the continental lithosphere and the underlying sublithospheric mantle are both anisotropic and laterally heterogeneous at a range of scales. To fully exploit modern three-component broad-band array data sets requires the use of comprehensive modelling tools. In this work, we investigate the use of a wide-angle, one-way wave equation to model variations in teleseismic 3-D waveforms due to 2-D elastic heterogeneity and anisotropy. The one-way operators are derived based on a high-frequency approximation of the square-root operator and include the effects of wave propagation as well as multiple scattering. Computational cost is reduced through a number of physically motivated approximations. We present synthetic results from simple 1-D (layer over a half-space) and 2-D (subduction zone) models that are compared with reference solutions. The algorithm is then used to model data from an array of broad-band seismograph stations deployed in northwestern Canada as part of the IRIS-PASSCAL/LITHOPROBE CANOE experiment. In this region radial-component receiver functions show a clear continental Moho and the presence of crustal material dipping into the mantle at the suture of two Palaeo-Proterozoic terranes. The geometry of the suture is better defined on the transverse component where subduction is associated with a ∼10 km thick layer exhibiting strong elastic anisotropy. The modelling reproduces the main features of the receiver functions, including the effects of anisotropy, heterogeneity and finite-frequency scattering.     

Scattering characteristics in heterogeneously fractured reservoirs from waveform estimation

Offset-dependent characteristics of seismic scattering are useful for characterizing fractured reservoirs. We use two models that have different background medium properties and different azimuthal AVO responses to study elastic wave propagation and scattering in gas-saturated, heterogeneously fractured reservoirs. Heterogeneous fracture density distributions are built through stochastic modelling. Synthetic seismograms are generated by 3-D finite difference modelling, and waveforms along crack-normal and strike directions are considered in this paper. The multiple signal classification (MUSIC) frequency estimator is used in waveform estimation to provide frequency-domain attributes related to seismic wave scattering by fracture heterogeneity. Our results indicate that the strength of the scattering field is a function of the background medium. The strength also increases with increasing fracture scatterer density and with decreasing correlation length of spatial variations of fracture density. The scattering field is weak at the top of the fractured reservoir. The first-order results are dominated by velocity anisotropy of the mean fracture density field. However, the base of the fractured reservoir corresponds to a strong scattering field on which fracture heterogeneity has a larger effect and is characterized by the loss of coherence.     

Surface wave tomography for azimuthal anisotropy in a strongly reduced parameter space

Large scale seismic anisotropy in the Earth's mantle is likely dynamically supported by the mantle's deformation; therefore, tomographic imaging of 3-D anisotropic mantle seismic velocity structure is an important tool to understand the dynamics of the mantle. While many previous studies have focused on special cases of symmetry of the elastic properties, it would be desirable for evaluation of dynamic models to allow more general axis orientation. In this study, we derive 3-D finite-frequency surface wave sensitivity kernels based on the Born approximation using a general expression for a hexagonal medium with an arbitrarily oriented symmetry axis. This results in kernels for two isotropic elastic coefficients, three coefficients that define the strength of anisotropy, and two angles that define the symmetry axis. The particular parametrization is chosen to allow for a physically meaningful method for reducing the number of parameters considered in an inversion, while allowing for straightforward integration with existing approaches for modelling body wave splitting intensity measurements. Example kernels calculated with this method reveal physical interpretations of how surface waveforms are affected by 3-D velocity perturbations, while also demonstrating the non-linearity of the problem as a function of symmetry axis orientation. The expressions are numerically validated using the spectral element method. While challenges remain in determining the best inversion scheme to appropriately handle the non-linearity, the approach derived here has great promise in allowing large scale models with resolution of both the strength and orientation of anisotropy.     

Solution of the inverse problem of seismology for laterally inhomogeneous media

Summary. The Backus-Gilbert method has been extended to the estimation of the seismic wave velocity distribution in 2-D or 3-D inhomogeneous media from a finite set of travel-time data. The method may be applied to the inversion of body wave as well as surface wave data. The problem of determining a local average of the unknown velocity corrections may be reduced to a choice of a suitable δ-ness criterion for the averaging kernel. For 2-D and 3-D inhomogeneous media the simplest criterion is to minimize a sum of 'spreads' over all the coordinates. The use of this criterion requires the solution (the averaged velocity corrections) to be represented as a sum of functions, each of which depends only on one coordinate. This is a basic restriction of the method. In practice it is possible to achieve good agreement between the solution and a real velocity distribution by a reasonable choice of the coordinate system.
Numerical tests demonstrate the efficiency of the method. Some examples of the application of the method to the inversion of real seismological data for body and surface waves are given.     

Seismic reflection traveltimes in two-dimensional statistically anisotropic random media

Velocity estimation remains one of the main problems when imaging the subsurface with seismic reflection data. Traveltime inversion enables us to obtain large-scale structures of the velocity field and the position of seismic reflectors. However, as the media currently under study are becoming more and more complex, we need to know the finer-scale structures. The problem is that below a certain range of velocity heterogeneities, deterministic methods become difficult to use, so we turn to a probabilistic approach. With this in view, we characterize the velocity field as a random field defined by its first and second statistical moments. Usually, a seismic random medium is defined as a homogeneous velocity background perturbed by a small random field that is assumed to be stationary. Thus, we make a link between such a random velocity medium (together with a simple reflector) and seismic reflection traveltimes. Assuming that the traveltimes are ergodic, we use 2-D seismic reflection geometry to study the decrease in the statistical traveltime fluctuations as a function of the offset (the source–receiver distance). Our formulae are based on the Rytov approximation and the parabolic approximation for acoustic waves. The validity and the limits are established for both of these approximations in statistically anisotropic random media. Finally, theoretical inversion procedures are developed for the horizontal correlation structure of the velocity heterogeneities for the simplest case of a horizontal reflector. Synthetic seismograms are then computed (on particular realizations of random media) by simulating scalar wave propagation via finite difference algorithms. There is good agreement between the theoretical and experimental results.     

Estimating rupture scenario likelihood based on dynamic rupture simulations: the example of the segmented Middle Durance fault, southeastern France

The Middle Durance fault system, southeastern France, is a slow active fault that produced moderate-size historical seismic events and shows evidence of at least one   M _w ≳ 6.5  event in the last 29 000 yr. Based on dynamic rupture simulation, we propose earthquake scenarios that are constrained by knowledge of both the tectonic stress field and of the 3-D geometry of the Durance fault system. We simulate dynamic rupture interaction among several fault segmentations of different strikes, dips and rakes, using a 3-D boundary integral equation method. 50 combinations of reasonable stress field orientations, stress field amplitudes and hypocentre locations are tested. The probability of different rupture evolutions is then computed. Each segment ruptures mainly as a single event (44 per cent of the 50 simulations test in this paper). However, the probability that an event triggers simultaneously along three segments is high (26 per cent), leading to a potential rupture length of 45 km. Finally, 2 per cent of the simulations occur along four adjacent segments, producing the greatest total rupture length of 55 km. The simulation results show that the southernmost segment is most easily ruptured (40 per cent), because of its favourable orientation with respect to the tectonic stress and of its favourable location for interaction with the other segments. South-bound unilateral propagation is slightly preferable (41 per cent), compared to north-bound unilateral and bilateral propagation modes. Although, these rupture scenarios cannot be directly translated into probabilities of occurrence, they do provide a better insight as to which rupture scenarios are more likely, an important element to better estimate near-field strong ground motion and seismic hazard.     

Three-dimensional seismic structure beneath the Australasian region from refracted wave observations

The earthquakes in the seismicity belt extending through Indonesia, New Guinea, Vanuatu and Fiji to the Tonga–Kermadec subduction zone recorded at the 65 portable broad-band stations deployed during the Skippy experiment from 1993–1996 provide good coverage of the lithosphere and mantle under the Australian continent, Coral Sea and Tasman Sea.
The variation in structure in the upper part of the mantle is characterized by deter-mining a suite of 1-D structures from stacked record sections utilizing clear P and S arrivals, prepared for all propagation paths lying within a 10° azimuth band. The azimuth of these bands is rotated by 20° steps with four parallel corridors for each azimuth. This gives 26 separate azimuthal corridors for which 15 independent 1-D seismic velocity structures have been derived, which show significant variation in P and S structure.
The set of 1-D structures is combined to produce a 3-D representation by projecting the velocity values along the ray path using a turning point approximation and stacking into 3-D cells (5° by 50 km in depth). Even though this procedure will tend to underestimate wave-speed perturbations, S -velocity deviations from the ak135 reference model exceed 6 per cent in the lithosphere.
In the uppermost mantle the results display complex features and very high S -wave speeds beneath the Precambrian shields with a significant low-velocity zone beneath. High velocities are also found towards the base of the transition zone, with high S -wave speeds beneath the continent and high P -wave speeds beneath the ocean. The wave-speed patterns agree well with independent surface wave studies and delay time tomography studies in the zones of common coverage.     

2.5-dimensional crosshole acoustic response: a variable density approach

This paper presents the development of a 2.5-D simulation technique for acoustic wave propagation in media with variable density and velocity. A comparative study of the 2-D and 2.5-D responses of a model reveals the spatially and temporally damped nature of the 2.5-D acoustic wave equations. The simulated results for constant and variable density models show that the density variation affects only the reflectivity of the layer. The computational cost for variable density models is 2.17 and 2.26 times that for constant density models for the 2.5-D and 2-D cases, respectively. Furthermore, the 2.5-D computational cost in the time domain is only about 10–15 per cent more than that for two dimensions, so this modest increase in computational cost can avoid the exorbitant 3-D computational cost.
Snapshots for a crosshole geometry were computed at various times in order to study the effect of heterogeneity on the amplitude and shape of the wave front. Extensive analysis of an oil-bearing reservoir with and without the inclusion of a gas zone was performed using a point source as well as multiple sources. In addition, the effects of the thickness of a low-velocity layer (oil-bearing) and of the location of the source have been studied. It is concluded from the numerical response that the waveguide action of the low-velocity layer depends on its thickness in terms of the dominant wavelength. Trapping of waves was not observed when the source was outside the low-velocity layer. Furthermore, the presence of heterogeneity in the low-velocity layer contributes considerably to the leakage of energy in the adjacent layers due to scattering/diffraction. It was found that, in the 2.5-D numerical simulation, the stability condition and the requirement of the number of grid points per wavelength to avoid grid dispersion are the same as for the 2-D case.     

Crustal structure of the Newfoundland rifted continental margin from constrained 3-D gravity inversion

The rifting history of the Atlantic continental margin of Newfoundland is very complex and so far has been investigated at the crustal scale primarily with the use of 2-D seismic surveys. While informative, the results generated from these surveys cannot easily be interpreted in a regional sense due to their sparse sampling of the margin. A 3-D gravity inversion of the free air data over the Newfoundland margin allows us to generate a 3-D density anomaly model that can be compared with the seismic results and used to gain insight into regions lacking seismic coverage. Results of the gravity inversion show good correspondence with Moho depths from seismic results. A shallowing of the Moho to 12 km depth is resolved on the shelf at the northern edge of the Grand Banks, in a region poorly sampled by other methods. Comparisons between sediment thickness and crustal thickness show deviations from local isostatic compensation in locations which correlate with faults and rifting trends. Such insights must act as constraints for future palaeoreconstructions of North Atlantic rifting.     

Seismic waves in rocks with fluids and fractures

Seismic wave propagation through the earth is often strongly affected by the presence of fractures. When these fractures are filled with fluids (oil, gas, water, CO₂, etc.), the type and state of the fluid (liquid or gas) can make a large difference in the response of the seismic waves. This paper summarizes recent work on methods of deconstructing the effects of fractures, and any fluids within these fractures, on seismic wave propagation as observed in reflection seismic data. One method explored here is Thomsen's weak anisotropy approximation for wave moveout (since fractures often induce elastic anisotropy due to non-uniform crack-orientation statistics). Another method makes use of some very convenient crack/fracture parameters introduced previously that permit a relatively simple deconstruction of the elastic and wave propagation behaviour in terms of a small number of crack-influence parameters (whenever this is appropriate, as is certainly the case for small crack densities). Then, the quantitative effects of fluids on these crack-influence parameters are shown to be directly related to Skempton's coefficient B of undrained poroelasticity (where B typically ranges from 0 to 1). In particular, the rigorous result obtained for the low crack density limit is that the crack-influence parameters are multiplied by a factor  (1 − B )  for undrained systems. It is also shown how fracture anisotropy affects Rayleigh wave speed, and how measured Rayleigh wave speeds can be used to infer shear wave speed of the fractured medium in some cases. Higher crack density results are also presented by incorporating recent simulation data on such cracked systems.     

A data-mapping approach to 3-D wavefield redatuming

In case of a complex overburden, the seismic data can be greatly improved by applying a full wavefield redatuming procedure. In practice, the application of the redatuming process to 3-D data acquired by conventional acquisition designs is non-trivial. Because of the large amount of data involved in the 3-D redatuming process and because of the sparseness of these data, it is impossible to apply conventional wave equation datuming directly.
We present a data mapping approach to redatuming (DMR), which follows the concept of Kirchhoff data mapping. A simplified background medium where no ray bending occurs is assumed for the medium below the datum in order to map an input data set referenced to the acquisition surface to an output data set referenced to the new datum level. The DMR method can be interpreted as a simplified version of the Kirchhoff summation redatuming (KSR) method, where one of the 2-D integrals over the acquisition coordinates can be solved analytically. Consequently, in this approach fewer traces are involved in the computation of one time sample (a 2-D integral is computed instead of a 4-D integral), which makes it particularly attractive for the application to 3-D data sets.
In this paper the theory underlying data mapping redatuming is discussed and the proposed approach is tested on fully sampled 2-D and 3-D synthetic data from models with both simple and complex velocity distributions in the subsurface.
The tests clearly show that the objective of producing results that are comparable to the conventional KSR has been achieved. The redatumed traces are dynamically and kinematically correct. Furthermore, these results confirm that the dependency of the new approach on the assumed medium below the datum level is, indeed, weak because the assumption of a velocity medium where no ray bending occurs is already sufficient to produce correct results.     

Windowed f-k spectra of a three-dimensional wavefield excited by a point source in a two-dimensional multilayered elastic medium

When full 3-D modelling is too costly or cumbersome, computations of 3-D elastic wave propagation in laterally heterogeneous, multilayered 2-D geological structures may enhance considerably our ability to predict strong ground motion for seismological and engineering purposes. Towards this goal, we extend the method based on the combination of the thin-layer finite-element and boundary-element methods (TLFE-BEM) and calculate windowend f - k spectra of the 3-D wavefield. The windowed f - k spectra are spatially localized spectra from which the local properties of the wavefield can be extracted. The TLFE-BEM is particularly suited for calculating the complete wavefield where surface waves are dominant in multilayered media. The computations are performed in the frequency domain, providing the f - k spectra directly. From the results for the 3-D wavefield excited by a point source in a 2-D multilayered, sloped structure, it can be said that the phase velocity of the fundamental-mode Rayleigh wave in a laterally heterogeneous multilayered medium, estimated from the windowed f - k spectra, varies with the location of the point source. For the model calculated in this article, the phase velocity varies between the value for the flat layered structure of the thick-layer side and that for the structure just under the centre of the window. The exact subsurface structure just under the centre of an array in a laterally heterogeneous medium cannot be obtained if we use the f - k spectral analysis assuming a flat layered structure.     

Seismic imaging of the laterally varying D" region beneath the Cocos Plate

We use an axisymmetric, spherical Earth finite difference algorithm to model SH -wave propagation through cross-sections of laterally varying lower mantle models beneath the Cocos Plate derived from recent data analyses. Synthetic seismograms with dominant periods as short as 4 s are computed for several models: (1) a D" reflector 264 km above the core–mantle boundary with laterally varying S -wave velocity increases of 0.9–2.6 per cent, based on localized structures from a 1-D double-array stacking method; (2) an undulating D" reflector with large topography and uniform velocity increase obtained using a 3-D migration method and (3) cross-sections through the 3-D mantle S -wave velocity tomography model TXBW. We apply double-array stacking to assess model predictions of data. Of the models explored, the S -wave tomography model TXBW displays the best overall agreement with data. The undulating reflector produces a double Scd arrival that may be useful in future studies for distinguishing between D" volumetric heterogeneity and D" discontinuity topography. Synthetics for the laterally varying models show waveform variability not observed in 1-D model predictions. It is challenging to predict 3-D structure based on localized 1-D models when lateral structural variations are on the order of a few wavelengths of the energy used, particularly for the grazing geometry of our data. Iterative approaches of computing synthetic seismograms and adjusting model characteristics by considering path integral effects are necessary to accurately model fine-scale D" structure.     

Palaeoseismicity of the Oquirrh fault, Utah from shallow seismic tomography

The magnitude and frequency of normal-fault palaeoearthquakes are usually determined by trenching studies that ascertain the size and number of colluvial wedges along the fault. Such information can be invaluable in predicting the seismic hazard and potential for a future earthquake in that region. Digging trenches across normal faults, however, is environmentally intrusive, expensive and limited in the penetration depth. To overcome these problems we propose the use of 3-D seismic tomography as a means to identify the shapes and sizes of colluvial wedges along normal faults. As an example,2-D and 3-D seismic surveys were conducted across the Oquirrh fault, Utah with the purpose of imaging the normal-fault structure to a depth of about 10  m. Results show that the 3-D tomogram clearly delineates the fault zone and a colluvial wedge, both of which correlate extremely well with the geological cross-section interpreted from an adjacent trench. The thickness of the colluvial wedge image is used in conjunction with a seismic section to compute an estimate of a 6.8 moment magnitude earthquake for the most recent event on this fault, which is in close agreement with the 7.0 estimate based on a nearby trenching study. These tomographic results demonstrate, for the first time, that seismic imaging methods can be used in some cases to estimate unambiguously the shapes of colluvial wedges and the sizes of prehistoric earthquakes. Thus, seismic tomography has the possibility of providing cheaper, deeper and wider, but less resolved, images of fault systems than the intrusive excavation of trenches across faults.     

Sequential integrated inversion of refraction and wide-angle reflection traveltimes and gravity data for two-dimensional velocity structures

A new algorithm is presented for the integrated 2-D inversion of seismic traveltime and gravity data. The algorithm adopts the 'maximum likelihood' regularization scheme. We construct a 'probability density function' which includes three kinds of information: information derived from gravity measurements; information derived from the seismic traveltime inversion procedure applied to the model; and information on the physical correlation among the density and the velocity parameters. We assume a linear relation between density and velocity, which can be node-dependent; that is, we can choose different relationships for different parts of the velocity–density grid. In addition, our procedure allows us to consider a covariance matrix related to the error propagation in linking density to velocity. We use seismic data to estimate starting velocity values and the position of boundary nodes. Subsequently, the sequential integrated inversion (SII) optimizes the layer velocities and densities for our models. The procedure is applicable, as an additional step, to any type of seismic tomographic inversion.
We illustrate the method by comparing the velocity models recovered from a standard seismic traveltime inversion with those retrieved using our algorithm. The inversion of synthetic data calculated for a 2-D isotropic, laterally inhomogeneous model shows the stability and accuracy of this procedure, demonstrates the improvements to the recovery of true velocity anomalies, and proves that this technique can efficiently overcome some of the limitations of both gravity and seismic traveltime inversions, when they are used independently.
An interpretation of field data from the 1994 Vesuvius test experiment is also presented. At depths down to 4.5 km, the model retrieved after a SII shows a more detailed structure than the model obtained from an interpretation of seismic traveltime only, and yields additional information for a further study of the area.