Danilo R. Velis

Stochastic sparse-spike deconvolution

Sparse-spike deconvolution can be viewed as an inverse problem where the locations and amplitudes of a number of spikes reflectivity are estimated from noisy data seismic traces. The main objective is to find the least number of spikes that, when convolved with the available band-limited seismicwaveletestimate,fitthedatawithinagiventolerance errormisfit.Thedetectionofthespikes’timelagsisahighly nonlinearoptimizationproblemthatcanbesolvedusingvery fastsimulatedannealingSA.Amplitudesareeasilyestimated using linear least squares at each SA iteration. At this stage, quadratic regularization is used to stabilize the solution, to reduce its nonuniqueness, and to provide meaningful reflectivitysequences,thusavoidingtheneedtoconstrainthe spikes’time lags and/or amplitudes to force valid solutions. Impedanceconstraintsalsocanbeincludedatthisstage,providing the low frequencies required to recover the acoustic impedance.Oneadvantageoftheproposedmethodoverother sparse-spike deconvolution techniques is that the uncertaintyoftheobtainedsolutionscanbeestimatedstochastically.Further,errorsinthephaseofthewaveletestimatearetolerated,foranoptimumconstant-phaseshiftisobtainedtocalibratetheeffectivewaveletthatispresentinthedata.Results using synthetic data including simulated data for the Marmousi2 modelandfield 3D data show that physically meaningful high-resolution sparse-spike sections can be derived fromband-limitednoisydata,evenwhentheavailablewaveletestimateisinaccurate.

Automatic first-breaks picking: New strategies and algorithms

We have developed three methods for the automatic pickingoffirstbreaksthatcanbeusedformarine,dynamite,orvibroseisshotrecords:amodifiedCoppens’smethod,anentropy-based method, and a variogram fractal-dimension method.Thetechniquesarebasedonthefactthatthetransitionbetweennoiseandnoiseplussignalcanbeautomaticallyidentified by detecting rapid changes in a certain attribute energy ratio, entropy, or fractal dimension, which we calculate withinmovingwindowsalongtheseismictrace.Theapplication of appropriate edge-preserving smoothing operators to enhancethesetransitionsallowedustodevelopanautomated strategy that can be used to easily signal the precise location ofthefirst-arrivalonset.Furthermore,weproposeamispickcorrectingtechniquetoexploitthebenefitsofthedatapresent in the entire shot record, which allows us to adjust the traceby-trace picks and to discard picks associated with bad or dead traces. As a result, the consistency of the first-break picks is significantly improved. The methods are robust under noisy conditions, computationally efficient, and easy to apply. Results using dynamite and vibroseis field data show that accurate and consistent picks can be obtained in an automated manner even under the presence of correlated noise, badtraces,pulsechanges,andindistinctfirstbreaks.

Simulated annealing wavelet estimation via fourth-order cumulant matching

The fourth-order cumulant matching method has been developed recently for estimating a mixed-phase wavelet from a convolutional process. Matching between the trace cumulant and the wavelet moment is done in a minimum mean-squared error sense under the assumption of a non-Gaussian, stationary, and statistically independent reflectivity series. This leads to a highly nonlinear optimization problem, usually solved by techniques that require a certain degree of linearization, and that invariably converge to the minimum closest to the initial model. Alternatively, we propose a hybrid strategy that makes use of a simulated annealing algorithm to provide reliability of the numerical solutions by reducing the risk of being trapped in local minima. Beyond the numerical aspect, the reliability of the derived wavelets depends strongly on the amount of data available. However, by using a multidimensional taper to smooth the trace cumulant, we show that the method can be used even in a trace-by-trace implementation, which is very important from the point of view of stationarity and consistency. We demonstrate the viability of the method under several reflectivity models. Finally, we illustrate the hybrid strategy using marine and field real data examples. The consistency of the results is very encouraging because the improved cumulant matching strategy we describe can be effectively used with a limited amount of data.

Minimum entropy deconvolution with frequency-domain constraints

A method for reconstructing the reflectivity spectrum using the minimum entropy criterion is presented. The algorithm (FMED) described is compared with the classical minimum entropy deconvolution (MED) as well as with the linear programming (LP) and autoregressive (AR) approaches. The MED is performed by maximizing an entropy norm with respect to the coefficients of a linear operator that deconvolves the seismic trace. By comparison, the approach presented here maximizes the norm with respect to the missing frequencies of the reflectivity series spectrum. This procedure reduces to a nonlinear algorithm that is able to carry out the deconvolution of band-limited data, avoiding the inherent limitations of linear operators. The proposed method is illustrated under a variety of synthetic examples. Field data are also used to test the algorithm. The results show that the proposed method is an effective way to process band-limited data. The FMED and the LP arise from similar conceptions. Both methods seek an extremum of a particular norm subjected to frequency constraints. In the LP approach, the linear programming problem is solved using an adaptation of the simplex method, which is a very expensive procedure. The FMED uses only two fast Fourier transforms (FFTs) per iteration; hence, the computational cost of the inversion is reduced.

Wavelet estimation revisited

This article deals with the ubiquitous problem of wavelet estimation. In essence, the problem consists of estimating, at best, two unknowns from one equation. In the seismic case, as it is well known, the model of the recorded seismic trace, xt, is generally formulated as a convolution, xt=wt*rt, where wt is the seismic source signature, or wavelet, and rt is the reflectivity series which is comprised of all primary and multiple reflections. Since, in general, we can assume that rt is white (or nearly so), the amplitude signature of wt may be estimated with some degree of confidence. It is, of course, the phase of wt that presents the difficulty. It turns out, perhaps somewhat counterintuitively, that, except for rather special circumstances, the phase of a signal uniquely determines the signal to within a scale factor. The amplitude spectrum plays no part whatsoever (we explore this further below). No wonder, then, that the phase is the all‐important factor.

Simulated annealing two-point ray tracing

We present a new method for solving the two-point seismic ray tracing problem based on Fermats principle. The algorithm overcomes some well known difficulties that arise in standard ray shooting and bending methods. Problems related to: (1) the selection of new take-off angles, and (2) local minima in multipathing cases, are overcome by using an efficient simulated annealing (SA) algorithm. At each iteration, the ray is propagated from the source by solving a standard initial value problem. The last portion of the raypath is then forced to pass through the receiver. Using SA, the total traveltime is then globally minimized by obtaining the initial conditions that produce the absolute minimum path. The procedure is suitable for tracing rays through 2D complex structures, although it can be extended to deal with 3D velocity media. Not only direct waves, but also reflected and head-waves can be incorporated in the scheme. One important advantage is its simplicity, in as much as any available or user-preferred initial value solver system can be used. A number of clarifying examples of multipathing in 2D media are examined.

Seismic characterization of thin beds containing patchy carbon dioxide-brine distributions: A study based on numerical simulations

We studied the seismic attenuation and velocity dispersion effects produced by wave-induced fluid flow in weakly consolidated sandstones containing patchy carbon dioxide CO2 ðÞ -brine distributions. The analysis also focuses on the velocity pushdown because of the presence of this gas, as well as on the role of the wave-induced fluid flow (mesoscopic) effects on the amplitude variation with angle (AVA) seismic response of thin layers containing CO2, such as those found at the Utsira Sand, Sleipner field, offshore Norway. We found that this loss mechanism may play a key role on conventional surface seismic data, suggesting that further data analysis may provide useful information on the characteristics of the fluid distributions in these environments. Numerical experiments let us observe that although mesoscopic effects can be very significant in this kind of media, the seismic response of a given isolated thin layer computed considering such effects is very similar to that of a homogeneous elastic thin layer with a compressional velocity equal to that of the original porous rock averaged in the effective data bandwidth. This suggests that the thin-bed prestack spectral inversion method published by the authors could be used to estimate representative compressional velocities and layer thicknesses in these environments. In effect, results using realistic synthetic prestack seismic data show that isolated CO2bearing thin beds can be characterized in terms of their thicknesses and representative compressional velocities. This information can be qualitatively related to CO2 saturations and volumes; thus, the prestack spectral inversion method could find application in the monitoring of the evolution of CO2 plumes at injection sites similar to that at the Sleipner field.

Thin-bed prestack spectral inversion

Prestack seismic data has been used in a new method to fully determine thin-bed properties, including the estimation of its thickness, P- and S-wave velocities, and density. The approach requires neither phase information nor normal-moveout (NMO) corrections, and assumes that the prestack seismic response of the thin layer can be isolated using an offset-dependent time window. We obtained the amplitude-versus-angle (AVA) response of the thin bed considering converted P-waves, S-waves, and all the associated multiples. We carried out the estimation of the thin-bed parameters in the frequency (amplitude spectrum) domain using simulated annealing. In contrast to using zero-offset data, the use of AVA data contributes to increase the robustness of this inverse problem under noisy conditions, as well as to significantly reduce its inherent nonuniqueness. To further reduce the nonuniqueness, and as a means to incorporate a priori geologic or geophysical information (e.g., well-log data), we imposed appropriate bounding constraints to the parameters of the media lying above and below the thin bed, which need not be known accurately. We tested the method by inverting noisy synthetic gathers corresponding to simple wedge models. In addition, we stochastically estimated the uncertainty of the solutions by inverting different data sets that share the same model parameters but are contaminated with different noise realizations. The results suggest that thin beds can be characterized fully with a moderate to high degree of confidence below tuning, even when using an approximate wavelet spectrum.

Estimating the distribution of primary reflection coefficients

The distribution of primary reflection coefficients can be estimated by means of the maximum entropy method, giving rise to smooth nonparametric functions which are consistent with the data. Instead of using classical moments (e.g. skewness and kurtosis) to constraint the maximization, nonconventional sample statistics help to improve the quality of the estimates. Results using real log data from various wells located in the Neuquen Basin (Argentina) show the effectiveness of the method to estimate both robust and consistent distributions that may be used to simulate realistic sequences.

Traveltime inversion for 2‐D anomaly structures

This work presents a traveltime inversion method that uses parametric functions to represent 2-D anomaly structures. These functions are described by a small set of unknown parameters which in turn are obtained after solving a highly nonlinear optimization problem via simulated annealing (SA). The procedure favors neither smooth nor high contrasting anomalies and keeps the number of unknowns very small so as to make the problem tractable using SA. Yet the strategy allows one to accommodate a large class of velocity models. Results indicate that this new approach typically yields better images than a standard linearized inversion based on a cell parameterization scheme.