User:Ageary/Yilmaz 10

From SEG Wiki
Jump to: navigation, search
Seismic Data Analysis
Seismic-data-analysis.jpg
Series Investigations in Geophysics
Title Seismic Data Analysis
Author Öz Yilmaz
DOI http://dx.doi.org/10.1190/1.9781560801580
ISBN ISBN 978-1-56080-094-1
Store SEG Online Store


Contents

10.0 Introduction

A recorded seismic wavefield represented by a shot gather has two components — traveltimes and amplitudes. Direct inversion of a seismic wavefield to estimate elastic parameters of the earth demands numerically intensive computations. Instead, most practical methods of inversion are applied to seismic traveltimes and amplitudes, separately. In earth modeling in depth, we discussed methods of layer velocity estimation — Dix conversion, stacking velocity inversion, and coherency inversion, which essentially are based on inversion of traveltimes. In earth imaging in depth, we discussed prestack Kirchhoff migration, which again, essentially is based on computing diffraction traveltimes. Traveltime inversion thus yields a structural model of the earth represented by a set of layer velocities and reflector geometries, which can then be used to derive a structural image of the earth by depth migration. The term structural inversion may be appropriately used to describe the process of structural modeling and imaging by way of inversion of traveltimes.

In reservoir geophysics, we shall discuss poststack amplitude inversion to estimate an acoustic impedance model of the earth and prestack amplitude inversion to derive the amplitude variation with offset (AVO) attributes. Amplitude inversion thus yields a stratigraphic model of the earth represented by a combination of acoustic impedance and AVO attribute changes within the layers themselves. The term stratigraphic inversion may be appropriately used to describe the processes of estimating the acoustic impedance and AVO attributes by way of inversion of amplitudes.

In this chapter, we shall discuss case studies in structural inversion of 2-D and 3-D seismic data. These case studies relate to structural complexities caused by

  1. extensional tectonism as in the cases of salt diapirs of the North Sea and the Gulf of Mexico,
  2. compressional tectonism as in the cases of overthrust belts of the Middle East and Rocky Mountains, and
  3. wrench tectonism as in the cases of the pull-apart basins of offshore Indonesia and Venezuela.

Results from inversion of a 2-D seismic data set must be evaluated within the bounds of 2-D imaging, and 3-D effects must always be kept in mind (introduction to earth imaging in depth). A structural inversion project can be a futile exercise if the 2-D data set has not been recorded along a dominant dip direction with minimal 3-D effects.

The first 2-D case study is from the Southern Gas Basin of the North Sea. The North Sea Basin has been subjected to extensional tectonics, primarily in the northeast-southwest direction with some rotational component. As a result, salt diapiric structures were formed and the overlying strata were subjected to faulting. Continuing extensional tectonics and salt movements caused further faulting of the overlying strata, thus forming collapsed structures especially above the apexes of the salt diapirs. Common structural targets in the North Sea are Permian sands of the Rotliegendes and Carboniferous substrata below the Zechstein diapiric formation. We want to obtain an accurate depth structure map of the top Rotliegendes formation in the survey area. This requires an estimate of the velocity-depth model above the Zechstein diapiric formation and removal of its deleterious effect on the underlying Permian sands of Rotliegendes and deeper targets.

The second 2-D case study is from the Gulf of Mexico. The basin in the Gulf of Mexico has been subjected to extensional tectonics, primarily in the north-south direction. As a result, the Jurassic salt loaded by the overburden sand-shale sequence began to be deformed first in the edges of the basin, forming diapiric structures. The extensional tectonism also gave rise to a series of growth faults, often enveloped by major listric faults and sometimes accompanied by counter-regional faults. As the sand-shale deposition of Miocene and subsequent ages increased the overburden pressure, some diapirs took dike-like forms and some formed overhangs. Additionally, some of the salt diapirs were squeezed along the fault planes and some even were moved laterally as far up as the water bottom, eventually detaching themselves from the source salt layer situated within the deeper strata and forming tabular salt bodies. Some of these tabular bodies joined together to form salt canopies. The subsequent deposition on top then deformed the shape of the salt canopies, giving rise to a rugose geometry along the top and the base. The target zones are, for some cases, immediately below the base-salt, and for some other cases, somewhat deeper, they are below an overpressured zone.

The third 2-D case study deals with irregular water-bottom topography associated with a reef body. The reef causes severe distortions of reflection travel-times associated with underlying target horizons. The fourth 2-D case study deals with imaging beneath a volcanic layer. Finally, the fifth 2-D case study is from the Gulf of Thailand. Shallow gas anomalies cause difficulty in imaging multileveled reservoirs along faults that are abundant in the area.

The first 3-D case study is from the Southern Gas Basin of the North Sea. The objective is to delineate the geometry of the base-Zechstein layer which forms a three-dimensionally complex diapiric structure. The second 3-D case study is from the Central North Sea. We demonstrate the use of the combination of the inversion methods listed in Table 9-1 to estimate layer velocities and delineate reflector geometries within the overburden above the Zechstein formation. The third 3-D case study is from offshore Indonesia. The objective is to image the complex fault blocks caused by pull-apart tectonism. In this case study, we also demonstrate structural and stratigraphic interpretation of the image volume derived from 3-D prestack depth migration. Finally, the fourth 3-D case study is from Northeast China. With this case study, we demonstrate a complete time-and-depth sequence for processing, inversion, and interpretation that involves both 3-D prestack time and depth migrations.

10.1 Subsalt imaging in the north sea

The first 2-D case study for structural inversion is from the Southern Gas Basin of the North Sea. Figures 10.1-1 and 10.1-2 show selected CMP gathers from the two segments of the line accompanied by the corresponding CMP stacked sections. The left segment (Figure 10.1-1) contains a salt diapir with relatively steep flanks and the right segment (Figure 10.1-2) contains a salt diapir which has a broader base. The Zechstein formation comprises a halite unit with a velocity of 4400 m/s and anhydrite-dolomite rafts of various sizes and shapes with a velocity of 5900 m/s. The strong velocity contrast across the top-salt boundary and the presence of the anhydrite-dolomite rafts within the diapiric formation itself give rise to raypath distortions and thus constitute a complex overburden structure. The effect of the complex overburden on the underlying target zone — base Zechstein and the underlying Carboniferous sequence, is evident on the CMP gathers and the CMP stacks. Specifically, note events with complex moveout on CMP gathers 481-881 below 1.5 s (Figure 10.1-1) and the traveltime distortions along the base-Zechstein reflection on the CMP stack (events A, B, C, and D in Figures 10.1-1, and events E, F, G, H, and K in Figure 10.1-2).

The subsurface geology in the Southern Gas Basin can be represented by a layered earth model that consists of distinct layer boundaries with significant velocity contrast. Also, there exist vertical velocity gradients in layers above the salt formation. A typical velocity-depth model for the Southern Gas Basin may be considered in two parts:

  1. An overburden above the salt layer with some faulting where time migration is applicable, and ray-paths associated with nonzero-offset traveltimes yield a moveout behavior that makes it possible to estimate layer velocities with sufficient accuracy using Dix conversion, stacking velocity inversion or coherency inversion, and
  2. a substratum that includes the salt layer and the layers beneath, where depth migration is imperative.

The boundary between the overburden and the substratum is defined by the top of the salt layer where the most severe ray bending takes place. To estimate a velocity-depth model for a typical Southern Gas Basin subsalt target, the following procedure composed from the list of inversion methods in Table 9-1 can be used:

  1. Coherency inversion to estimate layer velocities and 2-D poststack depth migration to delineate reflector geometries within the overburden, and
  2. constant half-space velocity analysis of image gathers from prestack depth migration to estimate the substratum velocity and stacking of image gathers to delineate the reflector geometry of the base-Zechstein (top-Rotliegendes) target horizon.

This procedure is applied layer-by-layer starting from the surface to resolve layer velocity and reflector geometry of one layer before moving onto the next (model building). Such an approach enables us to establish the accuracy of the model one layer at a time and minimizes accumulation of errors as we proceed down to the target zone.

Estimation of the overburden model

Start the analysis by interpreting the time horizons from the unmigrated CMP stack that correspond to layer boundaries with significant velocity contrast. These horizons are denoted in Figure 10.1-3 — water bottom (H1), base Miocene Unconformity (H2), base Upper Tertiary (H3), base Lower Tertiary (H4), base Cretaceous Chalk (H5), base Upper Triassic (H6a), base Lower Triassic (H6b), and base Zechstein (top Rotliegendes) (H7). The time horizons are assumed to be equivalent to zero-offset reflection times that are used in coherency inversion (models with horizontal layers).

Horizon H6b is the top-salt boundary, which separates the overburden and substratum parts of the model. Although interpreted, horizon H7 — base-salt boundary, is not included in the analysis sequence for modeling the overburden. Instead, it is dealt with as part of modeling the substratum.

The depth horizon associated with the water bottom is obtained simply by normal-incidence depth conversion of the time horizon (H1 in Figure 10.1-3). Then, the following sequence was applied to horizons H2-H6b, one layer at a time, starting at the top. For the sake of the discussion here, assume that the velocity-depth model for the first n − 1 layers already have been determined, and that the nth layer is under consideration.

  1. Perform coherency inversion along the horizon under consideration, pick semblance maxima, and derive an interval velocity profile as a function of the midpoint location along the line. Available velocity gradient information is incorporated into the interval velocity estimation. As demonstrated in models with horizontal layers, the accuracy in velocity estimation degrades with shorter effective cable length, faster velocity and deeper horizon. Spurious peaks on semblance curves and rapid lateral variations in velocity should be avoided. Accordingly, smoothing is applied to the velocity profile as much as geologically plausable, but not excessively so as to retain lateral velocity variations that are realistic.
  2. Create a gridded velocity-depth model that consists of two parts — the known part on top, with the n − 1 layers already established, and the unknown part underneath, defined as a half-space with its velocity equal to the velocity of the nth layer derived in step (a).
  3. Perform poststack depth migration using the gridded velocity-depth model from step (b) down to a depth just below the layer of interest.
  4. Interpret the depth horizon associated with the base of the layer under consideration from the depth-migrated section.
  5. To verify the accuracy of velocity estimation from coherency inversion and, if needed, to update the layer velocity (model updating) derived in step (a), perform prestack depth migration to create image gathers at some interval along the line. The velocity-depth model used for prestack depth migration is the same as that used for poststack depth migration in step (d).

Repeat steps (a) through (e) for all the layers within the overburden down to top-Zechstein boundary (Horizon 6b), and thus establish a velocity-depth model for the overburden (Figure 10.1-4).

Estimation of the substratum model

We now want to estimate a velocity field for the substratum region — Zechstein and the underlying strata. Again, consider the velocity-depth model in two parts — the overburden, the known part, and the substratum, the unknown part, which is defined as a half-space.

  1. Shown in Figure 10.1-4 are four different velocity-depth models with the same overburden but different constant velocities assigned to the substratum. These constant velocities span the range that correspond to the velocity variations within Zechstein. The anhydrite-dolomite rafts (5900 m/s) floating within the halite formation (4400 m/s) have little thickness (less than 100 m), but can influence the velocity of the Zechstein unit as a whole depending on the spatial distribution of the rafts.
  2. Perform prestack depth migration using the four earth models in Figure 10.1-4, and generate a set of image gathers along the line (Figures 10.1-5 and 10.1-6). Examine the image gathers from the four velocity panels at the same location and note that they all represent the same image associated with the overburden down to a depth that corresponds the the top-salt boundary. Nevertheless, the event associated with the base Zechstein exhibits different moveout characteristics depending on the half-space velocity. For instance, at midpoint location 321, flatness for the base-Zechstein event at a depth of 3500 m is achieved with the 4700-m/s substratum velocity (event A in Figure 10.1-5). Similarly, at midpoint location 1681, flatness for the base-Zechstein event at a depth of 3950 m is achieved, again, with the 4700-m/s substratum velocity (event B in Figure 10.1-6). In principle, one may be able to determine velocity nodes from these image gathers based on the flatness of the event that corresponds to the base-Zechstein boundary. It is not just flatness, but event strength and continuity across the offset axis of the image gather, that we take into consideration when picking velocity nodes. If the event is weak, it could mean that the velocity assigned to the layer above is erroneously too low or too high, causing migration errors. To gain confidence in the velocity determination, the image gathers are used in combination with the image-gather stacks as in the next step.
  3. Stack the image gathers as from step (b) to generate the depth images (Figures 10.1-7 and 10.1-8) associated with the earth models in Figure 10.1-4. The stack power associated with the base-Zechstein event can be used as an additional criterion in combination with the flatness criterion for image gathers (Figures 10.1-5 and 10.1-6) to derive the velocity profile for the Zechstein formation.
  4. Assign the velocity field for the substratum derived from step (c) to the half-space below the overburden. Then, combine this velocity field for the half-space with the overburden model to construct a new velocity-depth model.
  5. Next, perform prestack depth migration using the new velocity-depth model from step (d).
  6. Interpret the depth image from step (e) to delineate the base-Zechstein horizon. Subsequently incorporate this horizon into the velocity depth-model, and assign a velocity of 5400 m/s to the subsalt region based on well data to construct a final velocity-depth model (Figure 10.1-9).
  7. Perform prestack depth migration once more using the final velocity-depth model from step (f) to obtain the image gathers and their stacks which represent the final image in depth.
  8. Convert the depth-migrated sections from depth to time using the final velocity-depth model from step (f) to apply poststack spiking deconvolution, band-pass filtering and AGC scaling. Finally, convert the sections back to depth. It usually helps improve the vertical resolution of depth-migrated data to apply poststack deconvolution. This also is needed for a fair comparison of images from prestack depth migration and poststack depth migration, since the latter normally would be performed using stack with poststack processing applied.

Selected image gathers and the depth images using the final velocity-depth model (Figure 10.1-9) from the two segments of the line are shown in Figures 10.1-10 and 10.1-11. Compare the depth image obtained from prestack depth migration with the images obtained from poststack depth and time migrations (Figures 10.1-12 and 10.1-13). We should not expect significant differences between the three images within the overburden region where time migration often yields acceptable results. Note, however, the differences in the substratum region. Specifically, the base-Zechstein horizon (event A in Figure 10.1-10b and events E and F in Figure 10.1-11b) cannot be delineated from poststack depth migration, nor can it be delineated accurately from poststack time migration (events G and H in Figure 10.1-12b, and event K in Figure 10.1-13). The anhydrite-dolomite rafts (such as event D in Figure 10.1-11b) and the 3-D behavior of the diapiric structure limit the accuracy of the final image obtained from prestack depth migration. Specifically, note the poor image of the base Zechstein represented by events B and C in Figure 10.1-10b.

Image rays associated with the base-Zechstein boundary (Figure 10.1-14), however, clearly indicate the need for prestack depth migration for accurate imaging of the substratum region. Note the significant lateral shifts on the image rays, especially beneath the salt diapirs — demonstrative evidence of the presence of strong lateral velocity variations associated with a complex overburden.

Model verification

The final stage in earth modeling and imaging is the verification of the accuracy of the model itself (model updating). For an earth model in depth to be acceptable, it has to pass the following two tests:

  1. Image gathers from prestack depth migration using the earth model in question must exhibit flat events (Figures 10.1-10a and 10.1-11a). Events associated with multiples and converted waves (4-C seismic method) are not expected to be flat. Nevertheless, even with good models, usually, there also are some primary events that do not meet the flatness criterion. Violation of this criterion may occur because of erroneous layer velocities or 3-D effects that are not accounted for by 2-D modeling and imaging. Other sources of departures from flatness include strong lateral velocity variations that are much less than a cable length and effect of anisotropy on layer velocities.
  2. For an earth model in depth to be an acceptable representation of the subsurface geological model, it must be consistent with the seismic data used to estimate the model in question. To check for consistency, perform ray-theoretical modeling of zero-offset traveltimes associated with the layer boundaries included in the model (Figure 10.1-14). Then, superimpose the modeled traveltimes on the CMP-stacked section and observe any discrepancy between the modeled and the actual traveltimes (Figure 10.1-15). Here, we are assuming that the reflection traveltimes observed on a CMP-stacked section can be closely approximated by two-way zero-offset traveltimes. Actual traveltimes interpreted from the unmigrated CMP-stacked section are shown in Figure 10.1-3.

Aside from these two criteria, the estimated earth model needs to be validated by examining it for consistency with the structural model applicable to the area of interest. For instance, the faulting and folding implied by the the model must have the same patterns as in the true subsurface situation. Calibration to well tops is also part of the model validation procedure (model building).

10.2 Subsalt imaging in the gulf of mexico

The second 2-D case study for structural inversion is from the Gulf of Mexico. Figure 10.2-1 shows the DMO-stacked section and Figure 10.2-2 shows the poststack time-migrated section of a 2-D data set from the Gulf of Mexico. The high-amplitude event (such as event A in Figure 10.2-1) with complex traveltime is the top-salt reflection. Conflicting dips associated with the fault blocks within the overburden and the rugose top-salt boundary are preserved by way of DMO correction (Figure 10.2-1), and accurate imaging of the suprasalt region can be achieved by poststack time migration (Figure 10.2-2). Nevertheless, accurate imaging of the base-salt boundary and the subsalt region is only possible by way of prestack depth migration.

Aside from the water layer, a Gulf of Mexico velocity-depth model is typically represented in two parts:

  1. A background velocity field with vertical velocity variations characterized by gentle variations in the gradient, the absence of distinct layer boundaries, and mild-to-moderate lateral velocity variations.
  2. Tabular and diapiric salt bodies with different shapes, but with a constant velocity of 4450 m/s, embedded into the background velocity field.

To estimate a velocity-depth model for a Gulf of Mexico structural target below the tabular salt bodies, the following procedure composed from the list of inversion methods in Table 9-1 is used:

  1. Dix conversion of stacking velocities to estimate the background velocity field,
  2. Model updating and verification of the velocity field within the suprasalt region (model updating),
  3. Poststack or prestack depth migration to delineate the top-salt boundary,
  4. Assignment of the salt velocity into the half-space below the top-salt boundary,
  5. Prestack depth migration to delineate the base-salt boundary,
  6. Assignment of the the background velocity into the half-space below the base-salt boundary (the subsalt region),
  7. Prestack depth migration to obtain and verify the final earth image in depth.

In this case study, we shall first test the procedure used for the Southern Gas Basin line to estimate a layered earth model in depth for comparison with the procedure outlined above.

Layered earth model estimation

Start the North-Sea style analysis (subsalt imaging in the North Sea) by interpreting a set of time horizons from the time-migrated CMP stack (Figure 10.2-3a). These time horizons are then used to obtain the unmigrated time horizons (Figure 10.2-3b) by way of a ray-theoretical forward modeling scheme, sometimes called demigration. Horizon A in Figure 10.2-3a represents the top-salt boundary. The modeled time horizons are assumed to be equivalent to zero-offset reflection times, which are used in coherency inversion. Alternatively, time horizons could have been interpreted from the unmigrated stacked section, directly.

To estimate a layered earth model, the following procedure was implemented:

  1. Starting from the top, for each layer within the suprasalt region, perform coherency inversion to estimate the layer velocity, and
  2. 2-D poststack depth migration to delineate the reflector geometry associated with the layer boundary.
  3. Repeat steps (a) and (b) for all the layers within the suprasalt region, and thus establish the velocity-depth model for that region.
  4. Assign the salt velocity into the half-space below the top-salt boundary, and
  5. perform 2-D prestack depth migration to delineate the base-salt boundary.
  6. Then, assign a vertically varying velocity function into the subsalt region (sometimes referred to as salt flooding), and
  7. perform 2-D prestack depth migration to obtain and verify the final earth image in depth.
  8. Convert the image-gather stack from step (g) from depth to time using the velocity-depth model from step (d), apply poststack spiking deconvolution, band-pass filtering, and AGC scaling, then convert back to depth.

Figure 10.2-4 shows the velocity-depth model based on the procedure outlined above. Note that there are three detached salt bodies. Also note the lateral velocity variations, detected by coherency inversion, within each layer in the suprasalt region. The velocity field in the subsalt region is referenced to the water bottom, since the salt bodies do not have much influence on the vertical variations in the background velocity field. To appreciate the raypath distortions caused by the salt bodies, examine the image rays down to the base-salt boundary in Figure 10.2-5. The image rays clearly indicate the need for prestack depth migration for accurate imaging of the subsalt region.

Figure 10.2-6 shows the depth image from prestack depth migration using the velocity-depth model in Figure 10.2-4. Note the rugose top-salt boundary, and relatively smoother base-salt boundary. The accompanying image gathers are shown in Figure 10.2-7. Note the high-amplitude events (A and B) associated with the top- and base-salt boundaries. Flatness of events both in the suprasalt and subsalt regions provides the evidence that the model in Figure 10.2-4 is geologically acceptable. Note the presence of intrasalt events (C) with moveout; they often are associated with shale intrusions. These events could also be associated with multiples or converted waves (4-C seismic method).

As part of a verification procedure, the estimated earth model (Figure 10.2-4) needs to be tested for consistency with the input seismic data. Figure 10.2-8a shows the modeled zero-offset traveltimes superimposed on the unmigrated stacked section as in Figure 10.2-1. These traveltimes are associated with the layer boundaries included in the velocity-depth model of Figure 10.2-4. They should be compared with the traveltimes derived from the forward modeling of the time horizons picked from the time-migrated section (Figure 10.2-8b). Note that the two sets of traveltimes are fairly consistent within the suprasalt region. However, some discrepancy exists for the top-salt event (red horizon) and the base-salt event (yellow horizon).

Compare the depth images derived from prestack depth migration (Figure 10.2-6) and poststack depth migration (Figure 10.2-9) using the same velocity-depth model (Figure 10.2-4). While the images within the suprasalt region are comparable, prestack depth migration certainly yields a superior image of the base-salt boundary and the subsalt region.

Structure-independent model estimation

We shall now implement the procedure outlined at the beginning of this section.

  1. Start with the stacking velocity field (Figure 10.2-10). Note that already a pattern is inferred about the velocity-depth model — an upper region where the velocity field shows a consistent increase with depth and moderate lateral variations, and a lower region where anomalous zones exist. Actually, the upper region coincides with the suprasalt region and the lower region coincides with the salt sills and the subsalt region.
  2. Apply smoothing to the stacking velocity field in Figure 10.2-10 to derive what may be considered an rms velocity field (Figure 10.2-11) for Dix conversion.
  3. Perform Dix conversion to derive an initial velocity-depth model in the form of an interval velocity field (Figure 10.2-12). The upper half of this model appears geologically plausable — a velocity field typical of a sedimentary sequence without any distinct layer boundaries with significant velocity contrast. The lower half exhibits anomalous behavior that geologically is not meaningful. It is this region where Dix conversion would fail because of the presence of salt sills with rugose boundary that causes severe raypath distortions.
  4. Using the initial velocity-depth model from step (c) (Figure 10.2-12), perform prestack depth migration and generate image gathers (Figure 10.2-13). Flatness of events indicates that the initial velocity-depth model (Figure 10.2-12) based on Dix conversion of rms velocities is remarkably accurate within the suprasalt region (above event A which corresponds to the top-salt boundary).
  5. Nevertheless, small adjustments may be required to the velocity-depth model within the suprasalt region before moving down to the salt and subsalt regions. To update the model in the suprasalt region, we shall perform residual moveout analysis of image gathers. Specifically, a variation of the procedure described in model updating will be followed here. First, convert the image gathers at selected locations along the line from depth (Figure 10.2-13) to time domain (Figure 10.2-14) using the initial velocity-depth model (Figure 10.2-12). The image gathers now are equivalent to moveout-corrected CMP gathers, except that they are in their migrated positions.
  6. The next step is the application of inverse moveout correction to the image gathers from step (e) as shown in Figure 10.2-14. The velocity field used for inverse moveout correction is the rms velocity field shown in Figure 10.2-11, which is consistent with the interval velocity field (Figure 10.2-12) used to create the image gathers themselves.
  7. Now perform conventional velocity analysis using the selected image gathers in time (Figure 10.2-14) and create vertical rms velocity functions at all analysis locations along the line.
  8. Create a new rms velocity field from the vertical functions as shown in Figure 10.2-15. Compare with the initial rms velocity field in Figure 10.2-11 and note the details introduced to the lateral variations in the upper region.
  9. Using the updated rms velocity field from step (h), perform Dix conversion and create an updated velocity-depth model (Figure 10.2-16). Note that the upper region now closely resembles the depositional sequence in the suprasalt region as seen in the time-migrated stacked section (Figure 10.2-2). Also note that the top-salt boundary is now more evident in this model. The lower region still has to be considered as geologically implausable.
  10. Using the updated velocity-depth model (Figure 10.2-16), perform prestack depth migration. Selected image gathers are shown in Figure 10.2-17 and the depth image derived from stacking of the image gathers is shown in Figure 10.2-18. From the flatness of events on image gathers, we may convince ourselves that the model in Figure 10.2-16 can now be considered as final for the suprasalt region. Any further iteration of steps (e) through (j) would only yield insignificant refinement to the model. The model updating described in steps (e) through (j) is only valid if the velocity errors in the initial model are fairly small (model updating).
  11. Interpret the top-salt boundary from the depth image (Figure 10.2-18) and insert it as a layer boundary into the velocity-depth model (Figure 10.2-19). Often, it is adequate to use the depth image from poststack depth migration, rather than prestack depth migration as in this case, to interpret the top-salt boundary.
  12. Assign the salt velocity into the region below the top-salt boundary (Figure 10.2-20). Assume that the velocity within the salt bodies is constant (4450 m/s). This assumption in some cases may not be valid — shale intrusions into halite crystalline rock can alter the velocity within the salt sills, substantially.
  13. Using the new velocity-depth model (Figure 10.2-20), perform prestack depth migration to get a new depth image (Figure 10.2-21).
  14. Interpret the base-salt boundary from the depth image (Figure 10.2-21) and insert it as a layer boundary into the velocity-depth model (Figure 10.2-22).
  15. Finally, introduce the background velocity field into the subsalt region in the form of a vertically varying velocity function (Figure 10.2-23). As for the layered earth model in Figure 10.2-4, the velocity field in the subsalt region is referenced to the water bottom, since the salt bodies have very little influence on the vertical variations in the background velocity field.
  16. Using the final velocity-depth model (Figure 10.2-23), perform prestack depth migration.
  17. Convert the image-gather stack from step (p) from depth to time using the final velocity-depth model, apply poststack spiking deconvolution, band-pass filtering, and AGC scaling, then convert back to depth. Figures 10.2-24 and 10.2-25 show the depth image and selected image gathers, respectively.

Examine the image gathers in Figure 10.2-25 associated with the final velocity-depth model (Figure 10.2-23) and observe that there still are some events that do not meet the flatness criterion. Violation of this criterion, in this case, largely is due to 3-D effects that are not accounted for by 2-D earth modeling and imaging. In the Gulf of Mexico, the salt sills can have complex shapes as shown in Figure 10.2-26. Only by 3-D structural inversion of 3-D data, we can hope to delineate such salt bodies and image the subsalt regions.

10.3 Imaging beneath irregular water bottom in the northwest shelf of australia

The third 2-D case study for structural inversion is from the Northwest Shelf of Australia. This case study deals with the deleterious effect of an irregular water-bottom topography, which is associated with a reef, on the geometry of the underlying reflectors. Shown in Figure 10.3-1 is the CMP-stacked section. The elevation difference between the top and base of the reef is nearly 500 m, and the lateral extent of the reef is about 7.5 km. Note the typical pull-up effect on the reflections below.

The water-bottom reef causes anomalous behavior in the stacking velocity field as seen in Figure 10.3-2. Time migration using a smoothed form of the stacking velocity field produces a distorted image of the subsurface (Figure 10.3-3). Note, for instance, the typical overmigration character below the reef edges.

Earth modeling and imaging in depth

To correct for the effect of the irregular water-bottom topography associated with the reef on the geometry of target reflectors within the substratum, we shall perform earth modeling and imaging in depth. The procedure for earth modeling includes coherency inversion combined with poststack depth migration. Start the analysis by interpreting a set of six time horizons (including the water bottom) within the 0-2.5 s time window from the unmigrated CMP stack (Figure 10.3-1). The time horizons are superimposed on the stacking velocity field in Figure 10.3-2. These time horizons were assumed to be equivalent to zero-offset reflection times and, as such, were used in coherency inversion.

The combination of coherency inversion to estimate the layer velocities and poststack depth migration to delineate the reflector geometries was applied to the six horizons, one layer at a time, starting at the top. Figure 10.3-4 shows the velocity-depth model based on this procedure. Note that the reef has a low-velocity thin cap and a high-velocity interior. The weight of the reef mass may have caused the sagging of the layer boundaries (H2 and H3) associated with the shallow, unconsolidated sediments.

Figure 10.3-5 shows the depth image from prestack depth migration using the velocity-depth model in Figure 10.3-4. Note that the pull-up effect below the reef seen in the time-migrated section (Figure 10.3-3) has been removed. The shelf slope dipping up from left to right is evident from the substratum reflector geometries.

Selected image gathers along the line shown in Figure 10.3-6 contain flat events associated with primaries and events with large moveout associated with multiples. No attempt was made in this case study to attenuate multiples. Conventional CMP stacking attenuates multiples based on velocity discrimination between primaries and multiples. The image gathers in Figure 10.3-6 can be considered similar to CMP gathers that have been moveout corrected using primary velocities. As a result, primaries are flattened and multiples are undercorrected. Therefore, stacking of image gathers to obtain a depth image has a bonus effect of attenuating multiples. However, multiples can make analysis of image gathers for residual moveout a difficult task.

Flatness of the primary events implies the accuracy of the velocity-depth model (Figure 10.3-4) and the accuracy of the depth image (Figure 10.3-5). Compare the depth images derived from prestack depth migration (Figure 10.3-5) and poststack depth migration (Figure 10.3-7) using the same velocity-depth model (Figure 10.3-4). The overmigration effect at the reef edges is quite evident on the poststack depth-migrated section. On the other hand, we remind ourselves the fundamental notion about earth modeling and imaging in depth from introduction to earth modeling in depth — the depth image in Figure 10.3-5 can only be considered a plausable representation of the subsurface geology, and it is not the only representation. Three-dimensional effects increase the ambiguity in the earth model estimated from 2-D seismic data.

10.4 Imaging beneath volcanics in the west of shetlands of the atlantic margin

The fourth 2-D case study for structural inversion is from the West of Shetlands of the Atlantic Margin. This case study deals with earth modeling and imaging beneath volcanics. Shown in Figure 10.4-1 is the CMP-stacked section in two parts. The stratigraphic column comprises shallow depositional sequences that are bounded at the base by a major unconformity (event A). Below the unconformity, is a volcanic layer (event B) with varying thickness. The exploration objective is the potential structural high below 2 s at midpoint location 3440 and the faulted zone below 3.5 s between midpoints 440-1240. Water-bottom and peg-leg multiples associated with the unconformity were largely attenuated using the Radon transform technique (the radon transform) before prestack depth migration. Poststack time migration yields an image of the subsurface that falls short of being acceptable for interpretation within the zone of interest below the unconformity (Figure 10.4-2). To obtain an improved image of this deeper zone, we shall perform earth modeling and imaging in depth. The procedure for earth modeling includes coherency inversion combined with normal-incidence depth conversion applied within the shallow depositional sequences (the overburden), and prestack depth migration and image-gather analysis applied within the zone below the unconformity (the substratum).

Earth modeling and imaging in depth

Start the analysis by interpreting a set of six time horizons, including the water bottom and the unconformity itself, within the overburden from the time-migrated CMP stack (Figure 10.4-2). Then, use these time horizons to obtain the unmigrated time horizons by way of a zero-offset traveltime modeling scheme based on image-ray tracing. The modeled time horizons are assumed to be equivalent to zero-offset reflection times, which are then used in coherency inversion. Except for the water bottom, the following sequence is applied to the five horizons, one layer at a time, starting at the top.

  1. Perform coherency inversion along the line to estimate the layer velocities, and
  2. Perform image-ray depth conversion of the time horizons to delineate the reflector geometries associated with the layer boundaries.
  3. Repeat steps (a) and (b) for all five layers within the overburden down to the unconformity and create a velocity-depth model for the overburden (the upper portion of the model that includes the top five layers in Figure 10.4-3b).
  4. Assign a range of constant velocities to the half-space that corresponds to the substratum (the region below the unconformity A in Figure 10.4-3b) while using the same overburden model (the region down to the unconformity A in Figure 10.4-3b) and perform prestack depth migration to generate image gathers at some interval along the line.
  5. Sort the image gathers into velocity panels as shown in Figures 10.4-4 through 10.4-8.
  6. At a specific midpoint location, examine the moveout of events for flatness. Identify the velocity for which there exists a group of flat events within a depth range. The velocity-depth picks are indicated by the left-arrows in Figures 10.4-4 through 10.4-8. Note that for a velocity-depth pick to be valid, not only the event at that depth must be flat but also all events above it must be flat. If we were to pick an rms velocity function at an analysis location as was demonstrated in Figures 5.4-16 through 5.4-19, for the present case, we would only need to create one velocity panel. As an example, the rms velocity picks are denoted by the right-arrows in Figure 10.4-5.
  7. Combine the picks for the velocity-depth pairs that meet the flatness criterion at each analysis location and create a layer velocity profile (labeled as VL1 in Figure 10.4-3a) and a layer boundary.
  8. Insert the “velocity layer” (labeled as VL1 in Figure 10.4-3b) into the velocity-depth model and define a new half-space region below (Figure 10.4-3b).
  9. Then, repeat steps (d) through (h) for the remaining velocity layers — VL2, VL3, and VL4.
  10. Perform prestack depth migration (Figure 10.4-9) using the final velocity-depth model shown in Figure 10.4-3b.

Compare the results of from prestack depth migration (Figure 10.4-9) and poststack time migration (Figure 10.4-2), and note that an event (labeled as C) within the substratum has been uncovered. This event is at a depth of 3 km at the left edge of the section and dips down to the right. Also note that the complex structure — which may have a 3-D character — below 4 km between midpoints 440-1240 is imaged within the accuracy of 2-D migration.

Selected image gathers along the line (Figure 10.4-10) largely contain flat events associated with primaries, since multiples were attenuated by way of the Radon transform prior to prestack depth migration. Flatness of the primary events demonstrates the accuracy of the velocity-depth model (Figure 10.4-3) and the accuracy of the depth image (Figure 10.4-9). Note that a spatially varying mute has been applied to image gathers (Figure 10.4-10) to attain an optimum depth image.

Figure 10.4-11 shows depth-to-time conversion of the depth image shown in Figure 10.4-9. The conversion was done using the same velocity field that was used for prestack depth migration (Figure 10.4-3). This section should be compared with the time image derived from poststack time migration shown in Figure 10.4-2. Again, note the event at 3 s (labeled as C) in Figure 10.4-11; this event is almost impossible to follow in Figure 10.4-2.

This case study also included a test of the effect of multiples on prestack depth migration. The image gathers shown in Figure 10.4-10 and the corresponding depth image shown in Figure 10.4-9 are based on prestack data with Radon-transform multiple attenuation. Using the same velocity-depth model (Figure 10.4-3), but prestack data without multiple attenuation, prestack depth migration yields the depth image shown in Figure 10.4-12. Note the overwhelming presence of multiples in this depth image; in contrast, the depth image in Figure 10.4-9 is largely free of multiples.

Now examine the image gathers with and without multiple attenuation (Figures 10.4-10 and 10.4-13). Primary events are flat or nearly flat, whereas events associated with multiples have significant moveout. The multiples dominate the image gathers so much that the moveout analyses shown in Figure 10.4-4 through 10.4-8 could not have been conducted using data without multiple attenuation.

Figure 10.4-14 shows depth-to-time conversion of the depth image shown in Figure 10.4-12. Again, as in Figure 10.4-11, the conversion was done using the same velocity field that was used for prestack depth migration (Figure 10.4-3). Note that the event that is distinctively identified in Figure 10.4-11 (labeled as C) is somewhat contaminated by the strong multiples in Figure 10.4-14.

10.5 Imaging beneath shallow gas anomalies in the gulf of thailand

The fifth 2-D case study for structural inversion is from the Gulf of Thailand. This case study deals with the deleterious effect of shallow gas anomalies on a multileveled reservoir zone along growth faults. The velocity-depth model is characterized by near-horizontal layers with structure-independent velocity variations. Shown in Figure 10.5-1 is the CMP-stacked section. The shallow gas anomalies are represented by the high-amplitude piecewise-continuous reflections within the first 500 ms of the section.

The shallow gas anomalies give rise to fluctuations in the stacking velocity field as seen in Figure 10.5-2. For a horizon below the shallow anomaly zone, CMP raypaths at some locations travel through an overburden with laterally varying velocities associated with the shallow gas anomalies. As a result, traveltimes associated with deeper reflections are distorted.

Time migration (Figure 10.5-3) using a smoothed form of the stacking velocity field (Figure 10.5-2) produces a seemingly accurate image of the subsurface. Note, for instance, the fault below CMP 2170 within 1-2 s time window, and the v-shaped structure below CMP 1370. The subsurface appears to be so simple that one may be tempted to convert this section to depth using vertical raypaths and Dix-converted stacking velocities. However, the anomalous behavior of the stacking velocities will cause Dix conversion to yield meaningless interval velocities.

Earth modeling and imaging in depth

To resolve the effect of shallow gas anomalies, we shall perform earth modeling and imaging in depth. The procedure for earth modeling includes coherency inversion to estimate layer velocities combined with image-ray depth conversion to delineate reflector geometries. Start the analysis by interpreting a set of eight time horizons within the 0-2.5 s time window from the time-migrated CMP stack (Figure 10.5-3). These time horizons are used to obtain the unmigrated time horizons by way of ray-theoretical zero-offset modeling. The modeled time horizons are assumed to be equivalent to zero-offset reflection times, which are then used in coherency inversion.

Figure 10.5-4 shows the velocity-depth model based on the procedure outlined above. Note the lateral velocity variations within the layers in the neighborhood of the faulted zones. Figure 10.5-5 shows the depth image from prestack depth migration using the velocity-depth model in Figure 10.5-4. Note the crisp image of the fault below midpoint 2170 and the structural closures along the fault line between 1.5-2.5 km. Another feature of interest is the v-shaped structure centered around midpoint 1370. The boundaries of the closure defined by the two faults that merge at a depth of approximately 2.4 km below midpoint 1370 are not coincident with the changes in the reflector geometries because of faulting. This most likely is caused by the 3-D behavior of the fault, which is not imaged properly by 2-D depth migration.

Selected image gathers along the line shown in Figure 10.5-6 convincingly demonstrate the accuracy of the velocity-depth model (Figure 10.5-4). The flatness of events on the image gathers also provides evidence of the accuracy of the depth image (Figure 10.5-5). Compare the depth images derived from prestack depth migration (Figure 10.5-5) and poststack depth migration (Figure 10.5-7) using the same velocity-depth model (Figure 10.5-4). Although both are better than time migration (Figure 10.5-3), prestack depth migration yields a better image of the faults and the deep, subtle structural closures.

This case study demonstrates that seismic inversion for earth modeling and imaging in depth is not required just for targets below obviously complex structures (subsalt imaging in the North Sea and subsalt imaging in the Gulf of Mexico), but can also be required for resolving low-relief structures in the presence of subtle, shallow velocity anomalies.

10.6 3-D structural inversion applied to seismic data from the southern north sea

The first 3-D structural inversion case study is from the Southern Gas Basin of the North Sea. It involves a complex overburden associated with salt tectonism. The survey area is approximately 90 km2. The surface area is roughly a square with dimensions of 10 km in the crossline direction and 9 km in the inline direction. Survey statistics include more than 9 million traces recorded and nearly 300 000 stacked traces. The nominal fold of coverage is 32, the inline trace spacing is 12.5 m, and the crossline trace spacing is 25 m before trace interpolation and 12.5 m after trace interpolation.

To estimate a 3-D velocity-depth model, the following procedure composed from the list of inversion methods in Table 9-1 is used:

  1. Coherency inversion to estimate layer velocities and 3-D poststack depth migration to delineate reflector geometries within the overburden, and
  2. constant half-space velocity analysis of image gathers from prestack depth migration to estimate the substratum velocity with stacking of image gathers to delineate the reflector geometry of the base-Zechstein (top-Rotliegendes) target horizon.

The procedure outlined above is applied layer-by-layer starting from the surface to resolve layer velocity and reflector geometry of one layer before moving onto the next. This minimizes the influence of any errors made in one layer on the parameter estimates for the next layer.

Figure 10.6-1 shows selected crosslines from the unmigrated 3-D volume of DMO-stacked data. Superimposed on the displays in Figure 10.6-2 are the time horizons that correspond to layer boundaries with significant velocity contrast. These are, starting from the top, base Upper and Lower Tertiary, base Cretaceous chalk, base Upper and Lower Triassic, and base Zechstein (the red horizon). The target zones are the Rotliegendes sands just beneath the base Zechstein, and the underlying Carboniferous Westfalien sequence. Note the traveltime distortions along the base-Zechstein event (the red horizon) — a classic case of the deleterious effect of a complex overburden with strong lateral velocity variations. The complex overburden includes the collapsed structure above the apex and to the left of the salt diapir, and the complex geometry of the top-salt boundary (the brown horizon) which is where the most severe ray bending takes place.

Figure 10.6-3 shows the time surfaces derived from the interpretation of the 3-D volume of unmigrated stack data. These time horizons are considered equivalent to zero-offset time horizons, and as such, are used in coherency inversion to estimate layer velocities. Horizons TH1-TH6, in ascending order, correspond to base Upper and Lower Tertiary, base Cretaceous chalk, base Upper and Lower Triassic, and base Zechstein. The water depth in the survey area is shallow and the Upper Tertiary velocities near the water bottom are very close to the water velocity. Therefore, the water layer was considered as part of the Upper Tertiary layer. In fact, we observe that there is no appreciable velocity contrast at the base Upper Tertiary, and hence, the whole Tertiary section can be considered as one single layer in earth modeling.

Estimation of the overburden model

Based on the volume of stacked data (Figure 10.6-2) and the well data from the area, we make the following characterization of the earth model in depth:

  1. We may consider the earth model in depth in two parts — the overburden above the salt diapir and the substratum that includes the salt mass and the underlying strata.
  2. Time migration may be acceptable within the overburden, but depth migration is imperative within the substratum.
  3. The top-salt boundary (the brown horizon in Figure 10.6-2), which also is the boundary between the overburden and the substratum, is where the most severe ray bending takes place.
  4. The overburden above the Zechstein formation comprises layers with significant velocity contrast — Tertiary, Cretaceous chalk, and Triassic units.
  5. As a result of extensional tectonics, the overburden is intensely faulted. Note the collapsed structures above the apex and to the left of the salt diapir.
  6. Vertical velocity gradients within the overburden layers are significant and vary spatially.
  7. The base-Zechstein reflection times (the red horizon in Figure 10.6-2) are severely distorted by the complex overburden.

Based on these characteristics, an appropriate procedure for estimating an earth model in depth involves layer-by-layer application of 3-D coherency inversion to estimate layer velocities and 3-D poststack depth migration to delineate reflector geometries within the overburden. This is then followed by 3-D prestack depth migration for estimating the velocity field within the substratum.

Model representation by tessellation

Assume that we already have estimated the velocity-depth model for the Tertiary section, and that we want estimate the layer velocity for the Cretaceous chalk layer and delineate the reflector geometry associated with the base of that layer. Hence, we may consider the subsurface velocity-depth model made up of two parts — the known part that includes the Tertiary section, and the unknown part that includes the chalk layer and the layers beneath.

By combining the layer velocities and reflector geometries associated with the known part of the model, create a 3-D velocity-depth model represented in the form of a tessellated volume (introduction to earth modeling in depth). In a tessellated model, the volume associated with each layer is divided into a set of tetrahedra, the size and shape of which depend on the geometry of layer boundaries (Figure 10.6-4). As in the present case study, the estimated velocity and vertical velocity gradient from sonic logs are assigned to each corner of the tetrahedra within the known part of the model. In the unknown part of the model, tetrahedra are populated with constant trial velocities used in the next step (coherency inversion). Tessellated velocity fields can be desirable for efficient ray tracing used in prestack traveltime inversion for velocity estimation. Shown in Figure 10.6-4 is the base-Lower Tertiary layer boundary above which the model is known, and below which is the half-space that includes the layers yet to be determined. Note that the volumes within the known and unknown parts of the model have been subdivided into a set of tetrahedra.

3-D coherency inversion

Perform 3-D coherency inversion to estimate velocity nodes for the layer under consideration. Accompanying the velocity-depth model (Figure 10.6-4), use the time horizon (Figure 10.6-3) associated with the layer boundary picked from the 3-D unmigrated CMP-stacked data volume (Figure 10.6-2) and CMP gathers at analysis locations (Figure 10.6-5). Three-dimensional coherency inversion involves the same steps as for 2-D coherency inversion (models with horizontal layers), except for 3-D ray tracing. For a range of trial constant velocities for the half-space defined by the unknown part of the velocity-depth model, follow the steps below:

  1. Perform 3-D normal-incidence traveltime inversion of the time horizon associated with the layer boundary (Figure 10.6-3) under consideration using the trial constant velocity assigned to the layer above.
  2. Compute 3-D CMP traveltimes using the known overburden 3-D velocity-depth model.
  3. Measure the discrepancy between the modeled and actual CMP traveltimes by way of semblance.
  4. Assign the trial constant velocity that yields the maximum semblance as the velocity of the layer above.

For each layer, velocity estimation using 3-D coherency inversion is done at grid locations 1 km apart in the inline and crossline directions. The panels in Figure 10.6-5 show the CMP gathers at six analysis locations and the semblance curves derived from coherency inversion. The modeled traveltime trajectories associated with the selected velocities that correspond to semblance maxima are plotted on the CMP gathers. For a given CMP location, note from Figure 10.6-6 that the center data window exhibits a flat event coincident with the velocity selection made in Figure 10.6-5.

Vertical velocity gradients, which were available from sonic logs, were accounted for in the ray tracing to model the CMP traveltimes that are required in coherency inversion (models with horizontal layers). Additionally, we want to use available sonic logs to guide coherency inversion and choose an optimum range of constant velocities in computing semblances. As a reminder, coherency inversion also honors ray bending at layer boundaries (models with horizontal layers). Nevertheless, application of coherency inversion should be confined to parts of a subsurface velocity-depth model where no severe raypath distortions may occur (model building).

We circumvent picking of prestack reflection travel-times in coherency inversion by measuring the discrepancy between modeled and actual traveltimes based on semblance. Although picking of time horizons is done from the 3-D DMO-stacked volume of data, note that the common-cell gathers used in the analysis are not processed for dip-moveout (DMO) correction. This is because we model 3-D prestack traveltimes from source locations at the surface to the actual reflection points along each layer boundary and back to receiver locations at the surface associated with traces in a CMP gather.

Results of coherency inversion are used to pick velocities at analysis locations. In a 3-D structural inversion project, quality control of the velocity nodes at analysis locations is imperative for deriving a spatially consistent layer velocity field. We select the layer velocity at a given location based on a combination of the following factors to make sure that the estimated velocities are geologically plausable:

  1. The maximum of the semblance curve derived from coherency inversion as a function of trial constant velocity (Figure 10.6-5),
  2. the flatness of the event in the data windows along the modeled traveltime trajectories (Figure 10.6-6),
  3. the fit of the modeled traveltime trajectory over-layed on top of the actual event on the common-cell gather (Figure 10.6-5),
  4. the lack of anomolous behavior in the raypath associated with the modeled traveltime trajectory, and
  5. the magnitude of the velocity node with respect to the neighboring velocity nodes.

We create the velocity field for the half-space by spatial interpolation between the nodes (Figure 10.6-7). This is followed by the application of some smoothing to the velocity field.

We then create a gridded velocity-depth model that includes the known and unknown part of the model. Gridded velocity fields are desirable for doing 3-D poststack depth migration based on finite-difference schemes. Gridding recognizes layer boundaries and vertical velocity gradients.

3-D poststack depth migration

By using the gridded velocity field, we perform 3-D poststack depth migration down to a depth just below the current layer under consideration (Figure 10.6-8). The algorithm used in the present case for 3-D poststack depth migration is based on a frequency-space explicit scheme and the McClellan transform for designing a 3-D extrapolation operator (3-D poststack migration). It can handle arbitrary vertical and lateral velocity variations and is accurate for dips up to 80 degrees. A 2-D stable explicit operator is converted into a 3-D operator by applying the McClellan transform coefficients. For each frequency component and a given velocity ratio, the explicit operator is stored in a table and fetched as needed at each step of downward continuation. The McClellan transform coefficients are optimized to attain a near-perfect circular symmetry for the impulse response of the 3-D operator. This means that the algorithm, unlike the conventional one-pass 3-D poststack migration algorithms (3-D poststack migration) based on the splitting of the 3-D operator into inline and crossline components, does not cause azimuthal positioning errors.

The image volume derived from 3-D poststack depth migration is used in an interpretation session to incorporate the layer under consideration into the earth model.

  1. We interpret the base of the layer under consideration — the Cretaceous chalk, to delineate the reflector geometry in depth (Figure 10.6-8). Interpretation of the depth horizon is done on crosslines from the 3-D volume of poststack depth-migrated data at every tenth line.
  2. By using the interpretation results from the crosslines, horizon strands are created (Figure 10.6-9a). These strands then are spatially interpolated to create the surface that represents the reflector geometry associated with the base of the layer (Figure 10.6-9b). A way to represent the surface in the computer is by a set of triangles, the size and shape of which vary depending on the complexity of the reflector geometry (Figure 10.6-9c). The reflector geometry is taken into account during ray tracing used in coherency inversion to honor ray bending at layer boundaries, and in 3-D poststack depth migration to account for lateral velocity variations.
  3. Finally, the base Cretaceous chalk surface is added to the model (Figure 10.6-10), and the procedure that includes 3-D coherency inversion and 3-D poststack depth migration is repeated for the next layer within the overburden.

Figure 10.6-11 shows a cross-section of the 3-D velocity field in the crossline direction at each iteration of the procedure described above. Starting from the top, note how the velocity field is updated as the velocity estimate for the next layer is included in the model. The horizons indicated in Figure 10.6-11 correspond to base Tertiary (TH2), Cretaceous chalk (TH3), and Upper Triassic (TH4).

Figure 10.6-12 shows selected crosslines from 3-D poststack depth migrated volumes of data after each iteration. The reflector geometries associated with the layer boundaries included in the overburden model are delineated from such cross-sections — base Lower Tertiary (TH2) from Figure 10.6-12a, base Cretaceous chalk (TH3) from Figure 10.6-12b, base Upper Triassic (TH4) from Figure 10.6-12c, and base Lower Triassic (TH5) from Figure 10.6-12d.

Figure 10.6-13 shows the results of earth modeling using the procedure described above. The left column shows the map view of the layer velocities and the right column shows the perspective view of the base-layer boundaries. Starting from the top, we see layer velocities and reflector geometries for Lower Tertiary, Cretaceous chalk, Upper and Lower Triassic. Note the collapsed zone (A) on the base-Tertiary surface and the auxiliary structural feature (B) that is oblique to the axis of the collapsed zone. The latter becomes more prominent on the surface associated with the base-Cretaceous chalk (C). The ridge of the salt diapiric structure (D) and the collapsed zone (E) to the left of the salt diapir are evident on the surface associated with the base-Upper Triassic. Finally, note that the base-Lower Triassic (equivalent to the top Zechstein) is represented by a multisegmented surface. This is because the Lower Triassic section is missing in the collapsed zone (E) to the left of the salt diapir.

Figure 10.6-10  Tessellated earth model in depth that includes the base-Tertiary (top surface) and base-Cretaceous layer boundary (bottom surface).

Estimation of the substratum model

After the completion of model estimation for the overburden and before moving down to the substratum region, the overburden model must be verified. Model verification usually is done by examining selected image gathers from prestack depth migration along selected crosslines. If the overburden model is correct, events from the overburden on image gathers should exhibit a flat character with negligible moveout. If the residual moveout is significant, then it must be corrected for by residual moveout analysis (model updating), and the overburden model must be updated prior to estimating the substratum model.

To estimate the substratum model, we perform a series of 3-D prestack depth migrations using velocity-depth models that have the same overburden but different constant velocities assigned to the half-space that includes Zechstein and the underlying strata. We examine the image gathers (Figure 10.6-14) to derive velocity nodes for the substratum. As shown in the close-up display in Figure 10.6-15, image gathers associated with different half-space velocities differ only at and below the base Zechstein (shown by the arrow) since they all have been derived using the same overburden velocity-depth model. The event associated with the base Zechstein exhibits a different moveout character from one half-space velocity to another. Based on the flatness criterion (Figure 10.6-15) and the behavior of the geometry of base-Zechstein that can be traced from the image-gather stacks (Figure 10.6-16), we determine velocity nodes at selected locations for the Zechstein formation.

We expect some spatial variation in Zechstein velocities because of the effect of the high-velocity anhydrite-dolomite rafts within the halite unit. By spatial interpolation between the nodes, we create a velocity field for the half-space defined by the substratum region of the model (Figure 10.6-17). We then combine the overburden and the substratum regions of the velocity-depth model. Figure 10.6-18 shows cross-sections of the complete 3-D velocity-depth model along selected crosslines. Finally, we perform 3-D prestack depth migration using this velocity-depth model to obtain crosssections of the depth image along selected crosslines (Figure 10.6-19).

The algorithm for prestack depth migration used in the present case study is based on the Kirchhoff summation (3-D prestack depth migration). Output can be a subsurface image in depth along an arbitrary traverse, inlines and crosslines, or within a specified volume. Also, the output can be image gathers at specified locations for updating layer velocities. The algorithm involves computing traveltimes and summation of amplitudes as implied by the Kirchhoff integral. Traveltimes for all shot and receiver locations at the surface and reflection points in the subsurface are computed by efficient ray tracing. Ray bending at layer boundaries with velocity contrast is honored in computing traveltimes. Migrated data are obtained by summation along maximum-amplitude traveltime trajectories. Prior to summation, data are treated for operator antialiasing (3-D prestack depth migration).

Figure 10.6-20 shows selected image gathers along the crosslines shown in Figure 10.6-19. Except for some events on a few of the image gathers, events exhibit a flat character, which indicates that the earth model (Figure 10.6-18) and the earth image (Figure 10.6-19) are fairly accurate. The base-Zechstein event (denoted by the arrow in Figure 10.6-19) exhibits an improved continuity below the complex zone associated with the salt diapir in the center as compared to the image from 3-D poststack depth migration (Figure 10.6-21).

Note that the structural interpretation of the overburden from the images obtained by 3-D post- and prestack depth migrations should not be different. It is the substratum region where 3-D prestack depth migration would yield an improved image. The differences in amplitude characteristics and the bandwidth between the images from 3-D post- and prestack depth migrations are due to the use of different algorithms. The algorithm for 3-D poststack depth migration is based on an explicit finite-difference scheme, whereas the algorithm for 3-D prestack depth migration is based on the Kirchhoff summation. Although not done in the present case, results of prestack depth migration should be treated with a signal processing sequence commonly applied to stacked data in the time domain, such as poststack deconvolution and frequency filtering.

10.7 3-D structural inversion applied to seismic data from the central north sea

The second 3-D structural inversion case study is from the Central North Sea; it involves a 3-D survey that covers an area of nearly 1800 km2. The surface area is roughly a rectangle with dimensions of approximately 70 km in the northwest-southeast direction, which is coincident with the crossline direction, and 25 km in the inline direction. Survey statistics include 456 million traces recorded, 228 million prestack traces processed, 5.7 million stacked traces, and 11.4 million traces input to 3-D poststack time migration after crossline trace interpolation. The nominal fold of coverage is 40, the inline trace spacing is 12.5 m, and the crossline trace spacing is 25 m before trace interpolation and 12.5 m after trace interpolation.

We want to conduct a layer-by-layer earth model estimation, and use the following combinations of inversion methods:

  1. Coherency inversion to estimate layer velocities and 3-D poststack depth migration to delineate reflector geometries.
  2. Stacking velocity inversion to estimate layer velocities and image-ray depth conversion of time horizons interpreted from the time-migrated volume of data to delineate reflector geometries.

We then want to compare the depth structure maps for the target horizon base Zechstein (top Rotliegendes) obtained from the two procedures outlined above.

Figure 10.7-1 shows selected inlines from the unmigrated 3-D volume of DMO-stacked data. Superimposed on these displays are the time horizons that correspond to layer boundaries with significant velocity contrast. These are, starting from the top, base Upper and Lower Tertiary, base Cretaceous chalk, base Upper and Lower Triassic, and base Zechstein (the red horizon). The target zone is the Rotliegendes sands just beneath base Zechstein, and the underlying Carboniferous Westfalien sequence.

The Zechstein formation comprises two units of anhydrite-dolomite embedded in the halite unit. One of the anhydrite units is very close to the top boundary of Zechstein (the green horizon) and is concordant with it. The deeper, second anhydrite unit has a very complex geometry as a result of the intensive salt tectonics in the area. The abundance of diffractions within the Zechstein formation is associated with this anhydrite-dolomite unit. One of the objectives of this case study is to investigate whether the complex geometry of the second anhyrite-dolomite unit has any subtle, but deleterious effect on the geometry of the target horizon — base Zechstein (top Rotliegendes).

Figure 10.7-2 shows selected inlines as in Figure 10.7-1 from the 3-D poststack time-migrated volume of data. The velocity field used for 3-D poststack time migration is a smoothed version of the 3-D stacking velocity field. Note the complex geometry of the second anhydrite-dolomite unit within Zechstein. The faults along the base Zechstein act as seals on potential reservoir margins. Accurate positioning of these faults, and determination of their displacements in depth, are therefore, central to exploration and development objectives in the area.

Table 10-1. Lateral variations in horizon-consistent stacking velocities shown in Figure 10.7-4.
Horizon Layer Boundary Velocity Range (m/s)
2 Base Upper Tertiary 1450 – 1950
3 Base Lower Tertiary 1650 – 2350
4 Base Cretaceous 2150 – 2500
5 Base Upper Triassic 2550 – 2800
6 Base Lower Triassic 2850 – 3100
7 Base Zechstein 3200 – 3550

Figure 10.7-3 shows the time surfaces derived from the interpretation of the 3-D volumes of unmigrated and time-migrated stack data. Horizons 2-7, in ascending order, correspond to base Upper and Lower Tertiary, base Cretaceous chalk, base Upper and Lower Triassic, and base Zechstein. The time horizons interpreted from the unmigrated DMO-stacked volume are considered equivalent to zero-offset time horizons. They are used in Dix conversion of stacking velocities, stacking velocity inversion, and coherency inversion to estimate layer velocities. The time horizons interpreted from the time-migrated volume are used in map migration by way of image-ray depth conversion to obtain the corresponding depth horizons.

Figure 10.7-4 shows the horizon-consistent stacking velocities that were extracted from the 3-D stacking velocity field along time horizons interpreted from the unmigrated volume of stacked data shown in Figure 10.7-3. Table 10-1 shows the range of lateral velocity variations within each layer.

3-D coherency inversion combined with 3-D poststack depth migration

Examination of the earth image volume in time (Figure 10.7-2) and the velocity field (Figure 10.7-4) leads to the following characterization of the earth model in depth:

  1. The complex geometry of the second anhydrite-dolomite unit within Zechstein may influence the geometry of the target level — base Zechstein, and therefore, it may have to be included in the earth model.
  2. Lateral velocity variations within the overburden layers are, in some parts of the survey area, require imaging in depth.
  3. Velocity contrast at layer boundaries causes significant ray bending.
  4. Vertical velocity gradients within the overburden layers are significant and vary spatially. For instance, the gradient within the Cretaceous chalk layer can be as much as 1.5 m/s/m.

Based on these characteristics, an appropriate procedure for estimating an earth model in depth involves a layer-by-layer application of 3-D coherency inversion to estimate layer velocities and 3-D poststack depth migration to delineate reflector geometries. Assume that we already have estimated the velocity-depth model for the first n − 1 layers, and that we want to estimate the layer velocity for the nth layer and delineate the reflector geometry associated with the base of that layer. Hence, we may consider the subsurface velocity-depth model made up of two parts — the known part that includes the first n − 1 layers, and the unknown part that includes the nth layer and the layers beneath. The following sequence is conducted for each of the layers in the unknown part of the model.

  1. Create a 3-D velocity-depth model represented in the form of a tessellated volume.
  2. Perform 3-D coherency inversion to estimate velocity nodes for the layer under consideration. Coherency inversion requires time horizons picked from a volume of 3-D unmigrated CMP-stacked data (Figure 10.7-3) and CMP gathers at analysis locations. For each layer, estimate layer velocities using 3-D coherency inversion at grid locations 1 km apart in the inline and crossline directions. Figure 10.7-5 shows results of coherency inversion at one specific analysis location. The bottom frames are copies of the CMP gather at the analysis location, while the top frames are copies of the semblance curves derived from coherency inversion (in blue) and stacking velocity inversion (in red). The 3-D coherency inversion semblance curve is a measure of the discrepancy between the modeled and actual CMP traveltimes associated with the reflection event from the layer boundary under consideration. The stacking velocity inversion semblance curve is a measure of the discrepancy between the modeled and actual stacking velocities associated with the reflection event from the layer boundary under consideration. The vertical bars in red are aligned with the selected trial velocities that correspond to the center data windows in the panels shown in Figure 10.7-6. The modeled traveltime trajectories associated with the three selected velocities are plotted on the CMP gather displays in Figure 10.7-5. Note from Figure 10.7-6a, the center data window exhibits residual moveout which suggests erroneously too low velocity, corresponding to the selection made in Figure 10.7-5a. Similarly, the center data window in Figure 10.7-6c exhibits residual moveout which suggests erroneously too high velocity, corresponding to the selection made in Figure 10.7-5c. The center data window in Figure 10.7-6b, however, exhibits a flat-event character, which suggests that the selected velocity in Figure 10.7-5b is optimum. This velocity coincides with the peak value of the 3-D coherency inversion semblance curve and the minimum of the stacking velocity inversion semblance curve.
  3. Use the results of 3-D coherency inversion — the semblance spectrum and the modeled CMP traveltimes (Figure 10.7-5), and the residual moveout panels (Figure 10.7-6) to pick velocities at analysis locations.
  4. Assign the velocity nodes to the half space defined by the unknown part of the model which includes the layer to be delineated. Then, by spatial interpolation between the nodes, a velocity field is created for the half-space.
  5. Next, create a gridded velocity-depth model that includes the known and unknown part of the model.
  6. By using the gridded velocity field, perform 3-D poststack depth migration down to a depth just below the current layer under consideration.
  7. The next step involves interpretation of the base of the layer under consideration to delineate the reflector geometry in depth. Interpretation of the depth horizon is done on cross-sections from the 3-D volume of poststack depth-migrated volume of data at an interval that accounts for the complexity of the reflector geometry. By using the interpretation results from the cross-sections, we perform a surface fit to obtain the 3-D representation of the reflector geometry associated with the base of the layer.
Figure 10.7-4  Horizon-consistent stacking velocities that correspond to unmigrated time horizons shown in Figure 10.7-3. Horizons 2-7, in the ascending order, correspond to base Upper and Lower Tertiary, base Cretaceous chalk, base Upper and Lower Triassic, and base Zechstein.

Figure 10.7-7 shows selected inlines as in Figure 10.7-1 from the volume of 3-D poststack depth-migrated data following the layer-by-layer application of the procedure described above. Compare with the results of 3-D poststack time migration (Figure 10.7-2) and note the subtle differences in the geometry of the base-Zechstein boundary.

The second anhydrite-dolomite unit of Zechstein was included in the earth modeling as a layer of constant thickness (100 m). Note from Figure 10.7-8 the complex geometry of this unit. Although it has a complex geometry and a significant velocity contrast, this unit does not have much effect on the geometry of the base-Zechstein boundary. This can be explained by noting that the thickness of the unit is sufficiently small making the traveltime within the unit insignificant despite its very fast velocity (5900 m/s).

Figure 10.7-9 shows the results of earth modeling using the procedure described above. Starting from the top, we see layer velocities and reflector geometries for Upper and Lower Tertiary, base Cretaceous chalk, base Upper and Lower Triassic, and base Zechstein. To obtain an accurate depth structure map for the base-Zechstein boundary that is of interest to the explorationist, the layer velocities and reflector geometries of the layers above are estimated one layer at a time starting from the surface. As a result, the accumulation of errors in the earth model parameters at the zone of interest (base Zechstein) is minimized.

3-D stacking velocity inversion combined with 3-D image-ray depth conversion

We now consider a traditional approach for obtaining structure maps in depth — stacking velocity inversion to estimate layer velocities and image-ray depth conversion of time horizons interpreted from the time-migrated volume of data to delineate reflector geometries. Image-ray depth conversion is the usual implementation of map migration (model building).

Application of stacking velocity inversion to 3-D data involves the same steps as for 2-D stacking velocity inversion (model with low-relief structure) except for 3-D ray tracing. As a reminder, stacking velocity inversion can be considered accurate for velocity-depth models with smoothly varying reflector geometries and lateral velocity variations much greater than a cable length.

Figures 10.7-10 and 10.7-11 show the results of this alternative procedure. The input to stacking velocity inversion is the set of horizon-consistent stacking velocities (Figure 10.7-4). The input to image-ray depth conversion is the set of time horizons picked from the time-migrated volume of stacked data (Figure 10.7-3). As in coherency inversion combined with poststack depth migration, this procedure was conducted layer-by-layer starting from the surface. Actually, it is conducted twice — without (Figure 10.7-10) and with (Figure 10.7-11) vertical velocity gradients from well data included in the analysis. In the latter case, vertical gradients are accounted for in ray tracing to compute the CMP traveltimes for coherency inversion (models with horizontal layers).

The depth structure map for the target horizon — base Zechstein, has the gross features (Figure 10.7-10) similar to those of the map obtained from coherency inversion combined with poststack depth migration (Figure 10.7-9). Nevertheless, it is the subtle differences that can influence the interpretation.

A way to judge whether imaging can be done by time migration or needs to be done by depth migration is by examining behavior of image rays. Specifically, lateral displacement between the point of departure of an image ray at the target horizon and the point of emergence of the image ray at the surface is an indication of the degree of lateral velocity variations (Figure 10.7-12a). In the case of 3-D image rays, aside from lateral displacement, the azimuthal direction of the displacement is an additional attribute that needs to be examined (Figure 10.7-12b). Note that the largest displacements occur along the major fault zone within the survey area.

We shall end the discussion on this case study by comparing the results of the three inversion procedures that we have applied to the data. Shown in Figure 10.7-13 are the depth structure map of the target horizon — base Zechstein derived from the following procedures:

  1. coherency inversion, in which vertical velocity gradients were accounted for, combined with poststack depth migration,
  2. stacking velocity inversion, in which vertical velocity gradients were accounted for, combined with image-ray depth conversion,
  3. stacking velocity inversion, in which vertical velocity gradients were ignored, combined with image-ray depth conversion.

Shown in Figure 10.7-14 are the difference maps, where the map of Figure 10.7-13a was taken as the reference map, assuming that it is the most accurate between the three maps shown in Figure 10.7-13. These difference maps are based on the subtraction of the maps of Figures 10.7-13b and 10.7-13c from the reference map. Note that there are differences in the results obtained from the different procedures outlined above. While procedures (a) and (b) appear to produce results that are comparable, procedure (c) produces results with significant discrepancies along the major fault zone.

10.8 3-D structural inversion applied to seismic data from offshore indonesia

The third 3-D structural inversion case study is from offshore Indonesia; it involves a 3-D survey that covers an area of nearly 150 km2. The surface area is roughly a rectangle with dimensions of approximately 15 km in the inline direction and 10 km in the crossline direction. Survey statistics include more than 15 million traces recorded and 480 000 stacked traces. The nominal fold of coverage is 32, the inline trace spacing is 12.5 m, and the crossline trace spacing is 25 m. Figure 7.4-22 shows three inline sections from the volume of 3-D prestack time-migrated data from this 3-D survey.

Model building

To build the model, we shall apply the structure-independent inversion strategy discussed in model building. The steps for model building are outlined below.

  1. As for the time-to-depth conversion sequence described in model building, the Dix equation (9-1) is used to derive a 3-D interval velocity field from a 3-D rms velocity field. The latter is preferably estimated from prestack time-migrated data. Recall the sequence for 3-D prestack time migration described in 3-D prestack time migration. Following the application of NMO and 3-D DMO correction, each common-offset volume of data is migrated using an initial velocity field derived from the DMO-corrected data. The migrated data are then used to perform velocity analysis, once again, over the survey area at specified grid locations (Figure 7.4-19). The picked rms velocity functions are then used to create a 3-D rms velocity field. Figure 10.8-1 shows selected time slices from the 3-D rms velocity volume.
  2. Perform Dix conversion to derive the 3-D interval velocity field (Figure 10.8-2).
  3. By using this structure-independent 3-D interval velocity field, create an image volume from 3-D poststack depth migration. Figure 10.8-3 shows three inline sections from the image volume. The same inline sections with reverse polarity for better identification of some of the fault planes are shown in Figure 10.8-4. Selected depth slices from the image volume shown in Figure 10.8-5 exhibit a complex fault pattern.
  4. Interpret a set of depth horizons from the image volume of 3-D poststack depth migration created in step (c). The resulting depth structure maps for six horizons that correspond to layer boundaries with significant velocity contrast are shown in Figure 10.8-6.
  5. Intersect the 3-D interval velocity field from step (b) with the depth horizons from step (d) and extract a horizon-consistent interval velocity map for each layer (Figure 10.8-6).
  6. Combine the interval velocity maps from step (e) with the depth horizon maps from step (d) to construct an initial 3-D velocity-depth model.

Model updating

The steps for model updating by reflection tomography are outlined below.

  1. Perform 3-D prestack depth migration using the initial 3-D velocity-depth model (Figure 10.8-6) and generate the image gathers along selected inline traverses. Note from the selected image gathers shown in Figure 10.8-7 that majority of the events are nearly flat. Nevertheless, residual moveout along some of the events calls for a modification to the initial velocity-depth model. Figure 10.8-8 shows three inline sections derived from stacking the image gathers. The same inline sections with reverse polarity for better identification of some of the fault planes are shown in Figure 10.8-9. Compare these image sections from 3-D prestack depth migration with those derived from 3-D poststack depth migration shown in Figures 10.8-3 and 10.8-4, and note that prestack depth migration has improved the image quality within the zone with intensive faulting.
  2. Assume that residual moveout on image gathers can be described by a parabolic moveout equation when analyzed in time, and compute the horizon-consistent residual moveout semblance spectra for events on image gathers that correspond to the layer boundaries included in the velocity-depth model. Shown in Figure 10.8-10 are the semblance spectra for three horizons H2, H3, and H4 as in Figure 10.8-6 along three inline traverses.
  3. For each horizon, pick the residual moveout profiles by tracking the maxima of the semblance spectra for each of the inline traverses and create a residual moveout map as shown in Figure 10.8-11 associated with a reference offset of the image gathers — usually the maximum offset.
  4. Use the residual moveout maps for the reference offset from step (c) and compute the residual move-out for all offsets in the image gathers based on the parabolic moveout assumption. Build the traveltime error vector Δt of equation (9-7) using the residual moveout times.
  5. By using the interval velocity and depth-horizon maps associated with the initial velocity-depth model (Figure 10.8-6), construct the coefficient matrix L in equation (J-88) as described in Section J.6.
  6. Estimate the change in parameters vector Δp by way of the GLI solution given by equation (9-7). The tomography matrix equation (9-7) can be set up such that the residual moveout times given by the traveltime error vector Δt are used to perturb interval velocities and reflector geometries for not necessarily all of the horizons, but any combination of them. In fact, if judged to be appropriate, only one set of the parameters — either layer velocities or reflector geometries may be perturbed.
  7. Update the parameter vector p + Δp. Figure 10.8-12 shows the updated interval velocity and depth-horizon maps. Compare with the interval velocity and depth horizon maps associated with the initial model (Figure 10.8-6) and note that tomographic update was actually applied to the interval velocities while the reflector geometries represented by the depth horizon maps were not altered.
  8. Now combine the new set of interval velocity and depth horizon maps after the tomographic update (Figure 10.8-12) to construct a new velocity-depth model.
  9. Perform 3-D poststack depth migration to create the image volume using the updated velocity-depth model. Shown in Figure 10.8-13 are three inline sections from the image volume. Again, the same inline sections with reverse polarity for better identification of some of the fault planes are shown in Figure 10.8-14. Compare these image sections with those derived from 3-D poststack depth migration before the tomographic update shown in Figures 10.8-3 and 10.8-4, and note that the tomographic update of the model has improved the image quality within the zone with intensive faulting. Selected depth slices from the image volume are shown in Figure 10.8-15.
  10. Re-interpret the depth horizons using the image volume from step (i) and examine discrepancies with the previous interpretation. If required, repeat steps (a) through (i) until residual moveouts for events associated with the layer boundaries included in the model have been reduced to negligible magnitudes.

Imaging in depth

Following the tomographic updating described above, perform 3-D prestack depth migration using the final velocity-depth model (Figure 10.8-12) to generate the image sections along the selected inline traverses. Shown in Figure 10.8-16 are selected image gathers along the three inline traverses. Stack of the image gathers yields the final image from 3-D prestack depth migration. Shown in Figure 10.8-17 are three inline sections from the image volume and in Figure 10.8-18 are the same inline sections with reverse polarity for better identification of some of the fault planes. Compare these image sections with those derived from 3-D prestack depth migration before the tomographic update shown in Figures 10.8-8 and 10.8-9, and note that the tomographic update of the model has improved the image quality within the zone with intensive faulting.

Once the image quality from the selected inlines is judged to be acceptable, perform 3-D prestack depth migration to produce the entire image volume. Shown in Figure 10.8-19 are selected inline sections from the image volume derived from 3-D prestack depth migration using the final velocity-depth model. While the selected inline sections from 3-D prestack depth migration as in Figure 10.8-18 were created using the Kirchhoff summation algorithm, the entire image volume represented by the inline sections in Figure 10.8-19 was created using the common-offset 3-D phase-shift-plus-correction algorithm (3-D prestack depth migration). This choice for creating the image volume was primarily for faster turnaround. Selected crossline sections from the image volume in depth are shown in Figure 10.8-20. Note from the inline and crossline sections in Figures 10.8-19 and 10.8-20 the intensive faulting and folding caused by the wrench tectonism.

Volume-based interpretation

Interpretation of 3-D seismic data is based on one of the following two strategies:

  1. Image in time, interpret in time to create a set of time structure maps, and then perform time-to-depth conversion to create the corresponding set of depth structure maps, or
  2. Image in depth and interpret directly in depth to create a set of depth structure maps.

Whatever the strategy, the resulting depth structure maps are calibrated to well data.

We demonstrated the first strategy for interpretation in interpretation of 3-D seismic data. We now demonstrate the second strategy for interpretation. Figure 10.8-21 shows selected 3-D views of the image volume from 3-D prestack depth migration of the data as in Figure 10.8-19 in the inline direction, and Figure 10.8-22 shows selected 3-D views of the image volume in the crossline direction as in Figure 10.8-20. Begin the structural interpretation session by 3-D visualization of the image volume in depth to understand and characterize the subsurface structural model. Combine the visualization of the image volume in the inline and crossline directions with the the scanning of the depth slices (Figure 10.8-23). Note that the structural model is based on wrench tectonism that has caused the intensive faulting and folding.

By using the seed detection technique (interpretation of 3-D seismic data), we capture a set of surface patches for each of the depth horizons to be interpreted. Where seed detection fails, as in the complex fault zones, the surface paths are complemented with the horizon strands along the selected inlines and crosslines derived from line-based interpretation. The combination of seed detection and line-based interpretation produce a set of control points for each depth horizon as shown in Figures 10.8-24a through 10.8-29a. The control points are then used to create a complete surface for each horizon by grid-fitting (Section J.5) as shown in Figures 10.8-24b through 10.8-29b. The sixth horizon, which is the deepest, shown in Figure 10.8-29b is deeper in some areas than the maximum depth (4 km) associated with the image volume from 3-D prestack depth migration (Figures 10.8-21 and 10.8-22). All six depth horizons are shown in the 3-D view of Figure 10.8-30.

Stratigraphic interpretation uses an amplitude manipulation technique based on removal of opacity (interpretation of 3-D seismic data). Application of this technique to the depositional unit bounded by surface DH3 at the top and DH4 at its base (Figure 10.8-31a) leads to the discovery of an ancient channel (Figure 10.8-31b). This buried channel within the depositional unit of Figure 10.8-31a resembles the recent channel system at the water bottom which is dramatically manifested by the depth slices shown in Figure 10.8-32.

A close-up view reveals the disruption of the buried channel by the faulting associated with the wrench tectonism in the area (Figure 10.8-33a). By using the subvolume detection (interpretation of 3-D seismic data), the channel can be isolated as shown in Figure 10.8-33b. The color-coded view of the channel shows that elevation changes exist along the channel caused by the faulting that has occurred after the formation of the channel itself (Figure 10.8-34a). After its subvolume detection, the channel can be extracted from the volume associated with the depositional unit (Figure 10.8-31a) and viewed in combination with the surface that corresponds to the base of the unit (Figure 10.8-34b). The objective is to gain a complete understanding of the stratigraphic prospect represented by the channel in relation to the subsurface structural model.

Finally, the results of the structural interpretation (Figures 10.8-24b through 10.8-29b) are combined to build a solid earth model. Figure 10.8-35b shows the layer boundaries represented by the depth horizons that were derived from the structural interpretation. Each layer is then represented by a solid volume bounded by these surfaces as shown in Figure 10.8-35. Explosion of the entire solid model facilitates viewing of the individual depositional units represented by the solid segments (Figure 10.8-36). The solid segment associated with the zone of interest can be subsequently populated by the seismically derived petrophysical properties, which are constrained by the well data, to derive a reservoir model.

10.9 3-D structural inversion applied to seismic data from the northeast china

The fourth 3-D structural inversion case study is from the onshore northeast China; it involves a 3-D survey that covers an area of nearly 150 km2. The surface area is roughly a rectangle with dimensions of approximately 13 km in the inline direction and 12 km in the crossline direction. Survey statistics include more than 1.8 million traces recorded and nearly 120 000 stacked traces. The average fold of coverage is 30, the inline trace spacing is 25 m, and the crossline trace spacing is 50 m.

The fold of coverage is fairly uniform over the survey area as seen from the map in Figure 10.9-1. There are 19 swaths of data recorded using four cables that are 100-m apart, each with 60 receiver groups at 50-m interval. The shots were placed in the direction perpendicular to the swaths. A sketch of the receiver cable layout is also shown in Figure 10.9-1.

Shown in Figure 10.9-2 are selected shot records from each of the 19 swaths with a display gain applied. Each shot record contains 240 traces that correspond to the four 60-channel receiver cables. Examine the shot records and note that the quality of the recorded data is remarkably good, except for some coherent, dispersive, low-frequency and high-amplitude noise at near offsets that may be associated with the ground-roll energy.

With this case study, we shall develop a comprehensive, time-and-depth workflow for processing, inversion and interpretation of 3-D seismic data. This workflow is applicable to many of the oil and gas provinces with low-relief structures and complex structures that involve folding and faulting caused by extensional or compressional tectonism, and moderate-to-strong lateral velocity variations. The workflow, however, is not suitable for imaging beneath complex overburden structures associated with salt tectonism or overthrust tectonism, and strong-to-severe lateral velocity variations (2-D poststack depth migration).

We shall execute a unified seismic workflow in five phases as illustrated in Figure 10.9-3 and outlined below:

  1. 3-D DMO processing with deliverables that include a set of 3-D DMO gathers, a 3-D DMO stack volume, a 3-D DMO velocity field and an image volume derived from 3-D poststack time migration.
  2. 3-D prestack time migration with deliverables that include a set of common-reflection-point (CRP) gathers, an image volume from 3-D prestack time migration, a 3-D rms velocity field and a 3-D zero-offset wavefield modeled from the image volume. The latter may be treated as the equivalent of an unmigrated stack volume associated with the 3-D prestack time-migrated data.
  3. Stratigraphic inversion with deliverables that include the AVO attribute volumes derived from prestack amplitude inversion and an acoustic impedance volume derived from poststack amplitude inversion. In this phase, we also derive a 3-D interval velocity field by Dix conversion of the rms velocity field from phase 3 and use it to create an image volume from 3-D poststack depth migration.
  4. Structural inversion that involves the estimation of an initial 3-D velocity-depth model and model updating with deliverables that include a final 3-D velocity-depth model and an image volume derived from 3-D prestack depth migration.
  5. Structural and stratigraphic interpretation, and calibration to well data with the ultimate objective of building a reservoir model that is fundamentally a petrophysical representation of an earth model in depth.

In the present case study, we shall execute the entire time-with-depth workflow outlined above, except for stratigraphic inversion in phase 3 and reservoir modeling in phase 5.

3-D DMO processing

Prestack signal processing of the data included geometric spreading correction using a t2-scaling function (Figure 10.9-4), spiking deconvolution and time-variant spectral whitening (Figure 10.9-5). This parsimonious sequence sufficed to attain a desired flat and broad spectrum within the passband (Figure 10.9-6). Note that the coherent linear noise that dominates the raw shot record (Figure 10.9-4a) does not stand out in the shot gather after spectral balancing (Figure 10.9-5b).

Following the signal processing, the data are ready for DMO processing:

  1. Sort the common-shot gathers to common-cell gathers and perform stacking velocity analysis over a sparse grid of 2 × 2 km.
  2. Apply NMO and 3-D DMO correction to remove source-receiver azimuth and dip effects from the stacking velocities.
  3. Apply inverse NMO correction and repeat the stacking velocity analysis over a finer grid of 0.5 × 0.5 km. Figure 10.9-7 shows a subset of four velocity analysis panels over the survey area. Each panel comprises the common-cell gather at the analysis location and the velocity spectrum computed from it. Generally, the velocity picking from the DMO-corrected data was consistently reliable over the survey area.
  4. Create a 3-D DMO velocity field from the vertical functions picked in step (c).
  5. Apply NMO correction to the 3-D DMO-corrected data from step (c) using the velocity field from step (d) and stack the data. Figures 10.9-8 and 10.9-9 show selected inline and crossline sections, respectively, from the volume of the 3-D DMO stack following poststack deconvolution and band-pass filter.
  6. Perform 3-D poststack time migration using a laterally smoothed form of the 3-D DMO velocity field from step (d). Figures 10.9-10 and 10.9-11 show selected inline and crossline sections, respectively, from the image volume derived from 3-D poststack time migration. To circumvent the adverse effect of spatial aliasing on migration, especially in the crossline direction, trace interpolation before migration was considered. However, this process deteriorated the stacked data in zones of intensive faulting associated with extensional tectonism. Therefore, in lieu of trace interpolation, the unmigrated data in step (e) were filtered down to a passband of 6-36 Hz; hence the difference in the frequency content between the unmigrated data in Figure 10.9-8 and the migrated data in Figure 10.9-10.

The deliverables from phase 1 — 3-D DMO processing, include a set of 3-D DMO-corrected gathers, a volume of 3-D DMO-stacked data, a volume of 3-D DMO velocity field, and an image volume derived from 3-D poststack time migration.

3-D prestack time migration

While the interpreter may have a sneak preview of the subsurface image using the volume of 3-D poststack time-migrated data from phase 1, we move on to phase 2 for 3-D prestack time migration.

  1. Sort the 3-D DMO-corrected and moveout-corrected data from step (e) of phase 1 to common-offset volumes.
  2. Assume that each of the common-offset volumes of data is a replica of a 3-D zero-offset wavefield and perform 3-D zero-offset time migration using a regionally averaged, but vertically varying velocity function derived from the 3-D DMO velocity field from step (d) of phase 1.
  3. Sort the migrated common-offset volumes of data back to common-cell gathers, apply inverse moveout correction using the 3-D DMO velocity field from step (d) of phase 1, and repeat the velocity analysis over a grid of 0.5 × 0.5 km.
  4. Create a 3-D rms velocity field associated with the migrated data from the vertical functions picked in step (c).
  5. Apply NMO correction to the 3-D common-offset-migrated data from step (c) using the rms velocity field from step (d) and stack the data. Figures 10.9-12 and 10.9-13 show selected inline and crossline sections, respectively, from the image volume derived from the stack of the 3-D common-offset migrations. Again, to circumvent the adverse effect of spatial aliasing, the common-offset volumes of data were filtered down to a passband of 6-36 Hz prior to migration.
  6. Model a 3-D zero-offset wavefield from the image volume derived from the 3-D common-offset migration and stack in step (e) using the same vertically varying velocity function as in step (b). Figures 10.9-14 and 10.9-15 show selected inline and crossline sections, respectively, from the modeled 3-D zero-offset wavefield volume.
  7. The final step in phase 2 involves remigration of the stacked data. Specifically, perform 3-D poststack time migration for which the input data volume is the 3-D zero-offset modeled wavefield from step (f) and the migration velocity field is the 3-D rms velocity field from step (d). Figure 10.9-16 shows selected inline cross-sections from the 3-D rms velocity volume, and Figures 10.9-17 and 10.9-18 show selected inline and crossline sections, respectively, from the image volume derived from 3-D prestack time migration based on the above sequence. The migration artifacts on the left-hand side of some of the sections are caused by the missing data zone within the survey area (bottom left corner of the fold map in Figure 10.9-1).

The deliverables from phase 2 — 3-D prestack time migration, include a set of CRP gathers, a volume of the 3-D zero-offset wavefield which may be considered as equivalent to an unmigrated stack volume, a volume of the 3-D DMO rms velocity field, and an image volume derived from 3-D prestack time migration.

From RMS to interval velocities

As we note from the list of deliverables from phase 2, we are not implementing phase 2 — 3-D prestack time migration, just for imaging the subsurface. We shall make use of the 3-D rms velocity field and the 3-D zero-offset wavefield from phase 2 to move from time to depth domain in the analysis. Although it will not be considered in the present case study, phase 3 of the generalized workflow also includes amplitude inversion that uses the CRP gathers from 3-D prestack time migration and the image volume itself.

  1. Apply prestack amplitude inversion (analysis of amplitude variation with offset) to the CRP gathers and derive the AVO attribute volumes. These include the intercept and gradient volumes, P-wave and S-wave reflectivity and impedance volumes, pseudo-Poisson and fluid-factor volumes (analysis of amplitude variation with offset). To satisfy the underlying assumptions for AVO analysis, and specifically the Zoeppritz equations that describe the partitioning of energy associated with an incident compressional plane wave at a horizontal layer boundary into its reflected and refracted compressional-and shear-wave components, we want to use the CRP gathers from 3-D prestack time migration that contain events in their migrated positions in prestack amplitude inversion rather than the DMO gathers that contain events in their unmigrated positions.
  2. Apply poststack amplitude inversion (acoustic impedance estimation) to the AVO intercept volume to derive the acoustic impedance volume. It is assumed that the AVO intercept attribute is a close representation of the reflection amplitudes at vertical incidence and zero offset. Alternatively, poststack amplitude inversion can be applied to the image volume from 3-D prestack time migration derived in step (g) of phase 2.
  3. Apply Dix conversion to the 3-D rms velocity field from step (d) of phase 2 to derive a 3-D interval velocity field. Figure 10.9-19 shows selected inline cross-sections from the 3-D interval velocity volume. To satisfy the underlying assumptions for Dix conversion (Section J.4), we want to use the rms velocity field derived in conjunction with 3-D prestack time migration of the data in phase 2 that is associated with events in their migrated positions rather than the DMO velocity field in conjunction with 3-D DMO correction of the data in phase 1 that is associated with events in their unmigrated positions.
  4. Perform 3-D poststack depth migration of the 3-D zero-offset wavefield from step (f) of phase 2 using the 3-D interval velocity field from step (c) of the present phase of the workflow.

The deliverables from phase 3 include a volume of 3-D interval velocity field and an image volume derived from 3-D poststack depth migration.

Structural inversion

We want to construct a structurally consistent 3-D velocity-depth model using the deliverables from phase 3 and update it to obtain a final 3-D velocity-depth model. We then want to use this earth model in depth to create an earth image in depth from 3-D prestack depth migration.

  1. Interpret a set of depth horizons from the image volume derived from 3-D poststack depth migration (Figure 10.9-20a). These depth horizons correspond to layer boundaries with significant velocity contrast.
  2. To preserve the vertical and lateral velocity variations in the 3-D interval velocity field derived from Dix conversion (Figure 10.9-19), create a set of phantom depth horizons by subdividing each of the layers bounded by the depth horizons (Figure 10.9-20a) interpreted from the 3-D poststack depth-migrated volume of data into a set of 10 sublayers. The phantom horizons can be conformed either to the top-boundary or the base-boundary, or, as in the present case, to both the top- and base-boundary of the layer under consideration. Figure 10.9-20b shows the principal depth horizons as in Figure 10.9-20a interleaved with the phantom horizons in green.
  3. Intersect the 3-D interval velocity volume (Figure 10.9-19) from step (c) of phase 3 with the principal and phantom depth horizons from step (b) of the present phase and extract horizon-consistent interval velocity surfaces. Then, combine these layer velocities with the reflector geometries represented by the depth horizons to create an initial structurally consistent 3-D velocity-depth model. Figure 10.9-21 shows selected inline cross-sections from the initial 3-D velocity-depth model.
  4. Perform 3-D prestack depth migration and generate a set of image gathers along every tenth inline, and compute the residual moveout semblance spectra for model updating (model updating). Figure 10.9-22 shows selected panels of residual-moveout analysis. Each panel comprises the image gather at the analysis location and the residual moveout semblance spectrum computed from it. While many of the events in the image gathers show very small or no residual moveout (panel I-51), note that some events exhibit significant residual moveout errors (panel I-41). Perform residual moveout corrections and update the initial velocity-depth model (model updating).
  5. Repeat the residual-moveout analysis a few times until the residual-moveout errors have been reduced significantly. Figures 10.9-23 and 10.9-24 show inline and crossline cross-sections, respectively, from the final 3-D velocity-depth model. The residual-moveout analysis from the last iteration of model updating shown in Figure 10.9-25 indicates significant reduction of the residual-moveout errors (compare, for instance, panel I-41 in Figure 10.9-25 with that in Figure 10.9-22).
  6. Perform 3-D prestack depth migration using the final 3-D velocity-depth model from step (e). The input to 3-D prestack depth migration is the prestack data set with signal processing applied as in Figure 10.9-5b. Figures 10.9-26 and 10.9-27 show selected inlines and crosslines, respectively, from the image volume in depth derived from 3-D prestack depth migration.
Figure 10.9-20  (a) Depth horizons interpreted from the volume of 3-D poststack depth migration of the 3-D zero-offset wavefield as in Figures 10.9-14 and 10.9-15 using the 3-D interval velocity field as in Figure 10.9-19; (b) the same depth horizons as in (a) with the additional phantom horizons shown in green used to sub-divide each of the layers in (a) into a set of 10 layers to preserve the lateral and vertical variations in the gradient of the 3-D interval velocity field as in Figure 10.9-19.

The deliverables from phase 4 — structural inversion, include a volume of structurally consistent 3-D interval velocity field and an image volume derived from 3-D prestack depth migration.

Structural and stratigraphic interpretation

By using the techniques we learned in interpretation of 3-D seismic data, we now interpret the image volume from 3-D prestack depth migration — structural intepretation based on picking a set of depth horizons and stratigraphic interpretation based on amplitude manipulation. Figure 10.9-28 shows six depth horizons that coincide with the geologic markers in the area. The top three horizons, DH1, DH2, and DH3, were delineated largely by using the seed detection technique that involves connecting neighboring voxels with amplitudes that are within a specified range (interpretation of 3-D seismic data). The bottom three horizons, DH4, DH5, and DH6, were delineated by line-based interpretation that involves creating a set of horizon strands along selected inlines and crosslines.

The effect of extensional tectonism is evident on the top three depth horizons. Note in particular the series of faults across the lowlands indicated in green. Examine the depth horizons in reverse order, starting with DH6 and ending with DH1, to reconstruct the structural evolution of the area. In the beginning, there was an erosional surface cut by a deep channel as seen on horizon DH6. There may have been a period of compressional tectonism that gave rise to the gently folded surfaces in the lowlands portion of horizons DH5 and DH4. Then, a reversal from compressional tectonism to extensional tectonism began to cause subsidence of the basin, which in turn gave rise to the formation of the series of faults parallel to the ancient shoreline. The highlands indicated by the orange on all horizons have remained as hiatus until the geologic times approximately coincident with the age of the overburden above horizon DH1.

Now split the image volume into subvolumes that represent individual depositional units as shown in Figure 10.9-29. The extensional fault patterns are observed on the top surfaces of these subvolumes down to horizon DH4. To investigate the interior of each of these depositional units, cut them into thin slices as shown in Figure 10.9-30a. Then, apply transparency (interpretation of 3-D seismic data) to each of the thin slices to identify depositional features of interest. Note the presence of a stream with its traverse orthogonal to the fault patterns. Also shown in Figure 10.9-30b is the map view of the thin slice that exhibits the stream channel disrupted by the faults.

As indicated in Figure 10.9-3, the seismic pathway we followed from phases 1 through 5 in this case study should lead us to constructing a model for the reservoir. Begin with a detailed delineation of the top and base of the deltaic sequence that represents the reservoir unit shown in Figure 10.9-31a. Structural interpretation of the image volume derived from 3-D prestack depth migration (Figures 10.9-26 and 10.9-27) yields the top and base surfaces shown in Figure 10.9-32. Then, define the interior geometry of the reservoir by splitting the deltaic sequence bounded by the two surfaces shown in Figure 10.9-32 into a set of thin layers as shown in Figure 10.9-31b. Cross-sections of the reservoir unit along selected inline traverses after sublayering are shown in Figure 10.9-33.

Next, extract the subvolume associated with the reservoir unit from within the image volume derived from 3-D prestack depth migration (Figure 10.9-34). You may want to use the subvolume to estimate a set of seismic attributes. The sublayers (Figure 10.9-33) then are populated, based on any available well data, and results of the analysis of seismic attributes, such as amplitude variation with offset (analysis of amplitude variation with offset) and acoustic impedance (acoustic impedance estimation), by the petrophysical properties of the reservoir unit, such as, porosity, permeability and fluid saturation, all of which are allowed to vary laterally within each of the sublayers. In the final chapter, we shall review prestack and poststack amplitude inversion methods to infer petrophysical properties of reservoir rocks from seismic data as an aid to reservoir modeling.

Exercises

Exercise 10-1. Does a zero-offset wavefield modeled from exploding reflectors contain multiples?

Exercise 10-2. Consider two image sections derived from (a) time migration, and (b) depth migration followed by scaling from depth to time. To perform the depth-to-time scaling, you may use an interval velocity field based on the rms velocity field to do the time migration or the velocity-depth model to do the depth migration. Which option would make the two image sections in time more compatible?

Appendix K: Seismic modeling

K.1 Zero-offset traveltime modeling

Seismic modeling essentially is a simulation of a recorded seismic wavefield, seismic amplitudes, or seismic traveltimes. The input to seismic modeling is a representation of the earth’s reflectivity and a velocity-depth model. Seismic migration is a process of estimating earth’s reflectivity from a recorded seismic wavefield using a velocity-depth model. Therefore, seismic wavefield modeling may be viewed as the reverse process of seismic migration. As such, both seismic migration and seismic wavefield modeling algorithms are based on the wave equation.

Given a seismic wavefield P(x, z = 0, t) recorded over time t, at the surface z = 0, and along the spatial axis x, seismic migration yields the earth’s reflectivity P(x, z, t = 0) based on a process of wavefield extrapolation in depth z and collecting the image at time t = 0 (migration principles). Conversely, the fundamental ingredient of the modeling process is wavefield extrapolation in time t and collecting the result at depth z = 0. Both processes use wave equation as the basis for wave extrapolation.

Seismic modeling is different from data modeling (Appendix J). The latter involves, given an observed data set d, an estimation of a set of parameters p that are used to construct a model d′ of the observed data set d, such that the difference between the observed data set d and the modeled data set d′ is minimum based on a specific mathematical norm.

Throughout this book, we have seen numerous examples of seismic modeling:

  1. to explain a process such as deconvolution (Figure 2.1-3) or migration (Figure 4.0-8),
  2. to test an algorithm such as predictive deconvolution (Figure 2.4-13) or implicit frequency-space 3-D poststack time migration (Figure 8.4-2), or
  3. to understand a structural or stratigraphic phenomenon that may be of interest in exploration.

Just as there are several approaches to solving the wave equation for migration, there also are several types of modeling techniques. There are modeling techniques based on the Kirchhoff integral [1], finite-difference [2], and fk domain [3] solutions to the wave equation. The algorithms based on the scalar (acoustic) wave equation (Section D.1), which describes P-wave propagation, are suitable for structural modeling in which amplitudes are not as important as traveltimes. The algorithms based on the elastic wave equation (Section L.2), which describes both P- and S-wave propagation, are suitable for detailed stratigraphic modeling in which amplitudes are as important as traveltimes. Modeling based on one-way wave equations does not include multiples, while modeling based on two-way wave equations includes multiples in the simulated wavefields. In this appendix, we shall not discuss the details of specific algorithms for seismic modeling which span a broad range of applications including 2-D and 3-D, zero-offset and nonzero-offset, acoustic and elastic simulation. Instead, we shall provide examples of the most common seismic modeling strategies — zero-offset traveltime modeling, zero-offset and nonzero-offset acoustic wavefield modeling, and elastic modeling.

Shown in Figure K-1a is a velocity-depth model for a salt diapir with an overhang structure. The zero-offset traveltime response shown in Figure K-1b is created by normal-incidence ray tracing. Note that the top-salt boundary has given rise to a complex and multivalued traveltime trajectory. Note also that the base-salt boundary is flat and continuous in the velocity-depth model (Figure K-1a), whereas the reflection traveltime follows a discontinuous and multivalued trajectory (Figure K-1b).

Figure K-1  (a) A velocity-depth model of a salt dome with overhang, (b) the zero-offset traveltime response.

Zero-offset traveltime modeling using normal-incidence rays is a very useful and trivially simple tool for understanding the complexity of a reflection traveltime in field data. The disruptive behavior in traveltime trajectory associated with a layer boundary below a complex overburden as in Figure K-1b is observed also in real data (Figure 10.1-1).

K.2 Zero-offset wavefield modeling

Recall from introduction to migration that a stacked section often is assumed to be a close representation of a zero-offset wavefield. A modeled zero-offset wavefield therefore can be used to test poststack migration algorithms. Zero-offset wavefields can be simulated very efficiently using the exploding reflectors, also discussed in introduction to migration.

Wave-equation datuming [4], which was described in layer replacement, can be used to perform the simulation based on exploding reflectors. In particular, the datuming approach can propagate a wavefield from one irregular interface to another. Consider zero-offset modeling using the datuming technique of the velocity-depth model shown in Figure K-2a. Horizons 2 and 3 are the top and base of a salt dome. Start with the receivers situated along horizon 3. The corresponding zero-offset section (Figure K-2b) contains the reflection from the bottom of the velocity-depth model at z = 4000 m (not shown in Figure K-2a). Take this wavefield and extrapolate it to a new datum, horizon 2, using the salt velocity (5000 m/s). The resulting zero-offset section (Figure K-2c) contains the reflection (the deeper one) from the bottom of the velocity-depth model (z = 4000 m) and the reflection (the shallow one) from the base of the salt (horizon 3). Finally, extrapolate this wavefield (Figure K-2c) from horizon 2 to the surface (horizon 1 at z = 0) using the overburden velocity (3000 m/s) to get the 2-D zero-offset section in Figure K-2d. This section contains reflections from both the top and base of the salt. (The reflection from the bottom of the model arrives after the latest time shown on this section.) Note the velocity pull-up along the reflection from the base of the salt dome. Proper imaging of the top of the salt dome can be achieved by time migration (introduction to migration), while proper imaging of the base of the salt requires depth migration (introduction to earth imaging in depth).

K.3 Nonzero-offset wavefield modeling

Understanding complexities of recorded wavefields clearly requires nonzero-offset wavefield modeling. A finite-difference technique for modeling acoustic and elastic wavefields is described by Kelly [2]. Figure K-3 shows an example of acoustic modeling of a complex structure associated with overthrust tectonics. A seismic line is simulated over a 2-D complex structure (Figure K-3a). Selected common-shot (Figure K-3b) and CMP gathers (Figure K-3c) from this simulation show the many complexities in the arrivals. Since this is a two-way acoustic simulation, the modeled gathers contain not only primaries but also multiples. The zero-offset and stacked sections associated with this nonzero-offset data are shown in Figure K-4. Note the broad traveltime trajectories associated with the tight imbricate structures in the velocity-depth model (Figure K-3a).

An example of a nonzero-offset modeling application of wave-equation datuming is provided in Figure K-5. The shot gathers in Figure K-3b are computed to a flat datum level z = 0. Better simulation of the actual field conditions requires that the gathers be computed using an irregular topography. To do this, we can upward continue the shots and receivers to the new irregular datum represented by the topography shown in Figure K-5a, then compute the shot gathers in Figure K-5b and sort them to the CMP gathers shown in Figure K-5c. Compare Figures K-3b and K-3c with Figures K-5b and K-5c, and note the traveltime distortions.

Figure K-6 shows a velocity-depth model associated with a salt sill structure caused by salt tectonics in the Gulf of Mexico. Note that velocity variations in some parts of the sedimentary section are structure independent and represent overpressured zones. Selected common-shot gathers shown in Figure K-7 have been created by two-way acoustic wavefield modeling [5]; therefore, they contain both primaries and multiples. Each shot gather represents a modeled wavefield. Note the complex events in the gathers above and in the vicinity of the salt sill shown in Figure K-6.

Shown in Figure K-8 are selected CMP gathers sorted from the modeled shot gathers as in Figure K-7. Observe the events with complex moveout in the gathers above and in the vicinity of the salt sill. The zero-offset section obtained by collecting the zero-offset traces from the modeled shot gathers is shown in Figure K-9a, and the stacked section obtained from the CMP gathers as in Figure K-8 is shown in Figure K-9b. A zero-offset wavefield simulated by exploding reflectors does not include multiples, because the exploding reflectors are associated with simulation based on one-way wave equation [6]. When the simulation is based on the two-way acoustic wave-equation as in Figure K-9, the zero-offset and stacked sections both include primary and multiple reflections. Since it is wavefield modeling, not just traveltime modeling, the simulated shot gathers (Figure K-7), the associated CMP gathers (Figure K-8), and sections (Figure K-9) all contain the diffractions caused by the reflector discontinuities in the velocity-depth model (Figure K-6).

K.4 Elastic wavefield modeling

Elastic wavefield modeling primarily is used to understand the effect of lithology and pore fluids on seismic amplitudes (seismic resolution and analysis of amplitude variation with offset). Sherwood [3] developed an f − k method for nonzero-offset modeling of elastic waves in a 2-D horizontally layered medium. Figure K-10 shows five shot gathers derived from an earth model represented by a vertically varying velocity function that includes a water layer with five different thicknesses. The water depths are 5, 10, 15, 20, and 50 m. Note the guided wave energy, which is especially prominent in gathers corresponding to water depths of 5, 10, and 15 m. These gathers contain all primaries, both P-waves and S-waves, as well as all possible multiples and converted modes. By examining such modeled data, we can better understand the nature of coherent noise (guided waves and multiples) in both land and marine environments.

Figure K-11  Elastic modeling examples: (a) the response of a clastic section with a low-velocity shallow layer on the left and the response of the same clastic section with a fast-velocity shallow layer on the right. What is the low-frequency low-group velocity dispersive energy? (b) Seismic response of an all-shale section, (c) seismic response of a sand-shale section see text for details; [3].

A more interpretive application of elastic modeling is shown in Figure K-11 [3]. The synthetic shot gather on the left in Figure K-11a is from a clastic section with a shallow layer of low P-wave velocity. This layer has been replaced with a fast-velocity limestone for the record on the right. Note the invasion of the large-offset primary reflection data with coherent noise that is associated with this limestone layer. In the field, limestone on the surface often generates a large amount of coherent noise.

Figures K-11b and K-11c show two synthetic shot gathers for which the upper part of the depth model consists of a 50-ft water layer underlain by a 1020-ft shale section with a P-wave velocity of 5500 ft/s. The deeper portion of the model for Figure K-11b is an all-shale section with P-wave velocity increasing from 5600 ft/s on the top to 7700 ft/s on the bottom of the layer. The primary reflections PP between 380 and 850 ms are associated with this all-shale sequence. Figure K-11c shows the effects of including a 30 percent sand layer between 660 and 850 ms. The P-wave reflections PP in Figure K-11c show stronger amplitudes at large offsets. An analysis of amplitudes as a function of offset can provide hints for determining the sand-shale ratio, as well as fluid content, in some cases. Because the effects are complex, this type of modeling can be helpful in analyzing amplitude variations with offset. Also note the relatively strong converted PS and SP waves on the record from the sand-shale model. Modeling of this type also is useful when analyzing converted waves in multicomponent reflection data (4-C seismic method).

See also

References

  1. Hilterman, 1970, Hilterman, F. J., 1970, Three-dimensional seismic modeling: Geophysics, 35, 1020–1037.
  2. 2.0 2.1 Kelly et al., 1976, Kelly, K. R., Ward, R. W., Treitel, S. and Alford, R. M., 1976, Synthetic seismograms: A finite-difference approach: Geophysics, 41, 2–27.
  3. 3.0 3.1 3.2 3.3 Sherwood et al., 1983, Sherwood, J. W. C., Hilterman, F. J., Neale, R. N. and Chen, K. C., 1983, Synthetic seismograms with offset for a layered elastic medium; Offshore Technology Conference, Paper 4508.
  4. Berryhill, 1979, Berryhill, J., 1979, Wave-equation datuming: Geophysics, 44, 1329–1344.
  5. 5.0 5.1 5.2 5.3 5.4 5.5 O’Brien and Gray, 1996, O’Brien, M. and Gray, S., 1996, Can we image beneath salt?: The Leading Edge, 17–22.
  6. Claerbout, 1985, Claerbout, J. F., 1985, Imaging the earth’s interior: Blackwell Scientific Publications.

External links

find literature about
Ageary
SEG button search.png Datapages button.png GeoScienceWorld button.png OnePetro button.png Schlumberger button.png Google button.png AGI button.png