Capturing the spatial variability in soil is crucial for ground response analyses in the context of seismic hazard mitigation. The lateral variability in thickness and properties of the different soil layers is one of the main factors that determines the variability of the ground motion spectrum from one location to another. The absence of such lateral variability information in the subsoil in between the locations of Cone Penetration Tests (CPTs) may be compensated by the use of more densely sampled seismic data. In this research we aim to derive a shear-wave velocity field through seismic full-waveform inversion that yields a model resolution approaching that of high-resolution seismic CPT surveys. Following this, a datadriven correlation between geophysical and geotechnical information is attempted through the application of new machine-learning-based approaches. @en

The Marchenko method is capable of estimating Green’s functions between the surface of the Earth and arbitrary locations in the subsurface. These Green’s functions are used to redatum wavefields to a deeper level in the subsurface. The Marchenko method enables the isolation of the response of a specific layer or package of layers, free from the influence of the overburden and underburden. In this study, we apply the Marchenko-based isolation technique to land S-wave seismic data acquired in the Groningen province, the Netherlands. We apply the technique for combined elimination of the overburden and underburden. Our results indicate that this approach enhances the resolution of reflection data. These enhanced reflections can be utilised for imaging and monitoring applications.@en

Seismic interferometry (SI) retrieves the Green function between two receiver locations using their recordings from a boundary of sources. When using sources and receivers only at the surface, the virtual-source gathers retrieved by SI contain pseudo-physical reflections as well as ghost (non-physical) reflections. These ghost reflections are the results of the cross-correlation or auto-correlation (AC) of primary reflections from two different depth levels, and they contain information about the seismic properties of specific layers in the subsurface. We investigated the application of ghost reflections for layer-specific characterization of the shallow subsurface using SI by AC. First, we showed the technique's potential using synthetic data for a subsurface model with a lateral change in velocity, a gradient in depth for velocity, a thickness change and a velocity change of the target layer. Then, we applied the technique to shallow subsurface field data. We also focused on improving the retrieval of ghost reflections by removing the free-surface multiples and muting undesired events in active-source gathers before applying SI. Our results demonstrate that the ghost reflections can be used advantageously to characterize the layer that causes them to appear in the results of SI. Consequently, they can also provide valuable information for imaging and monitoring shallow subsurface structures.

@en

The absence of information on lateral variability in the soil is detrimental to estimating accurately the local site response in the event of an earthquake. To address this problem, the use of densely sampled seismic data together with sparsely distributed but detailed vertical soil profiles obtained from cone penetration tests (CPTs) is advantageous. This study explores the adaptation of suitable machine learning (ML) approaches to derive reliable, site- and depth-specific correlations between seismic shear-wave velocity (Vs) and cone-tip resistance (qc). Such correlation could be successfully established by combining information from seismic CPT surveys with available borehole information for the Groningen region in the Netherlands. It is found that, even over substantial distances, ML-based techniques offer site- and depth-specific correlations between Vs and qc.@en

To investigate the physical processes behind induced seismicities due to, for example, production of hydrocarbons from a reservoir, most of the earlier studies performed geomechanical simulations on a simple reservoir geometry. The effect of fluid depletion is, in general, simulated for such a simple geometry. Neglecting the contribution of realistic 3-D reservoir geometries can lead to a wrong estimation of the incremental stress field. A reliable estimate of the induced stress field is key to producing meaningful simulation results. We perform geomechanical simulations on a simple fault model as well as a more realistic model based on the known geological structures at the earthquake source-region in Zeerijp region, the Netherlands. Our results demonstrate that the angle of the fault intersection affects the incremental stress field, including the effective normal stress, the shear stress, and hence, the Coulomb stress and the SCU value. Our results also show a shift in the rupture pattern and the location of the maximum slip on the fault plane. We conclude that, to properly evaluate the effects of production activities and to simulate precisely the in-situ stress field and the induced seismicity, the incorporation of a realistic reservoir structure in modelling is essential.@en

Intersecting faults are often ignored in the geomechanical simulation of induced seismicity. To investigate the effects of fault intersection and the resulting reservoir geometry on induced seismicity, caused, for instance, by gas extraction, we have developed 3D geomechanical models considering two intersecting normal faults and the surrounding horst structure. We simulate the stress field and the dynamic fault reactivation in a uniformly depleted reservoir. We observe that a smaller intersection angle increases the incremental Coulomb stress at the lower reservoir juxtaposition, thus changing the temporal rupture pattern of the seismic event. In our dynamic simulation, the rupture propagates from the main fault to the secondary fault. We conclude that the fault intersection has important effects on the induced seismicity and should be taken into account when evaluating the seismicity risk in a specific region.

@en

Bulk-density ( ρ) of soil is an important indicator of soil compaction and type. A knowledge of the spatial variability of in situ soil density is important in geotechnical engineering, hydrology and agriculture. Surface geophysical methods have so far shown limited success in providing an accurate and high-resolution image of 3-D soil-density distribution. In this pursuit, 3-D seismic full-waveform inversion (FWI) is promising, provided the robustness and accuracy of density inversion via this approach can be established in the near-surface scale. Ho wever , simultaneous reconstruction of ρand seismic wave velocities through multiparameter FWI remains a challenging task. Near-surface seismic data are commonly dominated by dispersive surface waves whose velocities are controlled by the value and distribution of shear-wave velocity ( V S ). One major difficulty in estimating reliab ly ρfrom near -surface seismic data is due to the relati vel y low sensiti vity of the seismic w av efield to ρcompared to seismic v elocities. Additionally, the accuracy of the estimated ρdecreases due to error in V S -an issue known as parameter coupling. Parameter coupling makes it difficult to estimate accurately ρwithin the framework of conventional gradient-based FWI. More sophisticated optimization approaches (e.g. truncated Newton) can reduce the effect of parameter coupling, but these approaches are commonl y not af fordab le in near -surface applications due to heavy computational burden. In this research, w e ha v e inv estigated how choosing correctly the force direction of the seismic source can contribute to a higher accuracy of ρestimates through 3-D FWI. Using scattered wavefields, the Hessian, and inversion tests, an in-depth and systematic investigation of data sets corresponding to different force directions has been carried out. A comparison of the scattered wavefields due to a point-localized ρperturbation for different force directions shows the robustness of the horizontal-force data set to noise compared to the vertical-force data set. Fur ther more, for a point-scatterer model, an analysis of the gradients of the misfit function using the Hessian shows that utilizing a horizontal-force source enables one to reconstruct the high-resolution gradient with relati vel y small parameter coupling. Finally, inversion tests for two different subsoil models demonstrate that 3-D FWI on a horizontal-force-source seismic data set is capable of providing a more accurate 3-D ρdistribution in soil compared to a vertical-force-source data set. Our results show that the use of a horizontal-force source might allow avoiding computationally demanding, costly optimization approaches in 3-D FWI.

@en

Using post-critical reflection data, it is possible to obtain useful information that allows more reliable geological characterization of the subsurface. However, the strong distortion caused by the phase shift in post-critical wavelets makes the use of post-critical reflections rather challenging. For this reason, an approach which is capable of estimating the phase shift of each wavelet of a reflection event in a data-driven manner is desirable. In this vein, in case the frequency spectrum of a wavelet can be correctly estimated, it is possible to estimate the instantaneous phase shift. In this work, we propose an approach which can perform such estimation based on spectral recomposition of seismic data. We design an inversion approach in order to reconstruct the seismic spectrum of the wavelets of a reflection event, which subsequently allows us to estimate the instantaneous phase of each wavelet of the near-surface reflection events without performing prior velocity analysis and/or critical-angle estimation. After finding the instantaneous phase for each wavelet of a reflection event, we show next how one can find the respective phase shifts that can then be corrected.

@en

Using the critical angle information of a reflection event, it is possible to calculate several essential physical parameters that are key to reliable geological characterization of the subsurface. However, estimation of the critical angle usually requires several steps of seismic processing. For this reason, an approach which is capable of estimating the critical angle directly from the data is of interest. Once the critical angle is estimated, it is possible to estimate further the Poisson's ratio and the seismic velocities. In this work, we propose an approach which can perform this estimation, based on spectral recomposition of seismic data. We design an inversion scheme in order to reconstruct the seismic spectrum of wavelets of a reflection event, which subsequently allows us to estimate the critical angle of near-surface reflection events without performing prior velocity analysis. After finding the critical angle, we show next how to estimate the Poisson's ratio and the compressional- and shear-wave velocities of the medium above the reflector. The approach leads to quite accurate values for Poisson's ratio even for noisy data, in case the number of layers is not large.

@en

The Azambuja fault is a NNE trending structure located 50 km north of Lisbon, the capital and most populous city of Portugal. The fault has been considered as a possible source for the historical, large earthquakes. Understanding this fault is a priority in seismic hazard evaluation of this region. The fault has a clear morphological signature. Miocene and Pliocene sediments are tilted eastward and cut by steeply dipping mesoscale fault segments, presenting reverse and normal offsets with a net downthrow to the east. Neotectonic studies indicate a Quaternary slip on the fault of 0.05–0.06 mm/year. However, no direct evidence of the Azambuja fault affecting the Pleistocene or Holocene sediments was found so far. Here, we present the findings from high-resolution seismic reflection studies using both P- and S-waves over the Holocene deposits. The detection of small-throw faulting in ductile sediments is a challenging task. We show that multiple signatures, like perturbations in the reflection hyperbolae visible in shot and CMP gathers, interruptions of reflectors in stacked sections, lateral seismic velocity variations obtained by horizon velocity analysis, all at coincident locations, strongly suggest that the activity of the Azambuja fault has affected the Holocene sediments in the study area. The lateral velocity variations are corroborated by wavepath eikonal traveltime tomography and velocity analysis supported by seismic modeling. By means of 2D viscoelastic modeling, we explain the absence of fault-related diffractions and negligible back-scattered energy from the fault. Using data from nearby boreholes, we find that the 15 ka old alluvium cover has indeed been disturbed by the presence of shallow fault strands. Considering the estimated vertical throws and the empirical relationships between fault length, co-seismic rupture and magnitude, a slip rate of 0.07 mm/y, slightly larger than previously thought, is expected for this fault.

@en

Geomechanical modelling is generally used to simulate the nucleation of induce d earthquakes in, for instance the Groningen gas field. We apply quasi static simulation to investigate the stress changes from gas production. When a fault reaches a critical state, dynamic simulation provides information on the dynamic rupture during ea rthqu ake nucleation and the resulting wavefield . With the use of geomechanical modelling, it is possible to investigate the effects of the model parameters, e.g., depletion pattern and friction parameters. I n the modelling, the dynamic rupture at a finite fault is simulated both in space and time. The generated seismic wavefield from such a finite source is considered to be more realistic than the resulting wavefield from a point source. T he latter is often assumed in previous studies on the inversion of in duced earthquake data in the Groningen area. To link the wavefield generated by a geomechanically simulated finite source to the field seismic data for an earlier earthquake, we apply the same full moment tensor inversion to the waveform of a finite and of a point source . The inverted moment tensor from the field seismic observation provides a constraint to our geomechanical simulation. This allows us to perform a more realistic simulation of an induced earthquake.@en

Dynamic geomechanical modeling can generate the seismic wavefield caused by a fault rupture. In dynamic fault-rupture modeling, the source is considered to be finite, with a limited extent both in space and in time. This contrasts with the definition of a point source, which is generally assumed to explain the seismic wavefield caused by an earthquake. Most earlier seismic inversion studies, including those of the induced earthquakes caused by depletion of the Groningen gas field, were performed assuming a point source. Still, finding a point-source reference from the seismic wavefield, even when generated by finite faulting, is important in order to calibrate the geomechanical simulation with field-seismic observations. To this end, we have developed a workflow that links geomechanical forward modeling to seismic moment-tensor inversion. We have tested this workflow for the dynamic rupture considering a realistic 3D layered earth model. At first, we simulate the triggering of dynamic fault slip at the center of a fault plane. Next, we invert the seismograms recorded by receivers located on or near the surface to obtain the full moment-tensor point-source representation and the location of the earthquake. The results of inversion show similar waveforms for both the point source and the finite source. The location of the inverted point source is within 400 m from the center of the slip patch. The double-couple components of the inverted moment tensor also match with the strike and the dip of the fault plane.@en

Traditional least-squares full-waveform inversion (FWI) suffers from severe local minima problems in case of the presence of strongly dispersive surface waves. Additionally, recorded wavefields are often characterized by amplitude errors due to varying source coupling and incorrect 3D-to-2D geometrical-spreading correction. Thus, least-squares FWI is considered less than suitable for near-surface applications. In this paper, we introduce an amplitude-unbiased coherency measure as a misfit function that can be incorporated into FWI. Such coherency was earlier used in phase-weighted stacking (PWS) to enhance weak but coherent signals. The benefit of this amplitude-unbiased misfit function is that it can extract information uniformly for all seismic signals (surface waves, reflections, and scattered waves). Using the adjoint-state method, we show how to calculate the gradient of this new misfit function. We validate the robustness of the new approach using checkerboard tests and synthetic data contaminated by random noise. We then apply the new FWI approach to a field dataset acquired at an archaeological site located in Ostia, Italy. The goal of this survey was to map the unexcavated archaeological remains with high-resolution. We identify a known tumulus in the FWI results. The instantaneous-phase coherency FWI results also establish that the shallow subsurface under the survey lines is quite heterogeneous. The instantaneous-phase coherency FWI of near-surface data can be a promising tool to image shallow small-scale objects buried under shallow soil covers, as found at archaeological sites.

@en

Seismic incoherent noise and waves scattered from objects in the crossline directions can cause 2D elastic full-waveform inversion (FWI) to produce artifacts in the resulting 2D models. We develop a complete workflow that can determine subsurface S-wave velocity (VS) models inverted from 2D near-surface seismic data more stably. We make use of a combination of supervirtual interferometry and a matched filter to accurately retrieve dominant surface waves from the field data, whereas the incoherent noise and 3D scattering events are significantly suppressed. The subsurface structures obtained from inverting the retrieved data can be interpreted together with the sections resulting from FWI of the original data to mitigate the potential misinterpretation of artifacts. Our results demonstrate that it is possible to invert 2D near-surface seismic data even when the data quality is lowered by the presence of strong noise and 3D scattered events caused by objects located in the crossline direction.

@en

Seismic interferometry (SI) is a method that retrieves new seismic traces from the cross-correlation of existing traces, where one of the receivers acts as a virtual seismic source whose response is retrieved at other receivers. When using sources only at the surface, and the so-called one-sided illumination of the receivers occurs, we will retrieve not only the desired physical reflections but also non-physical (ghost) reflections. These non-physical reflections appear due to waves that propagate inside a subsurface layer. Thus, they contain information about the seismic properties of the specific layer. We illustrate the technique’s potential using numerically modelled data for a subsurface model with a low-velocity layer, which is also pinching out, and near-surface field data. We apply SI by cross-correlation and auto-correlation. Both resulting non-physical reflections are sensitive to the physical rock properties of the layer that causes them to appear in the SI results. Moreover, non-physical reflections in zero-offset gathers that result from SI by auto-correlation show very good conformity with the geometry of the subsurface layers.@en

Difficulties in detecting and characterising shallow objects close the surface with seismic shear waves are often problematic because of dominant surface waves. By sequencing a specific combination of two data driven processing steps followed by diffraction tomography can overcome these problems. Small scattering objects become visible in the final image that can have importance of the understanding of subsurface locations, such as areas of archaeological interest. On the other hand, deep changes in the electric resistivity on land are often problematic to detect and especially to monitor time-lapse change over long periods of time. The usual electrodes slowly erode and vanish. Geothermal heat production environments often lead to changes in the resistivity between in-situ water-filled formations and cooler injected water-filled formations of less than one order of magnitude. A dedicated set of capacitively coupled electrode could overcome to erosion problem. When placed in a well with composite casing, these could be used in measurements of much enhanced detectability. In that case it is necessary to have electrodes in a zone from below to above the target layer. By changing the source offset at the surface, optimal measurements can be done to detect the small and deep changes in resistivity.@en

Induced seismicity from a gas-producing region such as Groningen is believed to be caused by reservoir depletion due to long-term gas production. However, because of the complexity and uncertainty regarding the underground structure and composition, it is difficult to quantify the effect on induced seismicity due to gas production. Here we use finite-element modelling to investigate the seismogenic potential of a pre-existing fault reactivated due to fluid depletion, considering different model settings. By applying quasi-static poroelastic loading representing reservoir depletion, the stress and strain fields are derived from the resulting displacement field. The equilibrium of the fault is then evaluated using either rate-and-state or slip-weakening behaviour for friction. When the critical state is reached on the fault, where the shear stress is greater than the friction, the reactivation of the fault takes place. This reactivation is simulated by using a dynamic solver to observe the propagation and the arrest of the dynamic faulting, as well as the resultant wavefield due to seismic slip. By comparing the depletion value at both aseismic and seismic ruptures, and looking at the stress distribution on the fault, the pattern of rupture nucleation, and the resulting seismic wavefield, we are able to evaluate separately the effect of different model settings, including the geometry and material property of both caprock and reservoir, reservoir depletion pattern, and the friction law. Furthermore, by combining our study with the observed seismic wavefield, it is possible to obtain useful insights on the spatial variation in the source region.@en

Induced seismicity from a gas-producing region such as Groningen is believed to be caused by reservoir depletion due to long-term gas production. However, because of the complexity and uncertainty regarding the underground structure and composition, it is difficult to quantify the effect on induced seismicity due to gas production. Here we use finite-element modelling to investigate the seismogenic potential of a pre-existing fault reactivated due to fluid depletion, considering different model settings. By applying quasi-static poroelastic loading representing reservoir depletion, the stress and strain fields are derived from the resulting displacement field. The equilibrium of the fault is then evaluated using either rate-and-state or slip-weakening behaviour for friction. When the critical state is reached on the fault, where the shear stress is greater than the friction, the reactivation of the fault takes place. This reactivation is simulated by using a dynamic solver to observe the propagation and the arrest of the dynamic faulting, as well as the resultant wavefield due to seismic slip. By comparing the depletion value at both aseismic and seismic ruptures, and looking at the stress distribution on the fault, the pattern of rupture nucleation, and the resulting seismic wavefield, we are able to evaluate separately the effect of different model settings, including the geometry and material property of both caprock and reservoir, reservoir depletion pattern, and the friction law. Furthermore, by combining our study with the observed seismic wavefield, it is possible to obtain useful insights on the spatial variation in the source region.@en

Detecting small-size objects is a primary challenge at archaeological sites due to the high degree of heterogeneity present in the near surface. Although high-resolution reflection seismic imaging often delivers the target resolution of the subsurface in different near-surface settings, the standard processing for obtaining an image of the subsurface is not suitable to map local diffractors. This happens because shallow seismic-reflection data are often dominated by strong surface waves that might cover weaker diffractions and because traditional common-midpoint moveout corrections are only optimal for reflection events. We propose an approach for imaging subsurface objects using masked diffractions. These masked diffractions are first revealed by a combination of seismic interferometry and nonstationary adaptive subtraction, and then further enhanced through crosscoherence-based supervirtual interferometry. A diffraction image is then computed by a spatial summation of the revealed diffractions. We use the phase-weighted stack to enhance the coherent summation of weak diffraction signals. Using synthetic data, we show that our scheme is robust in locating diffractors from data dominated by strong Love waves. We test our method on field data acquired at an archaeological site. The resulting distribution of shallow diffractors agrees with the location of anomalous objects identified in the VS model obtained by elastic SH/Love full-waveform inversion using the same field data. The anomalous objects correspond to the position of a suspected burial, also detected in an independent magnetic survey and corings.@en