Prestack exploding reflector modeling: The crosstalk problem

Claudio Guerra and Biondo Biondi

Abstract

The recently introduced prestack exploding reflector modeling aims to model a small dataset comprised by areal shots, while keeping the correct kinematics to be used in iterations of migration velocity analysis. To achieve this goal, the modeled areal data must be combined into sets. This procedure generates data which is subjected to crosstalk during migration. Here, using a simple constant velocity model, we describe two different origins for the crosstalks and. We show that by applying the concept of random phase-encoding during the modeling of the areal shots we achieve an almost complete elimination of the crosstalks.

INTRODUCTION

Oil exploration is being conducted in increasingly more geologically complex areas. The imaging challenges, along with the continuous growth of computer power, have augmented the relative importance of depth migration – in particular, wave-equation migration. A variety of wave-equation migration algorithms have been developed to balance the competing constraints of cost and final image quality. Velocity model definition for depth migration, however, strongly depends on ray-based methods because of their flexibility and low cost when compared with wave-equation migration velocity analysis. In areas where ray theory fails and the wavefields become too complex, this dependence might result in inaccurate velocity models that yield poor quality or unreliable images even when using state-of-the-art wave-equation algorithms. Thus, in these areas it is desirable to use wave-equation migration to generate image gathers in order to evaluate the accuracy of the velocity model.

Every iteration of migration velocity analysis for depth migration includes: 1) migrating with the current velocity model; 2) measuring properties that diagnoses the current velocity model accuracy; and, 3) updating the velocity model. Conventionally, Kirchhoff migration is used in step 1 and steps 2 and 3 are accomplished by measuring the departures of reflectors from the horizontal in the common image gathers and by back-projecting them along rays onto the model space, respectively. In areas where ray theory fails, multi-pathing and shadow zones occurs, causing problems in updating the velocity model using ray methods.

Wave-equation methods for both migration and velocity update overcome the limitations of ray-based methods (??), but at the expenses of a high cost. The practical implementation of wave-equation methods, therefore, still depends on an efficient migration/modeling scheme.

Guerra and Biondi  2  Prestack exploding reflector modeling

Based on the linearity of wavefield extrapolators, different strategies have been formulated to reduce the computational cost of wave-equation migration. The synthesis of plane waves was first introduced by (Schultz and Claerbout, 1978). (Rietveld et al., 1992) used pre-defined wavefields of arbitrary shape to synthesize areal shots after propagation back to the surface. (Rietveld and Berkhout, 1992) proposed a depth migration method using the concept of controlled illumination, which can be considered as a generalization of Schultz and Claerbout’s method. Migration of phase-encoded shot-records (Romero et al., 2000) is a technique that provides the ability to migrated combined shots with minimal cross-talk. The combination of phase-encoded and plane-wave areal shots can be used to drastically decrease the computational effort (Sun et al., 2002). (Liu et al., 2006) provide a general framework to evaluate plane-wave composition in prestack source plane-wave migration and show the equivalence of shot profile and plane-wave migrations.

Recently, (Biondi, 2006; ?) introduced the concept of the prestack exploding reflector modeling. This method synthesizes source and receiver wavefields along the entire survey at the surface, in the form of areal data, starting from a prestack migrated image cube represented by subsurface-offset domain common-image gathers (SODCIGs). For the case of migration velocity analysis, the aim is to generate a smaller dataset than the one used in the initial migration while maintaining the kinematics necessary to perform migration velocity analysis.

If the migration velocity is correct and assuming sufficiently good illumination, the reflectors are focused around the zero subsurface offset and the method reduces to the conventional exploding reflector modeling. However, when the migration velocity is inaccurate and in the presence of geologic dip, a pre-processing step is needed before modeling. (?) describes how to generate dip-independent image gathers to be used as initial condition for the upward propagation of source and receiver wavefields. This is achieved by rotating the SODCIGs according to the geological dip in the wavenumber domain.

Conceptually, in the prestack exploding reflector modeling, the synthesized areal data is computed by upward propagating source and receiver wavefields using as initial condition one SODCIG. Depending on the relationship between shot and SODCIG sampling, this procedure actually can generate as many individual experiments as the existing in the original data. However, we take advantage of the linearity of the wave propagation to combine several experiments into a set composite records, therefore decreasing the amount of data to migrate.

Combining several experiments, though, gives rise to crosstalk during imaging. (Biondi, 2006) uses a decorrelation distance between SODCIGs sufficiently large to prevent crosstalk. Using this criterion, however, limits the maximum data reduction to the ratio between the number of SODCIGs and the decorrelation distance. Furthermore, he mentions that an additional saving can be achieved by using concepts similar to phase encoding (Romero et al., 2000) during the modeling step. This allows us to combine more shots into one set. From now on, upward propagated areal shots combined into one single areal shot will be called a set.

(Biondi, 2006; ?) uses the two-way wave equation to propagate the wavefields. Here, we use the one-way wave equation to achieve additional computational saving. Furthermore, under the one-way framework, we may to introduce a phase-encoding-like scheme during the modeling to reduce cross-talk when applying the imaging condition.

PRESTACK EXPLODING REFLECTOR MODELING

The modeling of 2D source and receiver wavefields at the surface, using the one-way wave equation and starting from a prestack image at a selected position, x_{ξ}, can be described by:

\begin{matrix}S(x,ω) D G(z_{ξ},x_{ξ} − h_{ξ};x,z D 0,ω) ∗ I_{s}(z_{ξ},x_{ξ},h_{ξ}), \\ R(x,ω) D G(z_{ξ},x_{ξ} C h_{ξ};x,z D 0,ω) ∗ I_{r}(z_{ξ},x_{ξ},h_{ξ}), \end{matrix}  (1)

where S(x,ω) is the source wavefield at z D 0; R(x,ω) is the receiver wavefield at z D 0; I_{s}(z_{ξ},x_{ξ},h_{ξ}) and I_{r}(z_{ξ},x_{ξ},h_{ξ}) are the prestack images used as initial condition for the source and receiver wavefield extrapolation, respectively; G(z_{ξ},x_{ξ} ± h_{ξ};x,z D 0,ω) represents the one-way operator which extrapolates the wavefields from the subsurface to the surface; h_{ξ} is the subsurface offset; z_{ξ} is depth; ω is the temporal frequency and x is the spatial coordinate in the data space coinciding with x_{ξ}. Notice that for a perfectly focused image around h_{ξ} D 0, equation 1 reduces to the conventional zero-offset exploding reflector modeling. The prestack images used as initial condition for the source and receiver wavefield extrapolation are supposed to be dip-independent gathers computed by changing the dip along the offset direction according to the apparent geological dip (?).

The number of individual experiments described in equation 1 equals the number of shots used in the initial migration if the CMP spacing is equal to the shot spacing. Therefore, to decrease the amount of input data to migration, areal shots need to be combined into sets after the upward propagation. We achieve this by regularly selecting individual experiments and adding them up into their set, after being upward propagated, according to:

S̃_{n}(x,ω) D ∑_{nD1}^{k} ∑_{iDn,k,N}S_{i}(x,ω),

R̃_{n}(x,ω) D ∑_{nD1}^{k} ∑_{iDn,k,N}R_{i}(x,ω),  (2)

where S̃_{n}(x,ω) and R̃_{n}(x,ω) contains k−sets of summed areal sources and areal receivers, respectively; N is the number of SODCIGs. Every other k areal data are selected to compose one set . Pairs of S̃_{k}(x,ω) and R̃_{k}(x,ω) are to be used as the areal source function and the areal receiver wavefield, respectively, in areal shot migration.

As we show in the next section, unrelated events pertaining to the same pair of sets of areal source and areal receiver correlate together during imaging by cross-correlation, generating crosstalk. The crosstalk degrades the image, obscuring the kinematics we are interested in to perform migration velocity analysis.

CROSSTALK GENERATION AND ATTENUATION

Wavefield propagation is a linear process. This allows us to linearly combine wavefields before the propagation effects are removed by migration. Several migration methods explore this linearity to decrease the amount of data to migration (Rietveld and Berkhout, 1992; Romero et al., 2000; Liu et al., 2006). However, while applying the imaging condition, energy from unrelated source and receiver wavefields crosscorrelate generating crosstalk. For the present case, let us consider the migration of a set, m, of added areal shots containing a single pair of individual experiments. The prestack image is formed by crosscorrelating the source and receiver wavefields

Ĩ_{m}(z_{ξ},x_{ξ},h_{ξ}) D ∑_{ω} S̃_{m}^{⋆}(z_{ξ},x_{ξ} − h_{ξ},ω)R̃_{m}(z_{ξ},x_{ξ} C h_{ξ},ω),  (3)

where^{⋆} represents complex conjugation. If S̃_{m}(x,ω) and R̃_{m}(x,ω) in equation 2 are comprised by two summed areal shots, the image Ĩ_{m}(z_{ξ},x_{ξ},h_{ξ}) will be given by:

Ĩ_{m}(z_{ξ},x_{ξ},h_{ξ}) D I_{1}(z_{ξ},x_{ξ},h_{ξ}) C I_{2}(z_{ξ},x_{ξ},h_{ξ})C

∑_{ω} S_{1}^{⋆}(z_{ξ},x_{ξ} − h_{ξ},ω)R_{2}(z_{ξ},x_{ξ} C h_{ξ},ω) C

∑_{ω} S_{2}^{⋆}(z_{ξ},x_{ξ} − h_{ξ},ω)R_{1}(z_{ξ},x_{ξ} C h_{ξ},ω).  (4)

In equation 4, the last two summation terms represent the crosstalk. (Romero et al., 2000) attenuate these migration artifacts by multiplying the source and the receiver wavefields by a function of frequency and space. By doing so, the crosstalks are dispersed (random phase encoding) throughout the image after stacking over frequency, during imaging, or shifted outside of the image space (linear phase encoding). Here, we use a similar approach while modeling the areal shots.

To illustrate crosstalk problem in the prestack exploding reflector strategy let us make use of a simple constant velocity model of 2 km/s with two intersecting reflectors, one horizontal and the other dipping 15^{∘}. The dataset is comprised of 200 split-spread shots spaced every 20 m with a maximum offset of 2000 m. Figure 1 shows the shot-profile migration results using the correct velocity. The trace spacing in the migrated result is 20 m and the number of subsurface offsets in the prestack image is 41 sampled every 20 m. The front panel corresponds to the zero-subsurface offset section, the side panel is a SODCIG selected at x=1.5 km and the upper panel is a constant depth slice. We use this prestack image to model the areal data.

Because we performed migration with the correct velocity and the energy is concentrated at zero-offset, the pre-processing step of rotating the SODCIGs according to the geological dip was not necessary before modeling. Therefore, to model the areal source and receiver, SODCIGs are upward propagated using the one-way wave equation (equation 1) without any pre-processing step. The result of one modeling experiment is collected at the surface as an areal shot and assigned to its specific set.

Figures 2 and 3 show sets of time-domain areal sources and areal receivers computed with different distances between SODCIGs. Since the areal source wavefield is backward propagated it exists for negative times. For both figures, for the set of areal data labeled a) the distance between SODCIGs is 51 traces (1000 m), for the set labeled b) is 41 traces (800 m) and for the set labeled c) is 11 (200 m). These distances in number of traces defines the

Figure 1: ]

Prestack cube obtained with shot profile migration with the correct velocity.[ER] mig08 [ER]

Figure 2: ]

Sets of areal sources computed with different distances between SODCIGs. a) 51 traces; b) 41 traces; and c) 11 traces.[ER] arsource [ER]

Figure 3: ]

Sets of areal receivers computed with different distances between SODCIGs. a) 51 traces; b) 41 traces; and c) 11 traces.[ER] arrecev [ER]

number of sets to be migrated and, consequently, the migration efficiency. As we show next, this also influences the intensity of the crosstalk.

The crosstalk intensity on a subsequent migration depends on the distance between SOD-CIGs making up a set of areal data. (Biondi, 2006) shows that, using a decorrelation distance greater than twice the maximum subsurface offset, the crosstalks can be eliminated for the case of non-intersecting reflectors. For the case of reflectors with different geological dip and, moreover, intersecting, this procedure does not guarantee the absence of crosstalk.

Figures 4, 5 and 6 show the areal shot migration using cross-correlation imaging condition of the datasets containing the sets shown in Figures 2 and 3. The number of subsurface offsets is 41. Figure 4 shows the migration of 51 sets (Figure 2a and 3a) , Figure 5 shows the migration of 41 sets (Figure 2b and 3b) and Figure 6 shows the migration of 11 sets (Figure 2c and 3c). The migration of the data with 51 sets shows no crosstalk in the SODCIG because the distance of 51 traces is bigger than the decorrelation distance. The migration of the data with 41 sets shows an acceptable level of crosstalk in the SODCIG. However, the crosstalk is strong in the SODCIG of the data with 11 sets. Notice that for all of the three results a strong crosstalk occurs at zero-subsurface offset (x,z-domain).

The crosstalk described in the previous paragraph have two distinct origins. The one occurring in the zero-subsurface offset section as a reflector with an intermediate dip results from the croscorrelation of reflections in the source wavefield with reflections in the receiver wavefield from a single upward propagation experiment. The cross-correlation occurs at times different from zero propagation time. Therefore, one possible strategy is to use a time-windowed imaging condition about zero propagation time. The crosstalk present in non-zero offsets in the SODCIGs has similar origin as the one described by (?). In the case of the prestack exploding reflector modeling, they are related to reflections in the source and receiver wavefields pertain-

Figure 4: ]

Areal shot migration of 51 sets, each set containing the summation of areal shots initiated at every 51^{th}SODCIG. Original figure was substituted.[ER] armig51W [ER]

Figure 5: ]

Areal shot migration of 41 sets, each set containing the summation of areal shots initiated at every 41^{th}SODCIG. Original figure was substituted.[ER] armig41W [ER]

Figure 6: ]

Areal shot migration of 11 sets, each set containing the summation of areal shots initiated at every 11^{th}SODCIG. Original figure was substituted.[ER] armig11W [ER] ing to different areal data summed together into the same set. (?) uses an imaging condition that crosscorrelates decomposed source and receiver wavefields as a function of local slope at every position and time. This procedure yields good results, but at the expense of a higher computational cost.

The fact that the areal source wavefield also contains reflections originates crosstalk. According to the exploding reflector model (Loewenthal et al., 1976), reflectors explode at time zero. Therefore, the reflections in the areal source wavefield should be focused at zero time of wavefield propagation. We can use this information to restrict the cross-correlation in a timewindow about zero time of wavefield propagation. As (?) shows, the cross-correlation of areal wavefields within a small window around zero time of wavefield propagation can diminish the crosstalk at zero-subsurface offset. Equation 5 shows that, using the time-windowed imaging condition, crosstalk at the zero-subsurface offset only occurs between events having similar traveltimes.

Ĩ(z_{ξ},x_{ξ},h_{ξ} D 0) D I_{1}(z_{ξ},x_{ξ},h_{ξ} D 0) C I_{2}(z_{ξ},x_{ξ},h_{ξ} D 0)C

∑_{−1t≤tD0≤1t} S_{1}^{⋆}(z_{ξ},x_{ξ},t)R_{2}(z_{ξ},x_{ξ},t) C

∑_{−1t≤tD0≤1t} S_{2}^{⋆}(z_{ξ},x_{ξ},t)R_{1}(z_{ξ},x_{ξ},t).  (5)

Applying this imaging condition using one-way propagators in frequency domain, requires storing both wavefields in memory for all the frequencies, which increases the memory requirements. On the other hand, applying the imaging condition for a small number of time samples saves much computational effort. Figures 7, 8 and 9 show the areal shot migration using the time-windowed cross-correlation imaging condition from equation 4. Again, the number of subsurface offsets is 41. Figure 7 is the migration of 51 sets, Figure 8 is the migration of 41 sets and Figure 9 is the migration of 11 sets. Notice how the crosstalk in the zero-subsurface offset sections have been largely attenuated in each migration result.

Figure 7: ]

Areal shot migration using timewindowed imaging condition of 51 sets, each set containing the summation of areal shots initiated at every 51^{th} SODCIG. Compare with Figure 4. Original figure was substituted.[ER] armig51T [ER]

Figure 8: ]

Areal shot migration using timewindowed imaging condition of 41 sets, each set containing the summation of areal shots initiated at every 41^{th} SODCIG. Compare with Figure 5. Original figure was substituted.[ER] armig41T [ER]

Figure 9: ]

Areal shot migration using timewindowed imaging condition of 11 sets, each set containing the summation of areal shots initiated at every 11^{th} SODCIG. Compare with Figure 6. Original figure was substituted.[ER] armig11T [ER]

For the migration of Figure 9, the time-windowed imaging condition, however, does not attenuate the severe crosstalk in the SODCIG. Ideally, the areal shots pertaining to the same set should be uncorrelated. A possible way to decrease the correlation between areal shots is by encoding them with a random phase function, a(z_{ξ},x_{ξ},ω), in such a way that the crosstalk is dispersed throughout the image. As (Romero et al., 2000) show, for the random phase encoding we can choose a(z_{ξ},x_{ξ},ω) D e^{if (z_{ξ},x_{ξ},ω)} , where f (z_{ξ},x_{ξ},ω) is a random sequence as a function of z_{ξ},x_{ξ} and ω.

The modeling with phase encoding synthesizes data according to

\begin{matrix}S_{a}(x,ω) D a(z_{ξ},x_{ξ},ω)G(z_{ξ},x_{ξ} − h_{ξ};x,z D 0,ω) ∗ I_{s}(z_{ξ},x_{ξ},h_{ξ}), \\ R_{a}(x,ω) D a(z_{ξ},x_{ξ},ω)G(z_{ξ},x_{ξ} C h_{ξ};x,z D 0,ω) ∗ I_{r}(z_{ξ},x_{ξ},h_{ξ}), \end{matrix}  (6)

and the sets of random phase encoded areal sources and receivers, S̃_{n}^{a}(x,ω) and R̃_{n}^{a}(x,ω), respectively, are given by:

S̃_{n}^{a}(x,ω) D ∑_{nD1}^{k} ∑_{iDn,k,N}S_{i}^{a}(x,ω),

R̃_{n}^{a}(x,ω) D ∑_{nD1}^{k} ∑_{iDn,k,N}R_{i}^{a}(x,ω),  (7)

Figures 10a and 10b show a set of random phase-encoded areal source and receiver wavefields, respectively. Alternatively, we can encode just in x_{ξ} and ω. Figures 11a and 11b show a set of (x_{ξ},ω)-random phase-encoded areal source and receiver wavefields, respectively.

Figure 10: ] Areal source (a) and receiver (b) wavefields computed with random phase encoding, a(z_{ξ},x_{ξ},ω).[ER] arrphz [ER]

Figure 11: ] Areal source (a) and receiver (b) wavefields computed with random phase encoding, a(x_{ξ},ω).[ER] arrphw [ER]

The image, Ĩ_{m}(z_{ξ},x_{ξ},h_{ξ}), after areal shot migration of one set, m, comprised by two randomly phase-encoded areal sources and receivers is

Ĩ_{m}(z_{ξ},x_{ξ},h_{ξ}) D I_{1}(z_{ξ},x_{ξ},h_{ξ}) C I_{2}(z_{ξ},x_{ξ},h_{ξ})C

∑_{ω} a_{1}^{⋆}a_{2}S_{m}^{a_{1}^{⋆}}(z_{ξ},x_{ξ} − h_{ξ},ω)R_{m}^{a_{2}}(z_{ξ},x_{ξ} C h_{ξ},ω) C

∑_{ω} a_{2}^{⋆}a_{1}S_{m}^{a_{2}^{⋆}}(z_{ξ},x_{ξ} − h_{ξ},ω)R_{m}^{a_{1}}(z_{ξ},x_{ξ} C h_{ξ},ω).  (8)

The cross-correlation between the phase-encoding functions a_{1}(z_{ξ},x_{ξ},ω) and a_{2}(z_{ξ},x_{ξ},ω) must provide random values that disperse the crosstalk throughout the image.

Figures 12 and 13 show the areal shot migration results of 11 sets of randomly encoded areal data using encoding functions varying in (x_{ξ},ω) and (z_{ξ},x_{ξ},ω), respectively. The comparison between Figure 12 and Figure 6 shows that encoding only in (x_{ξ},ω) yields a significant reduction of the SODCIG crosstalk. The crosstalk in the zero-subsurface offset, however, is much less affected. The comparison between Figure 13 and Figure 12 shows that encoding in (z_{ξ},x_{ξ},ω) performs as efficiently as encoding in (x_{ξ},ω) in dispersing the SODCIG crosstalk. However, for the crosstalk in the zero-subsurface offset section, caused by cross-correlation of reflections in the areal source and receiver wavefields originated from the upward propagation of the same SODCIG, encoding in (z_{ξ},x_{ξ},ω) performs much better than encoding in (x_{ξ},ω). This is because the additional encoding in (z_{ξ} decorrelates the exploding reflectors.

Since the number of phase-encoded areal data is small, the noise dispersed through the image is not negligible. By migrating more sets of (z_{ξ},x_{ξ},ω)-randomly phase-encoded areal data, these amplitudes can be largely attenuated. Figure 14 shows the superior result when compared to Figures 6, 12 and 13, obtained by migrating 4 groups of 11 sets of (z_{ξ},x_{ξ},ω)- randomly encoded areal data. Besides virtually eliminating the crosstalk, the speckled noise has lower amplitudes than in Figure 13.

Figure 12: ]

Areal shot migration of randomly encoded data. 11 sets of (x_{ξ},ω)– randomly encoded areal data were migrated. Compare with Figure 6. Original figure was substituted.[ER] arm11phw [ER]

Figure 13: ]

Areal shot migration of randomly encoded data. 11 sets of (z_{ξ},x_{ξ},ω)– randomly encoded areal data were migrated. Compare with Figure 6. Original figure was substituted.[ER] arm11phz [ER]

One drawback of applying random phase encoding is that it makes the reflectors randomly “explode” at times different from the zero time of the wavefields. Consequently, the timewindowed imaging condition can not be applied. Perhaps, the linear phase encoding will enable us to use the time-windowed imaging condition. This deserves future research.

CONCLUSIONS

The prestack exploding reflector method potentially can decrease data volume required to perform migration velocity analysis. To be cost efficient, the method relies on the linearity of

Figure 14: ]

Areal shot migration of randomly encoded data. 4 groups of 11 sets of (z_{ξ},x_{ξ},ω)–randomly encoded areal data were migrated. Compare with Figures 6, 12 and 13. Original figure was substituted.[ER] arm11phz4 [ER]

the wave-propagation to combine modeled areal shots into areal data super-sets . This process gives rise to crosstalk when migrating the areal data.

We show that applying a time-windowed imaging condition efficiently mitigates the crosstalk in the zero-subsurface offset section, but does little to attenuate SODCIG crosstalk. We present a strategy to perform random phase encoding during the modeling. The results are promising and much better than the time-windowed imaging condition with respect to attenuating SOD-CIG crosstalk. The next steps will be the application of this strategy on data migrated with an inaccurate velocity model as well as on data from complex geology.

Because of the random encoding, the reflectors do not focus anymore about the zero time of the wave propagation, disabling the use of the time-windowed imaging condition. Further research will investigate the use of linear phase encoding associated to the application of this imaging condition.

REFERENCES

Biondi, B., 2006, Prestack exploding-reflectors modeling for migration velocity analyis: SEP– 124, 45–60.

Liu, F., D. W. Hanson, et al., 2006, Toward a unified analysis for source plane-wave migration: Geophysics, 71, 129–139.

Loewenthal, D., L. Lu, R. Roberson, and J. Sherwood, 1976, The wave equation applied to migration: Geophys.Prosp., 24, 380–399.

Rietveld, W. E. A. and A. J. Berkhout, 1992, Prestack depth migration by means of controlled illumination: Geophysics, 59, 801–809.

Rietveld, W. E. A., A. J. Berkhout, and C. P. A. Wapenaar, 1992, Optimum seismic illumination of hydrocarbon reservoirs: Geophysics, 57, 1334–1345.

Romero, L., D. Ghiglia, C. Ober, and S. Morton, 2000, Phase encoding of shot records in prestack migration: Geophysics, 65, 426–436.

Schultz, P. S. and J. Claerbout, 1978, Velocity estimation and downward continuation by wavefront synthesis: Geophysics, 43, 691 714.

Sun, P., S. Zhang, and F. Liu, 2002, Prestack migration of areal shot records with phase encoding: 72nd Ann. Internat. Mtg, Soc. Expl. Geophys., Expanded Abstracts, pages 1172–1175.