Subtle fault detection plays a vital role in reservoir development studies because faults may form baffles or conduits that significantly control how a petroleum reservoir is swept. Small-throw faults are often overlooked in interpreting seismic amplitude data. However, seismic attributes can aid in mapping small faults. Over the years, dozens of seismic attributes have been developed that offer additional features for interpreters with as-sociated caveats. Using the Maui 3D seismic data acquired in the Offshore Taranaki Basin, New Zealand, we have generated seismic attributes that are typically useful for fault detection. We find that multiattribute analysis provides greater geologic information than would be obtained by the analysis of individual attribute volumes. We extract the geologic content of multiple attributes in two ways: interactive corendering of different seismic attributes and the unsupervised machine learning algorithm self-organizing maps (SOM). Corendering seismic attributes that are mathematically independent but geologically interrelated provides a well-integrated structural image. We suggest eight combinations of 16 various attributes useful for a human interpreter with interesting fault and fracture detection. Current interpretation display capabilities constrain corendering to only four attribute volumes. Therefore, we use principal component analysis and SOM techniques to efficiently integrate the geologic information contained within many attributes. This approach gathers the data into one classification volume based on the interrelationships between seismic attributes. We show that our resulting SOMclassification volume better highlights small faults that are difficult to image using conventional seismic interpretation techniques. We find that SOM works best when a fault exhibits anomalous features for multiple attributes within the same voxel. However, human interpreters are more adept at recognizing spatial patterns within various attributes and can place them in an appropriate geologic context.
Subtle faults and fractures as well as their clusters and zones are structural discontinuities within a geologic formation that show little or no relative movement along the discontinuity. The identification of small faults (i.e., those with small offsets near or below seismic tuning thicknesses) and fractures plays a crucial role in reservoir characterization, reservoir modelling, and identifying hydrocarbon migration pathways and bypassed oil accumulations, but it can be overlooked in the conventional interpretation of seismic amplitude data (Al-Dossary and Marfurt, 2006; Li and Lu, 2014). In this paper, we do not differentiate among small faults, fractures, or their zones; instead, we refer to these subtle structural features interchangeably as “faults” or “fractures.”
You must log in to read the rest of this article. Please log in or register as a user.
Over the past few decades, geoscientists have de- voted significant effort to develop seismic attributes to better map faults. A prominent example is the development of coherence algorithms that quantify the lateral change in the seismic waveform, and/or amplitude (Bahorich and Farmer, 1995; Luo et al., 1996; Gersztenkorn and Marfurt, 1999; Randen et al., 2000). However, coherence attributes display segmented anomalies or perhaps no anomalies in the presence of similar seismic reflectors across faults (Mai et al., 2009; Libak et al., 2017; Marfurt, 2018). Partyka et al. (1999) emphasize the spectral phase components to help to in- interpret faults. Henderson et al. (2008) use red-green-blue (RGB) blending of spectral magnitudes to detect faults. Taking advantage of spectral decomposition to enhance the coherence image, Li and Lu (2014) combine spectral decomposition and coherence to detect faults at different temporal scales. Alternatively, a multispectral coherence volume can be computed by constructing a multispectral covariance matrix along the structure dip calculated from the spectral voices (Mar-furt, 2017; Li et al., 2018). Chopra and Marfurt (2019) extend the multispectral coherence concept to multi offset and multi azimuth coherence analysis and show superior results.
Despite these efforts, it remains challenging to use coherence attributes to detect faults below the limits of seismic resolution (Gao, 2013). More geometric attributes have been developed to quantify the lateral variation of the seismic reflector geometry such as dip, curvature, and aberrancy attributes, which aid in illuminating small faults and fractures (Roberts, 2001; Al-Dossary and Marfurt, 2006; Di and Gao, 2016; Qi and Marfurt, 2018). The grey-level co-occurrence matrix (GLCM) is a texture attribute that measures how often various combinations of amplitude sample values occur within an analysis window (Gao, 2011; Matos et al., 2011). For the human interpreter, GLCM attributes may provide less useful images compared to other geometric attributes, which provide clear and easily understood geomorphic features such as channel edges, faults, and folds. However, these texture attributes work well for supervised and unsupervised classification machine learning algorithms (Marfurt, 2018).
With continuous advancements in computer capabilities and graphics, many computer-aided fault detection algorithms have been developed to enhance fault edges and interpret faults. For instance, Pedersen et al. (2002, 2003) introduce the application of the anti-tracking algorithm to variance volumes to enhance and sharpen fault edges. Barnes (2006) develops a discontinuity filter to enhance the steeply dipping discontinuities associated with faults and cancel all other discontinuities. Aqrawi and Boe (2011) compare the application of the ant-colony optimization algorithm on variance and Sobel filter volumes and find that the Sobel filter showed better re- results compared to the variance attribute, but ant- colony optimization is highly sensitive to noise. Dewett and Henza (2016) use a swarm intelligence algorithm for lineament connection and interpolation between discontinuous events within multiple frequency-based attributes and integrate the results by using normalization, addition, and rescaling besides the self-organizing maps (SOM) clustering technique. They conclude that SOM can retain some information of the original frequency volumes. Laudon et al. (2019) present a case study on the application of SOM to four geometric attributes to detect faults and fractures within the main reservoir intervals of the Denver-Julesburg Basin. They select the most positive curvature, energy ratio similarity, GLCM entropy, and GLCM homogeneity based on principal component analysis (PCA) analysis applied to various geometric attributes (Laudon et al., 2019).
In another approach, convolutional neural networks (CNNs) using pattern recognition have been trained on synthetic seismic data and applied to seismic data sets for fault detection (Pochet et al., 2019; Wu et al., 2019; Zhao, 2019). Their synthetic data generation included faults with random dip azimuth, dip magnitude, dis- placement along with a dip in addition to a certain degree of normal fault drags. Although the CNN showed sharp fault edges, it failed to detect small dip angle (approximately 9o) faults (Zhao, 2019). A limitation of the CNN workflow from Wu et al. (2019) is that the synthetic data did not include low-angle faults; thus, the trained CNN might not detect thrust or listric faults. In addition, when the CNN-trained model is applied to unfamiliar data with slightly different characteristics compared to the training data, the accuracy of the prediction drops dramatically (Lowell and Erdogan, 2019).
Laudon et al. (2019) show that a less tedious and more efficient way to detect small faults is to incorporate multiple seismic attributes in a machine learning technique instead of using a single seismic amplitude volume. However, the selection of the most suitable attributes to be used as inputs to machine learning algorithms is an open question that highly controls the results. In this paper, we investigate most of the available seismic attributes for fault interpretation starting from simple high spatial residual maps (e.g., horizon-based attributes) and moving to more advanced geometric seismic attributes such as dip, coherence, amplitude gradients, curvature, amplitude curvature, aberrancy, and GLCM attributes. We use a 3D seismic data set from the Maui field in Offshore Taranaki Basin, New Zealand. We in- investigate the various seismic attribute volumes generated for the C Sand reservoir of the Maui field and discuss the most useful attributes that can be used qualitatively by a human interpreter and quantitatively by PCA and SOM techniques and compare to a previous study (Laudon et al., 2019). We integrate the geologic information derived from multiple seismic attributes by using a multiattribute display and PCA and SOM techniques, and we compare it to workflow following a previous study (Laudon et al., 2019). We suggest a workflow that enables interpreters to apply PCA and SOM on the most appropriate seismic attributes to obtain a single classification volume that best shows small faults.
The Taranaki basin covers 200,000 km2 and is the only major hydrocarbon-producing basin of New Zealand (Reilly et al., 2016). It is located mostly offshore of New Zealand’s west coast (Figure 1) (Baur et al., 2010; Strogen et al., 2015). The structural development of the Taranaki basin was highly influenced by Austral- in-Pacific plate divergence during the Late Cretaceous- Paleocene period. The structural history of the Taranaki basin is complex with various tectonic elements affecting the basin which include rift transform, platform subsidence, passive margin, subduction, convergent transform, volcanic-arc, fold-thrust, and back-arc rift (Haque et al., 2016).
This study focuses on the Maui field, which is one of the largest New Zealand gas fields. OMV New Zealand owns and operates the Maui field, having acquired it from Shell Exploration NZ at the end of 2018. Maui covers approximately 1000 km2 and is located 40 km off of the west coast of the North Island at approximately
100 m water depth (Figure 1a) (Pannett et al., 2004). Based on the latest reserve estimation report released by the Ministry of Business, Innovation, and Employment (MBIE) in 2019, the estimated remaining gas re- serves of Maui field are 104.4 petajoules (99 bcf) and the estimated remaining oil reserves are 4.5 million barrels. The field produces gas from three main sand reservoirs called the Mangahewa, Kaimiro, and Farewell Formations of the Kapuni groups (Figure 1b). In the Maui field, Eocene-age shore-face to marginal marine Mangahewa and Kaimiro reservoirs are named the C Sand and D Sand, respectively, whereas the Paleocene Farewell sand is named the F Sand. The C Sand reservoir, which is the focus of this research, contains a shallower and larger gas column compared to the D Sand. However, it has a very thin oil-rim in the Maui B region. The upper D Sand reservoir has a smaller gas accumulation. Oil is produced from the lower D Sand and the F Sand reservoirs in the Maui B region (Pannett et al., 2004). Figure 2a displays the seismic-to-well tie for one of Maui wells that shows the seismic picks for the top of the C Shale, C Sand, and D Sand.
The Maui field is composed of two gentle anticlinal features: the Maui A region in the northeast and the Maui B region in the southwest (Figure 2b) that form a closure of 150 km2. The field is structurally complex.
The Maui A region is bounded to the east by the Plio- Pleistocene Cape Egmont normal fault (CEF) and the Pungawerewere normal fault (PF). The Maui B region is bounded to the southwest by the Late Miocene Whi- tiki reverse fault (WF) and is located within a hanging- wall anticline caused by inversion of a Cretaceous normal fault (Figures 1a, 2b, 2c, and 2e) (King and Thrasher, 1996). Three main tectonic phases affected the Maui field and formed the present-day Maui structure. The Late Cretaceous to Paleocene rifting phase caused by the Australian-Pacific plate divergence (Ha-que et al., 2016) formed graben and/or half-graben-like basins that are controlled by the Whitiki (WF) and Cape Egmont (CEF) faults (Figure 2c). The Late Eocene to Miocene compressive tectonic phase resulted from the subduction of the Pacific plate underneath the Australian plate. The Plio-Pleistocene differentiation phase resulted in extensional forces in Maui North and compressive forces in Maui South (Holt and Stern, 1994; Ha-que et al., 2016). These later tectonic movements led to the activation of preexisting structural fabrics and formed new faults. The Maui structural closure commenced in the Late Miocene due to the WF contraction movement and continues to develop at the present day because of the normal faulting activity (Nodder, 1993; Nicol et al., 2005).
A seismic line extracted from seismic volume over the Maui A region shows that the reservoir is influenced by normal faults (Figure 2d). These faults have a northeast-southwest structural trend with a few faults having a northwest-southeast structural trend. The reservoir layers dip (10°–15°) toward the east and west. A seismic line crossing the Maui B region shows slightly different structural features caused by the influence of structural mechanisms (Figure 2e). Although this part of the field is highly affected by a regional reverse fault (WF), a folding mechanism seems to have formed stratigraphic monocline units that have 15°–20° dip angles to the west. This indicates that the northern part of the field has a greater intensity of faulting compared to the southern part of the field. Some listric faults with rollover anticlines have been observed in the central part of the field (Haque et al., 2016).
Two 3D seismic surveys were acquired over the Maui field and made available for research by the MBIE, New Zealand. The first 3D seismic survey was acquired in 1991 and covers 1000 km2 and the second 3D seismic data cover 480 km2 over the Maui A and Maui B structures (Figure 1a). The 1991 was acquired as a quad-quad with four streamers and two shots. The shot interval was 18.75 m, and the receiver group interval was 25 m. The 2002 seismic data were acquired by using two shots and six streamers with 18.75 and 12.5 m shot and receiver group intervals, respectively. Both data sets were merged in the dip moveout processing step, pre-stack depth migration (PSDM) was applied to the data, and the final processed merged seismic volume was delivered in 2006. This volume displays the most coherent seismic image along with the CEF, better resolution of small faults, and a more accurate velocity model (Van der Veeken and Lutz, 2008). In this study, the final PSDM processed merged seismic volume is used for seismic interpretation and seismic attribute analysis. Shell Geoscience Solutions processed the seismic data in the Netherlands. Eight vertical wells (Maui-1, Maui-2, Maui-3, Maui-5, Maui-6, Maui-7, MB-Z11, and Rahi-1) are located within the 3D seismic survey area. Each well has conventional well logs (gamma-ray [GR], density, neutron, sonic, etc.) and completion reports. Some core reports and core photos are available for the Maui-5, Maui-6, and Maui-7 wells.
The final processed PSDM merged seismic volume covers a total area of 1000 km2 and was recorded to 3.5 s two-way traveltime (TWT) at a 4 ms sample interval. The data were cropped laterally to cover an area of 410 km2 and vertically to cover a time range of 1.8–2.8 s TWT. We cropped the original seismic volume to a smaller one to focus on areas with diverse fault styles within the C Sand reservoir and to reduce the computational time and hard-drive space. The seismic data were spectrally balanced during the 1991 and 2001 seismic volume merge (Van der Veeken and Lutz, 2008). How- ever, the data are contaminated by random noise and acquisition artefacts. Thus, we followed Chopra and Marfurt’s (2007) approach and applied principal component structure-oriented filtering (PC-SOF) along with the structural dip. This method enhances the signal that is aligned with the estimated dip and reduces the noise that crosses in other directions.
The seismic-to-well tie for the MB-Z11 well shows the C Shale, C Sand, and D Sand seismic events (Figure 2a). We map the top of the C Shale, C Sand, and D Sand as zero-crossings on quadrature component data, which is calculated from the preconditioned seismic volume. We consider the time structure map, horizon amplitude ex- traction, and high spatial residual maps, which is calculated by subtracting a high-precision manually picked seismic horizon from a smoothed seismic horizon (Brown, 2011), as horizon-based attributes. We calculate 30 typical geometric attributes available for a seismic interpreter that include volumetric dip, coherence, amplitude gradients, curvature, curvedness, shape in- dex, reflector rotation, Euler curvature, aberrancy, and GLCMs attributes. We then use a workflow to analyze and select the ones that are most useful for a human interpreter and those that work best for SOM analysis (Figure 3).
Seismic attribute analysis
We begin by analyzing 30 widely used horizon-based and geometric attributes for fault delineation (Figure 4). A description of each attribute is shown in Table S1 (supplementary information can be accessed through the following link: S1). As with most subsurface studies, there are limited well penetrations to independently verify the analyzed faults within the study area; as such, “real” faults will be evaluated from our analyses by verifying cross-sections for systematically offset reflectors (i.e., fault throws) and by assessing their map view likeness to known fault and fracture network patterns from geologic studies (Peacock et al., 2016) and structural analogue models (Wu et al., 2015). These interpreted subtle faults and fractures should be calibrated to well data (e.g., image logs or cores) when it becomes available. The goal of this analysis is to compare the relative strengths, weaknesses, and redundancies of the individual analyses to inform our later multiattribute analyses. The high spatial residual map and amplitude extraction (Figure 4a and 4b) along the C Sand reservoir shows the WF on the west side of the Maui B region. Some of the faults have clear seismic offsets (the white arrows), but in some areas, the attributes seem to be sensitive to noise (the red arrow). We evaluate the other 28 geometric attributes by using horizon slices extracted along the top of the C Sand reservoir. The inline dip, crossline dip, and dip magnitude attributes (Figure 4c–4e) show similar edges. Thus, the dip magnitude can be used as an edge detection attribute and the inline and crossline dips can be used to generate other attributes such as coherence and curvature. The multispectral energy
ratio similarity (Figure 4f) shows cleaner and sharper edges associated with the small faults in the central and northern parts of the Maui field compared to the broadband (Figure 4g) and corendered spectral energy ratio similarity (Figure 4h) slices. Inline and crossline root-mean-square (RMS) amplitude weighted amplitude gradients (Figure 4i and 4j) show similar edges. How- ever, the crossline RMS amplitude weighted amplitude gradient is affected by acquisition artefacts associated with the north–south-shooting orientation. As a result, these attributes do not provide additional geologic information compared to other attributes.
The long- and short-wavelength most-positive and most-negative principal curvature attributes are useful for mapping faults (Figure 4k–4m). However, the most- negative principal curvature attributes show better lineaments compared to those depicted by the most positive principal curvature. The amplitude curvature (Figure 4n) shows acquisition artefacts that mask the geology. These north–south-trending artefacts may result from the slight differences in the acquisition geometries between the 1991 and 2002 seismic surveys. As a result, this attribute was not as useful for fault detection in this study; however, it can be used to check the seismic data quality. The most-positive and most-negative principal curvatures can be combined to calculate the curvedness attribute, which describes the total de- formation and the shape index that describes the shape of the seismic reflector (Marfurt, 2018) (Figure 4o and 4p). Reflector rotation (Figure 4q) is one of the curvature attributes that shows the rotation of the faulted blocks at the faults’ tips. Euler curvature attributes (Figure 4r and 4s) show similar fault trends as seen in other attributes; thus, we consider these attributes redundant in this study.
The total aberrancy magnitude (Figure 4t) shows edges similar to those detected by the maximum aberrancy magnitude (Figure 4u). The intermediate and minimum aberrancy magnitudes (Figure 4v and 4w) do not provide useful information. As a result, we use the total aberrancy in further analyses and neglect the rest of the aberrancy attributes. Visual evaluation of the individual cyclic color-based attributes such as the dip azimuth (Figure 4x), the strike of curvature (Figure 4y), and total aberrancy azimuth (Figure 4z) suggests that these are not useful. Later, we will see the value of these attributes when they are corendered with other attributes. GLCM attributes (Figure 4aa–4dd) show the same faults seen in other attributes and do not add much additional geologic information on their own. However, the GLCM attributes work well with classification machine learning algorithms, as we will show later. In summary, some individual attributes (e.g., dip magnitude, multispectral energy ratio similarity, broadband energy ratio similarity, and others) seem more useful for detecting small faults whereas other attributes seem re- redundant or have limited value. In the next section, we will show that combining the structural information from multiple seismic attribute volumes can further enhance the subtle fault detection.
Interactive seismic attribute interpretation
Based on our analyses (Figure 4), we select the individual attributes of dip magnitude, dip azimuth, multispectral energy ratio similarity, broadband energy ratio similarity, long- and short- wavelength most-positive and most-negative principal curvatures, curvedness, shape index, reflector rotation, and total aberrancy attributes for further analysis. Here, we attempt to combine attributes to extract more information from the data. We combine the attributes in two ways: First, we use color blending of multiple seismic attributes to corender the geologic information captured from a maximum of four seismic attribute volumes; second, we use unsupervised machine learning algorithms, which enable integrating the structural content from as many attributes as needed.
Attribute selection for a human interpreter — Multiattribute display
Figure 5 shows a comparison between multiattribute maps, and a more detailed comparison of the northeast map corner is available in Table S2 (supplementary information can be accessed through the following link: S2). Coren- dering the dip magnitude and dip azimuth (Figure 5a) shows the edges of the main faults affecting the reservoir and the changes in dip of the seismic events that result from tectonic deformation. Corendering the long- and short-wavelength principal structure curvatures with multispectral energy ratio similarity (Figure 5b and 5c) shows the main faults affecting the area. The short-wavelength curvature attributes show a localized intense zone of fractures close to the main faults (Figure 5c). These fractures are below the seismic resolution with 19 m tuning thickness and cannot be detected by the long-wavelength curvature attributes (Figure 5b). The most positive and most- negative principal curvature attributes show positive and negative curvature anomalies that bracket the multispectral energy ratio similarity anomalies of the reverse fault WF (the red arrow) on the west side of Maui B (Figure 5b).
The lineaments seen in the most positive and most negative curvature attributes occur in which there are no observed multispectral energy ratio similarity
anomalies (the green arrows in Figure 5b). This indicates that there are small faults with no clear offsets in this part of the field. In the central and eastern parts of the Maui A region, the multispectral energy ratio similarity attribute shows segmented anomalies associated with Ihi faults (IF) (the cyan arrow in Figure 5b) and small faults shown by the white arrow. However, the most-positive and most-negative principal curvatures and the corendered dip azimuth and dip magnitude show continuous fault alignments of the IF and intersecting faults in the eastern part of Maui A region (Figure 5a and 5b). In the far northern part of the Maui field, the multispectral energy ratio similarity attribute shows a set of en echelon faults (the purple arrow in Figure 5b). In contrast, the curvature and dip attributes show two intersecting continuous faults.
The most-positive and most-negative principal curvatures modulated with their strikes (Figure 5d and 5e) show lineaments that can be interpreted as small faults or fractures that are below the seismic resolution. These are not well depicted by similarity or dip attributes. Corendered reflector rotation with multispectral energy ratio similarity and curvedness (Figure 5f) shows the faults along with the relative rotation of the faulted blocks about the faults’ tips in the central and northern parts of the field. This display indicates that the northern part of the area is affected by a set of en echelon normal faults, not two faults intersecting each other. The shape index modulated with curved- ness (Figure 5g) shows the total deformation of the seismic events and resolves not only the main faults affecting in the reservoir but it also shows lineaments parallel to the main fault trend that could be interpreted as fractures that are below the seismic resolution.
Total aberrancy modulated with multispectral energy ratio similarity (Figure 5h) shows sharp aberrancy anomalies at the same spatial position as the multispectral energy ratio similarity anomalies, illuminating the main faults affecting in the area. Aberrancy also shows continuous anomalies that clearly depict the IF and the faults in Maui south and east.
We now discuss vertical slices extracted from the seismic attribute volumes. Vertical slices A-A′ (the green line in Figure 5b) cross the far northern part of Maui and are extracted from the preconditioned seismic amplitude (Figure 6a); the corendered preconditioned seismic amplitude, broadband energy ratio similarity, dip magnitude and dip azimuth volumes (Figure 6b); the corendered preconditioned seismic amplitude, multispectral energy ratio similarity, dip magnitude and dip azimuth volumes (Figure 6c); and the corendered preconditioned seismic amplitude, long- wavelength most-positive and most-negative principal curvature attributes (Figure 6d). Faults that show a clear offset on seismic amplitude data can be detected by energy ratio similarity attributes (Figure 6b and 6c). However, the multispectral energy ratio similarity shows sharper and more continuous faults (the cyan ar- row in Figure 6a) compared to those depicted by broadband energy ratio similarity attributes. We notice a fault with clear offset (the purple arrow in Figure 6a) and deformation of the seismic event on the hanging wall as evidenced by changes in the dip azimuth and bent reflectors. This fault exhibits a similarity anomaly at the inferred fault position, whereas it shows dip and negative curvature anomalies on the hanging wall of the fault.
The green arrow in Figure 6a shows one of the en echelon faults interpreted on the multispectral energy ratio data (Figure 6b). There is a clear offset but no change in the dip of the seismic reflectors on either side of the fault and a little upward bend in the footwall re- flectors. This gives rise to a similarity anomaly at the inferred fault position and a positive curvature anomaly on the footwall of the fault (Figure 6d). Therefore, we prefer to use the multispectral energy ratio similarity attribute to map these faults because it provides sharper and more continuous similarity anomalies at the fault planes (Figure 6c) compared to those depicted by the curvature and dip attributes. We also observe some low-similarity anomalies (the orange arrow in Figure 6c) across the line that might indicate channelized features.
Vertical slice B-B′ (the yellow line in Figure 5b) crosses the IF and is extracted from the preconditioned seismic volume (Figure 7a); the corendered preconditioned seismic amplitude, broadband energy ratio similarity, dip magnitude and dip azimuth (Figure 7b); the corendered preconditioned seismic amplitude, multispectral energy ratio similarity, dip magnitude and dip azimuth (Figure 7c); and the corendered preconditioned seismic, long-wavelength most-positive and most- negative principal curvature attributes (Figure 7d).
In Figure 7a, the fault denoted by the purple arrow has a clear offset and no change in the dip of the seismic re- flectors on either side of the fault. However, these reflectors appear slightly folded near the fault (the red dashed line in Figure 7a). This fault shows clear energy ratio similarity anomalies (Figure 7b and 7c). The most-positive and most-negative curvature attributes show anomalies on the footwall and hanging wall of this normal fault (Figure 7d). One fault (the green arrow of Figure 7a) of the IF does not display clear offsets, but there are changes in the dip of the faulted blocks on each side of the fault that show a dip anomaly. In addition, the seismic reflectors on each side of this fault are bent and dis-play positive and negative anomalies on both sides of the fault (Figure 7d). As a result, this fault is not well
detected by the energy ratio similarity attributes. The other IF (the blue arrow) displays a clear offset and changes in the dip of the seismic reflectors on each side of the fault in addition to bent reflectors along the hanging wall of the fault. This fault exhibits clear energy ratio similarity anomalies (Figure 7b and 7c), a small dip anomaly, and small negative curvature anomalies on the hanging wall of this fault. Thus, horizon slices extracted from the dip and curvature attributes show more continuous edges of IF (Figure 5a–5c) compared to the energy ratio similarity images. The arrows in Figure 7c show low-energy-ratio similarity features that might indicate deformation within the shales overlying the C Sand reservoir (the cyan arrow) and small channelized features (the orange arrow).
In this section, we have shown how the visual analysis of geometric attributes can aid in detecting faults. We have discussed the eight most useful combinations of 16 different geometric attributes used to map faults and fractures in the Maui field. Given the large data volumes inherent in most 3D seismic surveys, delineation of small faults and fractures for reservoir development in a timely manner can be challenging. Corendering multiple attributes are useful but are limited to four attributes using RGB, cyan-magenta-yellow, or hue-lightness- saturation. As a result, we seek methods to combine large numbers of attributes that are also less reliant on visual interpretation. Below, we explore other ways to combine attributes using unsupervised machine learning.
Attribute selection for unsupervised machine learning clustering algorithm (SOM) — Human experience and principal component analysis
Computer-aided clustering techniques are useful for integrating the geologic content of more than three seismic attributes (Zhao et al., 2015; Marfurt, 2018). Unsupervised machine learning algorithms such as SOM apply a nonlinear clustering to multiple seismic attribute volumes to understand how these attributes relate to each other and to reduce the dimensionality of the data (Kohonen, 1982). A SOM groups similar data and uses a 2D topology map to display these groups or clusters. The number of clusters or classes is controlled by a user-defined number of neurons (called prototype vectors in other publications) within the 2D map. These clusters gather the information from the seismic attribute volumes into one classification volume that can facilitate geologic interpretation (Coléou et al., 2003; Roy, 2013).
Selecting the right attributes to be used as inputs to a SOM influences the classification results (Barnes and Laughlin, 2002). PCA is a linear mathematical technique that aids in distilling a library of seismic attributes into a small set of attributes, which contain the largest variations within the data set (Roden et al., 2015). However, PCA does not provide any information about the spatial positions of anomalous zones within the seismic attribute volumes. In general, machine learning algorithms are more quantitative and can perform re- peated tasks with as many voxels as needed more quickly compared to humans. However, humans can recognize larger scale fault patterns (Wu et al., 2015; Peacock et al., 2016) and put them in a geologic context. For instance, a human interpreter can infer the presence of a fault by observing that the most-positive and most-negative curvature anomalies bracket a coherence anomaly. The computer algorithms cannot understand these spatial patterns. Therefore, not all of the attributes that human interpreters use for fault interpretation can be used as inputs to machine learn- ing algorithms. Machine learning is more ideal when attributes that highlight a feature of interest lie within the same seismic voxel (Infante-Paez and Marfurt, 2019). Based on our previous discussion, the dip magnitude, coherence, curvedness, total aberrancy magnitude, and GLCM attributes show anomalies at the inferred fault planes. Thus, these attributes are considered as candidates for SOM analysis.
The SOM workflow used in this study is shown in Figure 8. We use the multispectral energy ratio similarity, dip magnitude, curvedness, total aberrancy magnitude, GLCM homogeneity, and GLCM entropy as inputs to PCA. We feed these six attributes with their associated weights obtained from the first four PCA eigenvectors (i.e., principal components of the six attributes) to the SOM analysis. The PCA and SOM analyses are performed on the seismic samples (voxels) within the in- interval of interest, which is from the C Shale horizon to the D Sand horizon. The SOM classification uses an
8 × 8 topology. A horizon slice extracted from the SOM classification volume along the top of the C Sand reservoir (Figure 9a) clearly shows all the small faults in the reservoir. In the northern part of the field, the horizon slice shows the sharp edges of a northwest-southeast fault (the white arrow in Figure 9a; F1 in Figure 9b) that is cross-cut by a set of en echelon normal faults (F2 in Figure 9b). Imaging such fault patterns and understanding the fault mechanism is important for advanced structure and reservoir analyses. The horizon slice (Figure 9a) also shows clear and sharp edges for faults in the Maui east (the yellow arrows), Maui west (the blue arrow), and Maui south (the cyan arrow) areas. Mapping small faults and obtaining the sharpest fault edges is critical for reservoir development studies. In the Maui field, the SOM result shows the sharpest edges of the conjugate IFs (the black arrow) compared to other attributes (Figures 4 and 5). These faults comprise small segments that seem to touch each other at their tips (Figure 9c). Understanding this fault pattern enabled the operator (Shell in 2006) to drill MA-02A well on the west side of IF. The field was producing from the C Sand reservoir on the east side of IF for 29 years. In 2006, the MA-02A well (the black circle in Figure 9a) penetrated the same reservoir on the west side of IF and encountered gas with the original gas- water contact depth (i.e., unswept reservoir) (Telford and Murray, 2008).
Although many case studies have been published re- regarding the application of SOM to facies classification in different geologic settings (Roy et al., 2010; Matos et al., 2011; Zhao et al., 2015, 2018; Infante-Paez and Marfurt, 2019), there has been less focus on the application of SOM to fault detection (Dewett and Henza, 2016; Laudon et al., 2019). To demonstrate the importance of selecting the appropriate attributes for detecting small faults, we test the Laudon et al.’s (2019) workflow on the same data set (Figure 9). Following Laudon et al. (2019), we use the most positive principal curvature (K1), energy-ratio similarity, GLCM entropy, and GLCM homogeneity attributes as our SOM analysis inputs (right side of Figure 9). We combine geologic information from the most positive principal curvature attribute, which shows anomalous features at different spatial positions with respect to the inferred fault locations, with similarity and GLCM attributes that display anomalous features at the inferred fault planes. A horizon extraction along the top of the C Sand reservoir from the produced previous SOM (Figure 9d) shows the faults in the Maui east (the yellow arrow) and south
the cyan arrow) areas of the field are detected but appear blurred compared to the SOM slice obtained using principal components of the attributes (compare Figure 9a and 9d). Our simulation of Laudon et al.’s approach (Figure 9d) does not detect the faults that affect the far northern (the white arrow) and northwest parts of the field (the blue arrow), where the SOM result seems to be noisy. A comparison of Figure 9b and 9e shows a similar outcome. Our SOM approach more clearly shows the en echelon faults affecting Maui north (Figure 9b), whereas the results obtained using Laudon et al.’s approach look noisy (Figure 9e). In addition, the SOM classification volume (the black arrow in Figure 9d) failed to depict the segmented faults within the IFs. Comparing Figure 9c and 9f, we note that our method (Figure 9c) clearly shows sharper segmented faults within IF whereas the Loudon et al. approach shows blurred IFs (Figure 9f). Closer investigation of the SOM results and turning on the neurons related faults shows that the SOM from this study SOM (Figure 9g) produces cleaner and sharper faults relative to the simulated Laudon et al. (2019) SOM (Figure 9h). Therefore, our proposed workflow incorporates the principal components of six geometric seismic attributes (i.e., multispectral energy ratio similarity, dip magnitude, curvedness, total aberrancy magnitude, GLCM homogeneity, and GLCM entropy) may provide alternative seismic attributes for application to further SOM studies.
Small faults that are obscured on seismic amplitude data may form important seals or conduits that dramatically control hydrocarbon migration and trapping within a petroleum reservoir. They are thus important to delineate. We found that seismic attributes are useful for detecting subtle faults in the Maui field of New Zealand. However, interpreting individual attributes can be a tedious and time-consuming task. Multiple attribute analyses, such as corendering, PCA, and SOM algorithms, help to extract more insightful geologic information compared to individual volumetric attribute analyses. However, unlike human interpreters, who are capable of recognizing spatial patterns within data and assign geologic context, the SOM technique works best when the multiple seismic attribute volumes exhibit anomalous features within the same seismic voxel. From our investigation, we found incorporating appropriate seismic attributes that exhibit anomalous features at the same seismic voxel to PCA and SOM analyses aid in combining the geologic context from six seismic attributes into one enhanced classification volume. This clustered volume enables us to better detect subtle faults within the C Sand reservoir of Maui Field that were not easily revealed from conventional seismic interpretation techniques. In the future, given sufficient amounts of labelled data and selecting the appropriate attributes, we anticipate that supervised machine learning algorithms such as CNNs may be able to mimic the ability of a human interpreter to delineate small faults within hydrocarbon reservoirs.
We would like to thank the Ministry of Business, Innovation, and Employment (MBIE), New Zealand, for providing a data set containing 3D seismic, well logs, and associated reports required for this study. The generosity of Geophysical Insights and Schlumberger is greatly appreciated for supporting us with Paradise 3.2.3 and Petrel G&G software v.2019, respectively. We express our gratitude to K. Marfurt for useful discussions and his support of the Attribute Assisted Seismic Processing and Interpretation software. Discussions with S. Stanley, B. Taylor, and S. Fuston were useful and much appreciated. We would like to thank D. Johnston for commenting on an earlier version of the manuscript. We thank reviewers P. Lanzarone and D. Dewett, one anonymous reviewer, and the associate editor T. Zhao for the constructive and helpful comments. H. Hussein and J. Wu acknowledge the computer work- station equipment purchased by the Governor’s University Research Initiative provided by the State of Texas and University of Houston. Publication fees were supported by the UH Center for Tectonics and Tomography led by John Suppe.
Data and materials availability
We requested the data from New Zealand via the following website and they gave us permission to use this data for research: https://data.nzpam.govt.nz/GOLD/system/mainframe.asp
Al-Dossary, S., and K. J. Marfurt, 2006, 3D volumetric multi- spectral estimates of reflector curvature and rotation: Geophysics, 71, no. 5, P41–P51, doi: 10.1190/1.2242449. Aqrawi, A. A., and T. H. Boe, 2011, Improved fault segmentation using a dip guided and modified 3D Sobel filter: 81st Annual International Meeting, SEG, Expanded Abstracts, 999–1003, doi: 10.1190/1.3628241.
Barnes, A. E., 2006, A filter to improve seismic discontinuity data for fault interpretation: Geophysics, 71, no. 3, P1–P4, doi: 10.1190/1.2195988.
Barnes, A. E., and K. J. Laughlin, 2002, Investigation of methods for unsupervised classification of seismic data: 72nd Annual International Meeting, SEG, Expanded Abstracts, 2221–2224, doi: 10.1190/1.1817152.
Baur, J. R., P. R. King, T. Stern, and B. Leitner, 2010, Development and seismic geomorphology of a Miocene slope channel mega system, Offshore Taranaki Basin, New Zealand: Seismic imaging of depositional and geomorphic systems: SEPM Society for Sedimentary Geology, 30, 618–649, doi: 10.5724/gcs.10.30.0618.
Brown, A. A., 2011, Interpretation of three-dimensional seismic data: AAPG Memoir 42, SEG Investigations in Geophysics.
Chopra, S., and K. J. Marfurt, 2007, Seismic attributes for prospect identification and reservoir characterization: SEG Geophysical Developments Series.
Coléou, T., M. Poupon, and K. Azbel, 2003, Unsupervised seismic facies classification: A review and comparison of techniques and implementation: The Leading Edge, 22, 942–953, doi: 10.1190/1.1623635.
Dewett, D. T., and A. A. Henza, 2016, Spectral similarity fault enhancement: Interpretation, 4, no. 1, SB149– SB159, doi: 10.1190/INT-2015-0114.1.
Di, H., and D. Gao, 2016, Improved estimates of seismic curvature and flexure based on 3D surface rotation in the presence of structure dip: Geophysics, 81, no. 2, IM13–IM23, doi: 10.1190/geo2015-0258.1.
Franzel, M., and S. Back, 2019, Three-dimensional seismic sedimentology and stratigraphic architecture of pro- grading clinoforms, central Taranaki Basin, New Zealand: International Journal of Earth Sciences, 108, no. 2, 475–496, doi: 10.1007/s00531-018-1663-1.
Gao, D., 2011, Latest developments in seismic texture analysis for subsurface structure, facies, and reservoir characterization: A review: Geophysics, 76, no. 2, W1– W13, doi: 10.1190/1.3553479.
Gao, D., 2013, Integrating 3D seismic curvature and curvature gradient attributes for fracture characterization: Methodologies and interpretational implications: Geophysics, 78, no. 2, O21–O31, doi: 10.1190/geo2012-0190.1. Gersztenkorn, A., and K. J. Marfurt, 1999, Eigenstructure- based coherence computations as an aid to 3-D structural and stratigraphic mapping: Geophysics, 64, 1468–1479, doi: 10.1190/1.1444651.
Haque, A. E., M. A. Islam, and M. Ragab Shalaby, 2016, Structural modelling of the Maui Gas Field, Taranaki basin, New Zealand: Petroleum Exploration and Development, 43, 965–975, doi: 10.1016/S1876-3804(16)30114-8.
Henderson, J., S. J. Purves, and G. Fisher, 2008, Delineation of geological elements from RGB color blending of seismic attribute volumes: The Leading Edge, 27, 342–350, doi: 10.1190/1.2896625.
Holt, W. E., and T. A. Stern, 1994, Subduction, platform subsidence, and foreland thrust loading: The late Tertiary development of Taranaki Basin, New Zealand: Tectonics, 13, 1068–1092, doi: 10.1029/94TC00454.
Infante-Paez, L., and K. J. Marfurt, 2019, Using machine learning as an aid to seismic geomorphology, which attributes are the best input? Interpretation, 7, no. 3, SE1–SE18, doi: 10.1190/INT-2018-0096.1.
King, P. R., and G. P. Thrasher, 1996, Cretaceous-Cenozoic geology and petroleum systems of the Taranaki basin, New Zealand: Institute of Geological and Nuclear Sciences (New Zealand), Monograph 13.
Kohonen, T., 1982, Self-organized formation of topologically correct feature maps: Biological Cybernetics, 43, 59–69, doi: 10.1007/BF00337288.
Laudon, C., S. Stanley, and P. Santogrossi, 2019, Machine learning applied to 3D seismic data from the Denver- Julesburg basin improves stratigraphic resolution in the Niobrara: Unconventional Resources Technology Conference, SEG, Global Meeting Abstracts, 1–17, doi: 10.15530/URTEC-2019-337.
Li, F., and W. Lu, 2014, Coherence attribute at different spectral scales: Interpretation, 2, no. 1, SA99–SA106, doi: 10.1190/INT-2013-0089.1.
Libak, A., B. Alaei, and A. Torabi, 2017, Fault visualization and identification in fault seismic attribute volumes: Implications for fault geometric characterization: Interpretation, 5, no. 2, B1–B16, doi: 10.1190/INT-2016-0152.1.
Lowell, J., and V. Erdogan, 2019, Combining artificial intelligence with human reasoning for seismic interpretation James: 89th Annual International Meeting, SEG, Expanded Abstracts, 2353–2357, doi: 10.1190/segam2019- 3216650.1.
Luo, Y., W. G. Higgs, and W. S. Kowalik, 1996, Edge detection and stratigraphic analysis using 3D seismic data: 66th Annual International Meeting, SEG, Expanded Abstracts, 324–327, doi: 10.1190/1.1826632.
Mai, H. T., K. J. Marfurt, and S. Chávez-Pérez, 2009, Coherence and volumetric curvatures and their spatial relationship to faults and folds, an example from Chicontepec basin, Mexico: 79th Annual International Meeting, SEG, Expanded Abstracts, 1063–1067, doi: 10.1190/ 1.3255033.
Marfurt, K. J., 2018, Seismic attributes as the framework for data integration throughout the oilfield life cycle: SEG Distinguished Instructor Short Course.
Matos, M. C. De, M. M. Yenugu, S. M. Angelo, and K. J. Mar-furt, 2011, Integrated seismic texture segmentation and cluster analysis applied to channel delineation and chert reservoir characterization: Geophysics, 76, no. 5, P11–P21, doi: 10.1190/geo2010-0150.1.
Nicol, A., J. Walsh, K. Berryman, and S. Nodder, 2005, Growth of a normal fault by the accumulation of slip over millions of years: Journal of Structural Geology, 27, 327–342, doi: 10.1016/j.jsg.2004.09.002.
Nodder, S. D., 1993, Neotectonics of the offshore Cape Egmont fault zone, Taranaki basin, New Zealand: New Zealand Journal of Geology and Geophysics, 36, 167– 184, doi: 10.1080/00288306.1993.9514566.
Pannett, S., S. Slager, G. Stone, and S. Dekker, 2004, Constraining a complex gas-water dynamic model using 4D seismic: Annual Technical Conference and Exhibition, SPE, Extended Abstracts, 209–225.
Partyka, G., J. Gridley, and J. Lopez, 1999, Interpretational applications of spectral decomposition in reservoir characterization: The Leading Edge, 18, 353–360, doi: 10.1190/1.1438295.
Peacock, D. C. P., C. W. Nixon, A. Rotevatn, D. J. Sander- son, and L. F. Zuluaga, 2016, Glossary of the fault and other fracture networks: Journal of Structural Geology, 92, 12–29, doi: 10.1016/j.jsg.2016.09.008.
Pedersen, S. I., T. Randen, L. Sønneland, and Ø. Steen, 2002, Automatic fault extraction using artificial ants: 72nd Annual International Meeting, SEG, Expanded Abstracts, 512–515, doi: 10.1190/1.1817297.
Pedersen, S. I., T. Skov, A. Hetlelid, P. Fayemendy, T. Ran- den, and L. Sønneland, 2003, New paradigm of fault interpretation: 73rd Annual International Meeting, SEG, Expanded Abstracts, 350–353, doi: 10.1190/1.1817918.
Pochet, A., P. H. B. Diniz, H. Lopes, and M. Gattass, 2019, Seismic fault detection using convolutional neural networks trained on synthetic post stacked amplitude maps: IEEE Geoscience and Remote Sensing Letters, 16, 352–356, doi: 10.1109/LGRS.2018.2875836.
Qi, X., and K. Marfurt, 2018, Volumetric aberrancy to map subtle faults and flexures: Interpretation, 6, no. 2, T349– T365, doi: 10.1190/INT-2017-0114.1.
Randen, T., E. Monsen, C. Signer, A. Abrahamsen, J. O. Hansen, T. Sæter, J. Schlaf, and L. Sønneland, 2000, Three-dimensional texture attributes for seismic data analysis: 70th Annual International Meeting, SEG, Ex-wana: Two-phase Cretaceous-Paleocene rifting in the Taranaki Basin, New Zealand: International Conference & Exhibition, AAPG.
Telford, C., and D. Murray, 2008, MA-02A well completion report: NZP&M, Ministry of Business, Innovation & Employment (MBIE), New Zealand Unpublished Petroleum Report PR4491.
Van der Veeken, J. P. M., and J. B. M. Lutz, 2008, Maui PreSDM 2006 processing: NZP&M, Ministry of Business, Innovation & Employment (MBIE), New Zealand Unpublished Petroleum Report PR4491.
Wu, J. E., K. McClay, and E. Frankowicz, 2015, Niger Delta gravity-driven deformation above the relict Chain and Charcot oceanic fracture zones, Gulf of Guinea: Insights from analogue models: Marine and Petroleum Geology, 65, 43–62, doi: 10.1016/j.marpetgeo.2015.03.008.
Wu, X., Y. Shi, S. Fomel, L. Liang, Q. Zhang, and A. Z. Yu- sifov, 2019, FaultNet3D: Predicting fault probabilities, strikes, and dips with a single convolutional neural network: IEEE Transactions on Geoscience and Remote Sensing, 57, 9138–9155, doi: 10.1109/TGRS.2019.2925003.
Zhao, T., 2019, 3D convolutional neural networks for efficient fault detection and orientation estimation: 89th Annual International Meeting, SEG, Expanded Abstracts, 2418–2422, doi: 10.1190/segam2019-3216307.1.
Zhao, T., V. Jayaram, A. Roy, and K. J. Marfurt, 2015, A comparison of classification techniques for seismic facies recognition: Interpretation, 3, no. 4, SAE29–SAE58, doi: 10.1190/INT-2015-0044.1.
Zhao, T., F. Li, and K. J. Marfurt, 2018, Seismic attribute selection for unsupervised seismic facies analysis using user-guided data-adaptive weights: Geophysics, 83, no. 2, O31–O44, doi: 10.1190/geo2017-0192.1.
Most Popular Papers
Systematic Workflow for Reservoir Characterization in Northwestern Colombia using Multi-attribute Classification