01 August, 2014

Spectral Decomposition for 3-D Geomodeling

 

In recent years, a number of workflows have been developed to use spectral decomposition for 3-D geomodeling to extract architectural elements such as channels.

Spectral decomposition attributes use a localized time-frequency decomposition within individual traces to characterize seismic data by its spectral content.

Since localized spectral content can be linked to physical phenomena of channels such as fill velocity and thickness, such techniques have been popular in the past decade for characterizing and interpreting of fluvial systems including turbidites.

The original work using spectral decomposition focused upon qualitative interpretation upon time slices and/or phantom horizons with a goal of broad scale interpretation and reservoir description.

In recent years, however, a number of workflows have been developed to use spectral decomposition for 3-D geomodeling to extract architectural elements such as channels. Most of these methods are proprietary – and, to my knowledge, few have been documented or published.

A typical workflow might involve the use of the instantaneous frequency attribute to determine the dominant frequency associated with a channel. Once this information has been ascertained, an interpreter then might use some spatially associated thresholding (i.e. amplitude clipping) method to extract the channel.

These results then would undergo further processing, such as smoothing or denoising, to produce a final interpretation.


Figure 1 illustrates how a workflow like this would be constructed:

  • Figure 1a shows a time slice of seismic amplitude data presenting a portion of a channel system from the West Cameron block of the shelf region of Louisiana.
  • Figure 1b shows a vertical slice cutting across this system (line A-A’) perpendicular to the direction of paleo-flow. At this location, there is one distinct channel (channel a) and two closely associated channels (channels b and c) with a narrow inter-fluvial region.
  • Figures 1c and d show time and vertical slices through the instantaneous frequency attribute at locations matching those shown in figures 1a and b, respectively.

An examination of these images suggests that the channels might be associated with frequencies in the 40-60 Hz range. This relationship is much stronger for channels “a” and “c,” which have better defined tops than channel “b.”

Note that 50 Hz corresponds to a period of 20 ms, which is the approximate trough to trough thickness associated with all three channels.

One possible pitfall is related to the spatial uncertainty associated with spectral attributes. Specifically, since spectral content is calculated using a vertical window, the calculated spectral information might be too imprecisely positioned to allow for good modeling in a vertical sense.

This is further compounded by the inclusion of spectral phase, which might allow for the fitting of a rotated wavelet that is not properly centered upon the architectural element of interest.

The net result of these problems could be an extracted model that looks geological but that is not properly located vertically. Better understanding of the vertical localization of spectral decomposition attributes is necessary to better understand if this is a true concern.


To assess whether a 3-D geomodeling workflow was feasible and valid, I decided to look at the vertical characteristics of some spectral magnitude attributes for my channels of interest.

My primary question was, “Are these channels completely and properly isolated vertically using spectral decomposition attributes?”

One aspect of the vertical localization will depend upon the algorithm, as there are a number of methods to calculate spectral decomposition attributes. For the purposes of my study, I chose to use a matching pursuit (MP) algorithm, which is regarded as having superior localization.

MP spectral decomposition works by decomposing a seismic trace into a set of time-shifted wavelets selected with replacement from a library of pre-computed wavelets. The decomposition is done iteratively by fitting component wavelets and calculating residuals until the residuals are sufficiently small.

Figure 2 shows some spectral magnitude attributes as calculated using a MP spectral decomposition algorithm calculated using a library of Ricker wavelets.

♦ Figure 2a shows the same vertical slice of seismic amplitude data shown in figure 1b. In this slice, I have shown an interpretation of the bottom of the channels picked at a trough (green horizon) and a pick of the top based upon a peak (yellow horizon).

Note that the trough-to-peak distance for all three channels is approximately 30 ms, while the trough-to-trough distance in these channels, as previously noted, is approximately 20 ms. These correspond to frequencies of approximately 30 Hz and 50 Hz, respectively.

♦ Figure 2b shows the 30 Hz component. In this, channels “a” and “c” are relatively well defined by a low magnitude region. This region is especially well defined for channel “c.” Channel “b” is not as well defined.

A low magnitude region matches our intuition for a peak-to-trough region when decomposed using a unimodal wavelet with an internal peak and trough as the internal components will destructively combine when convolved with the two positive sides of the wavelet.

Channel “c,” while well defined, is not delineated to its proper top.

♦ Figure 2c shows the 50 Hz component. In this, all of the channels are relatively well defined, this time by a high magnitude region.

Again, channel “c” is particularly well defined. However, there is some less well defined regions along a portion of the top of channel “a.” Either channel b’s top or bottom would be somewhat misplaced spatially, depending upon the choice of cutoff for segmenting.

Additionally, the inter-fluvial region between channels b and c is sharply defined by a particularly strong low magnitude event.


In conclusion, these results show that there is a promise to algorithms using spectral decomposition attributes to do 3-D geomodeling of channels.

In the case of the 50 Hz decomposition, the theoretical analysis of the wavelet is based upon the consideration of the bottom of the channel and an internal reflector, while the highlighted region appears to cover the full range of the channel.

I believe this has to do with the phase of the fit wavelets. This is a conjecture that suggests future work to better understand the role of phase in modeling channels.

Finally, while the results do show vertical localization, all of the images of the attributes are noisy. This emphasizes the importance of good post-processing workflows to provide parsimonious and geologically sound models.

Additionally, as some boundaries are misplaced vertically, it is my belief that naïve algorithms will give less desirable results than one that incorporates some heuristical processing to ensure geologically and geophysically consistent models.


I would like to thank Victor Aarre and AAPG members Kurt Marfurt, and Roger Slatt (Honorary) for many insightful discussions concerning this and related topics. Additionally, I would like to thank Schlumberger for providing the data used in this work.

Finally, I would like to thank the sponsors of the Attribute Assisted Seismic Processing and Interpretation (AASPI) Consortium for their ongoing support.

(Editor’s note: Brad Wallet is director of the Crustal Imaging Facility at the University of Oklahoma, Norman, Okla.)

Figure 2 – Panel a) shows the same vertical slice of seismic data shown in figure 1b with an interpretation of the channels based upon a top peak and a bottom trough. Panels b) and c) show the 30 and 50 Hz spectral components, respectively. Note that the channels are generally well interpretable in the spectral components, though some ambiguity exists in various cases.

Figure 1 – Panel a) shows a time slice of seismic data containing channels, and b) shows a vertical slice of the same seismic data taken perpendicular to the paleo-flow of the channels. Three channels – a, b and c – are identified. Panels c) and d) show the instantaneous frequency attribute for the views shown in a) and b), respectively. These images suggest frequencies in the 40-60 Hz range are likely associated with the channels.