Prestack imaging of compressed seismic data: A Monte Carlo approach
Descrição do Produto
Prestack imaging of compressed seismic data: a Monte Carlo approach
E. Bonomi, CRS4 Geophysics Area, Italy, L. Cazzola, ENI-AGIP Exploration and Production Division, Italy
Summary We analyze the depth wave eld extrapolation of a single dataset obtained by compressing sources and shotgathers. Our conclusion is that, running a sequence of M migrations, for a large class of random seismic data compressions, the noise variance of the average imaging condition decreases proportionally to M ,1 . This Monte Carlo approach may render 3D wave eld prestack migration cost e ective.
Introduction Prestack depth migration algorithms are usually divided into two classes: wave eld extrapolation techniques based on the implementation of a one-way wave equation solver [1] and Kirchho techniques using ray tracing as the wave propagation model [2]. Kirchho techniques, based on the eikonal approximation of the scalar wave equation associated with source and receiver location, is recognized as the most exible method of imaging prestack 3D seismic data. However, Kirchho migration, using rst arrival traveltime solutions, may encounter severe imaging inaccuracies in geological situations with fast lateral velocity variations. The use of a more complete description of the eikonal propagation, including multipathing, would increase the quality of migrated images at the cost of cumbersome and expensive ray tracing algorithms, which then could hardly handle caustics. Wave eld extrapolation techniques are well adapted to common-shot seismic data because in this domain they naturally split into two independent identical propagation problems, one related to the source wave eld and the other to the recorded wave eld, both followed by the standard imaging condition. Common-shot, wave equation migration has proven to be a very accurate imaging tool, either in the (x; y; !)-domain or in the (kx ; ky ; !)domain. Unfortunately the cost of prestack 3D wave eld extrapolation, especially with tens of thousands of shots, is prohibitive for industrial data set applications. To obtain a drastic reduction of the computational e ort, a new imaging technique based on the simultaneous wave eld extrapolation of many compressed shot-gathers was suggested in ref.[3]: sources and shot-gathers are phase encoded with a di erent sequence of random numbers. Our analysis shows that random seismic data compression and the resulting processing lead to a sort of Monte Carlo quadrature. In other words, the resulting imaging condition, averaged on M runs, each one representing a depth migration, provides, in addition to the expected correct contribution, a large number of spurious random
SEG 1999 Expanded Abstracts
terms that vanish in the limit of large M . Their stochastic variance is proportional to M ,1 . The method is valid for any sequence of random numbers, all independently drawn from the same probability distribution function with mean value zero and variance one. In presence of thousands of shots, to achieve an acceptable noise level, the number of runs does not exceed a small fraction of the number of shots.
Theory Migration is an inversion process which reconstructs the map R(x; z) of local re ectivity from the only available information, the velocity eld c(x; z) and the prestack seismic volume Q(xs ; xg ; t) = p(xs ; xg ; 0; t), where p(xs ; xg ; z; t) denotes the pressure measured at time t by a receiver at position (xg ; z ) after an impulsive source at position (xs ; z ) has been initiated at t = 0. In the case of a sequence of N independent shots, the seismic volume is written:
Q(xs ; xg ; t) =
X( N
n=1
xs , sn )Tr (n) (xg ; t):
(1)
is the source position at the n-th shot and T r(n) represents the resulting eld of recorded traces. As stated by seismic imaging theory [4], moving pairs of sources and receivers downward along ray paths, we nd that re ectors along points of discontinuity in the medium are eventually located at (xs ; z)=(xg ; z)=(x; z). As source and receiver approach each other, the travel time t between them goes to zero. On these points, the re ection coecient and the wave eld p are related as follows: sn
R(x; z ) p(x; x; z; t = 0):
(2)
The basic strategy for prestack wave eld migration is to formulate the propagation of p in the space-frequency domain. Linearity of the one-way equation and its formulation in the eld coordinate system allow the uncoupling of the downward extrapolation of the seismic volume, Eq.(1), in two depth wave eld propagations, one for the source, ps , and the other for the shot-gather, pg . The only coupling between them is introduced through the imaging condition, Eq.(2), from which the map of the local re ectivity is derived. For the sequence of N independent shots, the expression of R that results for sources and receivers inside the earth at each depth multiple of z takes the canonical form:
R(x; j z ) =
XR N
n=1
n (x; j z ):
( )
(3)
Main Menu
Monte Carlo Imaging
R n represents the re ectivity image of the earth's crust, and then the correct re ectivity of the medium. ( )
X^
resulting from the migration of the n-th shot-gather:
R n (x; j z ) = ( )
n (x; j z; !)^p(n) (x; j z; !) : (4) g
( )
ps
Methodology
To obtain numerical estimates of Eq.(7), a Monte Carlo importance sampling approach can be used [5]. More (n) p^s denotes the acoustic eld produced by the impul- precisely, we run a sequence of M depth migrations, sive source wavelet p(sn) (xs ; 0; t) = (xs , sn )(t), while starting each time from a new randomly compressed p^g denotes the depth extrapolation of the shot gather seismic volume, Eq.(5), and we construct the sequence (n) R1 ; R2 ; : : : ; RM of images to form the arithmetic mean pg (xg ; 0; t) = Tr (n) (xg ; t). Considering that in marine surveys N is of the order of 1 M R (x; j z ) : tens of thousands, Eqs.(3) and (4) eloquently illustrate the (8) A ( x ; j z ) = m M computational cost of prestack imaging. To drastically reM m =1 duce the processing, a solution is possible by compressing, in the space-frequency domain, sources and seismic data. AM is an unbiased estimator for hri. Monte Carlo theSo, a unique source term is composed by a linear com- ory states that, as M increases, the variance of the mean bination of all source wavelets and, similarly, a unique of AM decreases pointwise as M ,1. In other words, gather of seismic traces is assembled by superposing all value when the number M of migrations increases, the following shot gathers. The resulting compressed seismic volume is pointwise limit limM !1 AM = h r i almost surely exists. N The noise reduction reported in ref.[3] has now a precise Q^ (xs ; xg ; !) = a!;n (xs , sn )a!;n T^r(n ) (xg ; !); (5) statistical meaning. Our analysis actually explains why, increasing the number of stacked migrations, the phase n;n encoding de ned by the random variable a = exp(i ), where each coecient a!;n is a complex number. Eq.(5) with uniformly distributed on [0; 2), provides a better is the initial condition of the depth extrapolation based image. A second possibility would be, for instance, to on the one-way equation. The linearity of the one-way sample each an;! from a gaussian distribution with mean propagation leads to the following spurious expression for 0 and variance 1. A third choice would be to assign to the re ectivity each a one of the values 1 with probability 1=2. !
X
X
0
0
0
X )= N
X )+ N
n;!
Computer experiments show that these three compression r (x; j z r (n) (x; j z r (n;n ) (x; j z ); (6) methods, practically equivalent, are reliable in the sense n=1 n6=n that they provide a subsurface image that can be interThe accuracy is only controlled by the number M where, omitting the coordinates (x; j z), we have that preted. of stacked migrations that, for industrial applications, is much smaller than N , the number of shot-gathers. The (n) r = j a!;n j2 p^(sn) (; !) p^(gn) (; !) ; cost of migrating a single shot or a random superposition ! of many shots being the same, prestack imaging of compressed seismic data by the Monte Carlo approach may (n) (n;n ) (n ) r = a!;n a!;n p^s (; !) p^g (; !) : constitute a breakthrough in seismic processing. 0
0
X X
0
0
0
!
Suppose now that the collection an;! forms a suite of in- Examples dependent, identical random variables, then the stochastic To illustrate the method, we rst studied the imaging of a mean values of r(n) and r(n;n ) take the form: single point scatterer, depth 1200 m, in a 2D constant velocity eld, v = 4000m/s. The synthetic dataset (marine hr(n) i = h j a j2 i p^(sn) (; !) p^(gn) (; !) ; acquisition) is composed by 90 shot-gathers, with a far ! o set equal to 1000m, simulated using a nite-di erence (n) (n ) (n;n ) 2 p^s (; !) p^g (; !) ; hr i = jh a ij solver. The migration domain is 1400m 3000m, the recording time 4s, x = z = 10m and t = 4ms. ! (1) displays the scatterer image reconstructed by the 2 with hj a j i denoting the mean value of the squared mod- Fig. standard shot-gather processing (left). A commercial 2 ule of the random variable a, and jh aij denoting the code, implementing nite-di erence frequency-domain squared module of the mean value. Let us choose any extrapolation kernel,a was used for all the tests run in probability distribution such that hai = 0 and hj a j2 i = 1, this work. In a second computer experiment, a unique then, by computing the expectation of Eq.(6), we recover source term and a unique shot-gather of seismic traces the correct imaging condition, Eqs.(3) and (4), were composed by summing, respectively, all source wavelets and the 90 shot-gathers setting a!;n = 1 in N (5). Fig.(1) shows the resulting image (right): the h r i(x; j z ) = R(n) (x; j z ) ; (7) Eq. scatterer is completely delocalized and unrecognizable.
X X 0
0
0
X n=1
SEG 1999 Expanded Abstracts
Main Menu
Monte Carlo Imaging
Fig. 1: Scatterer images resulting from the standard shot-gather processing (left) and from the migration of one source and one shot gather, both compressed with no encoding (right).
Fig. 2: Scattererimages obtained from the depth extrapolationof one source and one shot-gather, both compressedimplementing Eq.(5) with the 1{random encoding: result after one migration (left) and result after 100 Monte Carlo migrations (right). Observe the drastic noise reduction, see also Fig.(3)
SEG 1999 Expanded Abstracts
Main Menu
Monte Carlo Imaging
Fig. 3: The inverse squared L2 {norm of the relative error is proportional to M , the number of migrations. As expected, Fig. 4: A more realistic example of subsurface image: 100 simthe three sampling methods are practically equivalent. ulated shot gathers migrated only 40 times. The compression was obtained by random phase encoding.
In a third computer experiment, a sequence of migrations with compressed sources and shot-gathers was run by implementing Eq.(5). To each coecient a!;n the value +1 or ,1 was assigned with probability 1=2. Fig.(2) displays, on the left, the scatterer image after only one migration: the correct structure is barely discernable through the noise. By stacking the sequence of images, averaging the imaging condition on M runs, the level of noise decreases in agreement with the analysis. This is illustrated by the behaviour in M ,1 of the squared L2 -norm of the difference between the correctly migrated image (standard shot-gather processing) and the average image de ned by Eq.(8). In Fig.(3), the inverse squared L2 -norm is displayed in arbitrary units. Fig.(2) shows on the right the resulting scatterer image after 100 migrations: the image is now interpretable. The same study was performed with the random phase encoding and the gaussian distribution for the coecients a!;n . Although Fig.(3) displays a better result for the random phase encoding, the three sampling methods provide very similar images. Fig.(4) illustrates a more realistic 2D synthetic example (marine acquisition): 100 shot gathers, migration domain 13000m4000m, far-o set 4000 m, recording time 4 s, x = 50; z = 10 m and t = 4 ms. The velocity eld is inhomogeneous and complex. The resulting image was obtained by random phase encoding, after only 40 migrations.
Conclusions Our main contribution, con rmed by computer experiments, states that, running a sequence of M migrations with compressed source and shot-gather, for a large class of random encoding, the variance of the average image error descreases as M ,1 . In other words, as long as the collection a!;n forms a suite of independent and identical random variables with mean value 0 and variance 1, it is always possible to extract meaningful information from
SEG 1999 Expanded Abstracts
compressed seismic data. This result does not depend on the size of the dataset. Computer experiments on industrial datasets corroborates all the potentiality of prestack Monte Carlo imaging which may ultimately render 3D wave eld imaging highly competitive with respect to the less accurate but cost e ective Kirchho technique.
Acknowledgements We wish to thank ENI-AGIP Exploration and Production Division and the Sardinian Regional Authorities for permission to publish.
References
[1] Bonomi E., Brieger L. M., Nardone C. M. and Pieroni E., 1998, PSPI: a scheme for high performance echo-reconstructive imaging, Computers in Physics, 12:126. [2] Garavaglia F., Bernasconi G., Drufuca, G. , 1998, 3D prestack migration of the SEG/EAGE overthrust model, 68th Ann. Internat. Mtg., Soc. Expl. Geophys, Expanded Abstracts. [3] Morton S. A., Ober C. C., 1998, Faster shot-record depth migrations using phase encoding, 68th Ann. Internat. Mtg., Soc. Expl. Geophys, Expanded Abstracts. [4] Stolt R. H., 1976, Migration by Fourier transform, Geophysics, 43:23. [5] M. H. Kalos and P. A. Whitlock, 1986, Monte Carlo methods, Volume 1: basics, John Wiley and Sons.
Main Menu
Lihat lebih banyak...
Comentários
