Advances in 3D seismic tomography enhance subsurface interpretation

To accommodate different geological conditions (e.g. deepwater compaction, hard rock environments, obtrusive structures, and local low/high velocity anomalies) today's most advanced tomography offers both grid-based and model-based solutions. For large-scale models, especially those that require velocity parameter updates on high-resolution grids, the size of the set of equations to be solved (i.e. the tomography matrix) can be too large to be handled even by the largest super-computers. Consequently, different approaches are required. New "direct" tomography methods propose to overcome these inherent limitations by reducing the memory and disk space required through more intense computation.

Subsurface velocity model

In pursuit of the optimum seismic image and model properties, geoscientists have discovered the value of more rigorous velocity model parameterizations. The subsurface velocity model generally is defined as a set of layers where the velocity at each layer can vary both laterally and vertically. Anisotropy effects (the change of velocity with direction) also can be assigned to specific sedimentary layers, and generally are modeled by a set of mutual perpendicular axes. Vertical transverse isotropy (VTI) and tilted TI (TTI) axes of symmetry are the most common anisotropic systems, especially within compacted shale/sand layers. In these layers the axial velocity is slower than the perpendicular velocity as the seismic wavefronts respond to differences in the "fabric" or depositional layering in the geologic formation. Horizontal TI (HTI) model representations, accounting for vertical fractures (e.g., fractured carbonate or fractured shale layers), are considered a specific case of TTI while orthorhombic and tilted orthorhombic (TO) models simultaneously account for layering and fractures within the same layers.

The background starting velocity model is assumed to be somewhat smooth, enabling ray tracing. TTI (and TO) anisotropic models require reliable information about the orientation of the subsurface local reflectors. This can come from extractions of subsurface dip, azimuth and continuity, from seismic data as seismic attributes, or from interpretation data. Advanced tomography systems offer different tools to automatically compute, smooth, and condition these structural attributes.

The tomographic approach

Traveltime tomography is based on solving a large set of linear equations that relate the desired (unknown) subsurface velocity model update parameters, with the input data containing traveltime errors along different types of rays traveling across the model. These traveltime errors result from incorrect velocity model assumptions. Tomography is formulated as a geologically constrained inversion by using a priori information about the subsurface model, such as the depth of well markers and spatial variation characteristics of the velocity field within the different geological layers.

Using tomography, it is possible to simultaneously update all anisotropic velocity components, such as the axial velocity Vp, and anisotropic parameters (δ and ε Thomsen TI) needed to define the diversity of velocity models described previously. Different variances can be set for each of the material parameters at each layer, which makes the tomography versatile and allows full control over the resulting updated model.

For example, sea water velocity can be fixed (zero variance), the velocity variation within salt diapirs, basalt layers or hard carbonate rocks can be forced to be uniform, and velocity updates within sedimentary layers can have lateral and vertical variations, allowing local anomalies.

Traveltime errors

Reflection tomography is designed to handle reflection seismic data, where traveltime errors along the rays (the required input for tomography) are measured generally from residual moveouts (RMO) automatically picked along common image gathers (CIGs). The CIGs may be generated by different types of seismic depth migrations.

Reflection tomography supports surface offset domain (computed by Kirchhoff migrations), image offset domain (computed by wave equation migrations such as RTM), and subsurface angle domain (computed by ray-based or wave-based local angle domain imaging solutions) CIGs. Depending on the input seismic data and the type of CIGs, the solution can contain single-, multi- and even full-azimuth RMO distributions. The richer the information about the traveltime errors along the reflected rays (e.g., from all opening angles and all azimuths), the higher the accuracy, resolution, and certainty of the updated model.

Overall, the aim of tomography is to update the subsurface model with a velocity field that minimizes traveltime errors; in other words, running the migration with the updated velocity field will flatten events within CIGS and remove the RMO effect.

The RMO can be defined in both parametric form (up to five parameters describing the second and fourth order RMO variation along the angle/offset axis and the variability with azimuth) and non-parametric form, where the RMO values are automatically picked and stored for each trace. The quality and reliability of the RMO data are keyδ to the success of the tomography in converging into the "right" model.

Adding traveltime information from refraction data and long offset walk-away check-shots (if available) can be extremely beneficial. The refraction traveltimes can significantly increase the accuracy and resolution of the shallow velocity model, and the walk-away check-shot traveltimes can contribute to accurate determination of anisotropic parameters.

Advanced 3D visualization provides a comprehensive solution for quality control and editing of the automated RMO picked data, where the RMO curves can be displayed with a massive amount of gathers. The systems must be able to filter out outliers according to predefined attributes.

Tomography equations

Tomography requires an efficient and accurate isotropic/anisotropic ray tracing technique with a satisfactory level of robustness, even when applied within complex background velocity models. Ray tracing is performed from subsurface points along the major reflectors (horizon interfaces) and along coherent events within the layers (inter-layer horizons), which are automatically picked along the 3D subsurface image (migrated volume). Specular rays (ray pairs that obey Snell's law at the reflection surfaces) are then traced up to the surface, where each ray forms a linear relation between the velocity model perturbations in its vicinity and the known total traveltime errors. The ensemble of all rays traced from all subsurface points and for all angles/offsets, forms the global set of tomography equations (tomography matrix) to be solved (inverted).

Today's most advanced tomography systems offer both grid-based and model-based methods, depending on the problem to be solved. For example, Gulf of Mexico sediments affected by long period compaction are normally parameterized by a Cartesian grid, whereas North Sea subsurface models are usually characterized by the aging of layers using a layer-based representation. The tomography workflow also supports hybrid approaches where, for example, shallow velocity anomalies can be parameterized with grids, and deeper parts with geological layers. Unlike grid-based tomography, in which only the velocity parameters are updated on a predefined coarse 3D grid, model-based tomography updates both the velocity field along the geological layers and the location of subsurface horizons.

Implementation

Tomography is conventionally implemented in two independent stages:

Constructing the tomography equations (building the tomography matrix) and solving the tomography equations (inverting the tomography matrix). The first part includes ray tracing from all subsurface points and for all angles/offsets – the most time-consuming part. The second part reads the tomography matrix and inverts it with different input parameters (constraints): data and model variances, resolution, imposing the velocity updates to follow the background structure, and other geological constraints, avoiding recalculation of the tomography equations.

Despite the advantages of 3D grid-based tomography, the size of the tomography matrix for large-scale models, which is quadratic with the number of model parameters, can be too large to be handled even by the largest super-computers. One way to overcome this is a novel approach referred to as "direct" tomography. While the resulting updated subsurface velocities are identical to those of conventional 3D grid tomography, the difference is in the implementation. In direct tomography, the tomography matrix is never generated explicitly, dramatically reducing the memory and disk space required by the application, where inversion (including ray tracing and the construction of the tomography equations) is performed in a single stage using an iterative process.

The main and obvious advantage of direct tomography is that it can handle large 3D updated grids (either due to the size of the survey or the requirements for high resolution) that could not be solved in a practical way before. Moreover, this implementation is efficient when using clusters with massive amounts of nodes. The runtime of direct tomography mainly relates to the CPU power in globally running the ray tracing many times (for each iteration), whereas conventional 3D grid-based tomography suffers mainly from being I/O bounded, where most of the time is spent reading and manipulating the tomography matrix (many terabytes of data).

Workflows

Tomography supports different workflows, depending on the type and complexity of the subsurface geological models and the available input data. Overall, a "standard" workflow should start with a simple, long-wavelength, background velocity model that includes the main trends of the model. The tomography, after several iterations, should bring in the small anomalies (if any) and the high-resolution features. It is also recommended to divide the workflow to independently update the shallow, medium and deep parts of the model, where a final iteration can be performed over the whole model.

Anisotropy velocity model determination requires access to some well data (well markers or check-shots) in order to define the axial velocity Vp and the small offset (near axial direction) δ Thomsen parameter values. Long offset RMOs are needed to define the far offset (near horizontal) ε Thomsen parameter values.

When the input seismic data contains surface offset azimuthal information such as 3D land data or wide (rich) marine azimuth data, full-azimuth imaging in the local angle domain can provide extremely useful continuous full-azimuth angle gathers. These can be used to better define the different subsurface anisotropic model parameters, since the additional azimuth information minimizes the non-uniqueness of the velocity model parameterization. The full-azimuth angle domain RMO extracted along these gathers, and the amplitude inversion applied to them (AVAZ, or amplitude versus angle and azimuth) can be further used to detect velocity directivity and geomechanical properties. These may point to stress/fracture systems, and help estimate their orientation and intensity.

Tomography also can be used to simulate different model scenarios while preserving the total traveltime (traveltime errors along the rays are zero). In this case, image gathers should be flat for each simulated model. An example of this type of tomography implementation is when the input consists of only depth mistie maps computed from sparse selected well markers.

The ideal system

Velocity model determination workflows, especially when using tomography, comprise several key stages. Each stage needs to be optimized to enable fast and reliable QC that will ensure a convergence process. The ideal tomography system includes:

• A 3D geomodeling solution to build and update geologically plausible subsurface structure models
• A comprehensive anisotropic ray tracing system
• Efficient and accurate pre-stack depth migrations for generating the subsurface image and the offset/angle domain CIGs
• An optimized processing and interpretation system for automatic event picking in both post-stack and pre-stack domains
• A comprehensive visualization and interpretation system to analyze, QC, and edit the picked events
• High-performance computing to use multi-node clusters to solve tomography equations.

By using a fully integrated system, the velocity model determination process reduces uncertainty and enables the geoscientist to define more reliable and accurate subsurface models.