A local binary patterns/variance operator based on guided filtering for seismic fault detection

Aiming at suppressing noise interference, improving the fault detection ability of seismic data, fully excavating the effective information in seismic data, and further improving the accuracy of fault detection, this study proposes a seismic fault detection method that combines the local binary pattern/variance (LBP/VAR) operator with guided filtering. The proposed method combines the advantages of LBP/VAR and guided filtering to remove noise from seismic data, and can simultaneously smooth the data and preserve linear features. When compared with several existing methods (coherent operator, LBP/VAR operator, LBP/VAR operator based on median filtering, and Canny operator based on guided filtering), the proposed method exhibits a better SNR, a better ability to identify small faults, and robustness to noise. This novel algorithm can control the balance between noise attenuation and effective signal preservation as well as effectively detect faults in seismic data. Therefore, the proposed method effectively improves the fault identification accuracy, facilitates the gas-bearing analysis of the structure, provides guidance for the actual well location deployment of the project, and has important practical significance for oil and gas exploration and development.


Introduction
Fault detection is a basic and important aspect of seismic data interpretation, the reliability of which is greatly affected by the signal-to-noise ratio (SNR). Due to the influence of various complex factors, such as the vibration of the external environment, the interference of industrial alternating currents, and complex shallow surfaces, seismic data inevitably contain random noise, which results in a relatively low SNR [1][2][3]. Thus, reducing random noise during seismic data processing, while preserving complex structural information has important practical significance. Many filtering methods are used to remove noise and improve the SNR, such as average filtering, median filtering [4], alpha mean filtering, and PC filtering [5]. However, these filtering methods also blur transverse discontinuities. Dossary et al. [6] reviewed filtering methods that can simultaneously smooth seismic data and preserve linear features and indicated the need to develop new seismic signal denoising technology.
Guided filtering has exhibited good image processing effects in other fields [7][8][9]. This algorithm employs a local linear model as an anisotropic filter [10], which can effectively remove seismic noise and preserve edge characteristics. Li et al. [11] added a first-order edge perception constraint and boundary perception weight to the filter, resulting in a weighted guided filtering algorithm that avoids the halo effect as much as possible. Dai et al. | https://doi.org/10.1007/s42452-021-04866-0 [12] introduced an additional measure of spatial similarity, developed full-connection guided filtering, and improved the self-adaptability of the local area of the image. Moreover, Kou et al. [13] achieved good results in image detail enhancement and significance detection using gradient domain-guided filtering.
According to the fault characteristics of seismic data, the local binary pattern (LBP) can be used for data processing and detection. The LBP was proposed by Ojala [14] as an algorithm for describing image texture features. Ojala et al. [15] first extended the basic LBP operator to any circular neighborhood from the perspective of a fixed topological structure. Others have improved the LBP by reducing the impact of noise and extending the greyscale between pixels by integrating the local binary mode into the local multi-valued mode [16]. Recently, many scholars have made improvements to the LBP and proposed various deformation LBP operators applicable to different fields [17][18][19]. Guo et al. [20] extended the completed LBP (CLBP) operator from the perspective of the coding mode. Tadi-Bani et al. [21] propose a new content based image retrieval approach using combination of color and texture information in spatial and transform domains, which provides higher precision than many existing methods. Heikkila et al. [22] proposed a centrosymmetric LBP (CS-LBP), which effectively reduces the dimensions of the basic LBP operator to 16 dimensions. Fekri-Ershad [23] proposed titled multi threshold uniform-based local ternary patterns with notation MTULTP, which is a skillful combination of LTP and MLBP with novelty in feature extracting and local pattern selecting. Furthermore, to integrate local spatial structure and contrast information, Lei et al. [24] proposed the LBP/VAR operator, which can detect faults and perform threshold optimization processing, resulting in significant advantages in detecting seismic faults.
In this study, the advantages of guided filtering and the LBP/VAR operator are combined to the extract local texture information of images for seismic fault detection. This method is robust to light and rotation and can compare textures between images. Experiments are conducted using model data, seismic data, and fault pictures, which confirm the ability of the proposed operator to effectively describe faults.

Guided filtering
The guided image filter is a local linear image filter [7,8] with good edge preservation characteristics that requires a guiding graph. This guiding image can be a single image or an input image; if an input image is used, the guided filtering will preserve the edges, which can be used for image reconstruction [25].
Guided filtering is used to process 2D seismic images in the time domain; it is assumed that the input 2D seismic image is p i , the guiding image is I i , and the output image is q i . I i and p i can be the same or different seismic images. The filter linear transformation relationship between I i and q i is where a k , b k are the coefficients of the linear function when the center of the window is located at k, r (k) is the center of the pixel k, and r is the radius of r (k) . We then find the gradient for both sides of equation: where a is a constant, which means that the edge information of q i and I i has a linear relationship. In order to determine the linear coefficient, the difference between q i and p i should be minimized. The linear coefficient can be obtained by minimizing the cost function, which is defined as: where ε is the regularization factor to adjust coefficient a k and prevent it from becoming too large. The coefficients a k and b k are obtained by solving the minimum extremum of Eq. 3; for details refer to Zhang [26].

LBP/VAR operator
The LBP is an algorithm proposed by Ojala et al. [14] to assist in describing the local contrast of images and to extract the local texture information of images. The traditional LBP has many limitations. In recent years, many scholars have improved the LBP in concrete applications and put forward various deformation LBP operators.
The rotation invariant unified local binary pattern (LBP riu2 ) is in a pixel region with a radius of R and P is the number of neighboring points. Fekri-Ershad studied circular symmetric neighbors for various radius (R) and number of neighbors (P) [27]. The gray value of the center point (g c ) from the gray value of the surrounding pixels of the circularly symmetric neighborhood g p . LBP riu2 of the center point in the local region is defined as The values of s(x) , g p (x, y),and U(LBP P,R ) are computed as is the sign function, and U(LBP P,R ) determines the number of jumps of binary encoding (bit 0/1 change) and sets the local binary pattern of U < 2 as the unified pattern. The LBP riu2 operator is a gray-scale rotation invariant operator, which can well describe the spatial structure of images and omits the contrast between pixels.
In fault detection using seismic data, it is important to distinguish the contrast of data points (pixels). Therefore, to retain the contrast between pixels, VAR is defined as: where µ is the expectation of the gray value of the sampling points in the neighborhood, and VAR can provide rich texture comparison information. The larger the VAR value is, the greater the contribution to the region's differentiation will be, and therefore, the coding weight corresponding to the region will be large.
The LBP/VAR operator is defined as: where: The LBP/VAR operator integrates the local spatial structure and contrast information, and has the advantages of both. It is not only robust to illumination and rotation, but also includes the texture contrast between images. Fault detection experiments using this operator have a good fault detection rate.

LBP/VAR Operator based on guided filtering
In this study, seismic data were optimized and de-noised by guided filtering, and then further enhanced images. Then, LBP/VAR operator technology was used to extract local texture information to identify seismic faults. A seismic fault detection method based on guided filtering of the LBP/VAR operator is proposed in this paper.
To further illustrate the LBP/VAR operator based seismic fault detection method, we demonstrate the selection of radius r and regularization parameter ε in the guided filtering and discuss the algorithm steps in detail.

Parameter selection
The smoothing effect of guided filtering on an image is determined by the window radius, r, and regularization parameter, ε. To calculate the matrix mean square error (MSE) of the theoretical model, a seismic section is adopted in this study. According to Fig. 1, the smaller the values of r and ε, the closer the mean square deviation to the original data; thus, fewer details of image are filtered out, resulting in higher similarity between the filtered data and the original data. On the other hand, filtering out more detail of image in a progressively smoother image. According to Fig. 1a, ε exhibits only a small change above 0.16; this change becomes increasingly small once r is larger than 5 (Fig. 1b). In this study, ε = 0.01; thus, when r = 1, the dominant effect is denoising, and when r = 16, the dominant effect is smoothing.

Algorithm steps
In this study, the input image is used as the guide graph for conducting self-guided image filtering to maintain the edge characteristics. The different steps of the algorithm are as follows: 1. Conduct self-guided image filtering on the input seismic profile, using r = 1 and ε = 0.01 for the dominant denoising function. The result of guided filtering, q 1 , loses few texture details. 2. Take q 1 as the input diagram of guided filtering and conduct self-guided image filtering. At this time, r = 16 and ε = 0.01; thus, the filtered image loses more details, and the filtering result q 2 is obtained. The filtered details are q 1 -q 2 . 3. Obtain the enhanced image. We use q 3 = t × (q 1q 2 ) + q 1 , t is the detail enhancement coefficient. Experiments show that when t = 3, the enhanced image has great advantages in computational efficiency and fault recognition effect in this paper. 4. Grayscale enhanced image. The amplitude value of the guided filtering image is converted into gray 256 levels. That is, gray level 0 corresponds to the lowest amplitudes and gray level 255 corresponds to the highest amplitudes on the guided filtering image. 5. According to Eq. (4), the pixel value of each 3 × 3 region in the slice of seismic data is calculated, and a new value is obtained, until each pixel in the image is LBP processed, and the local binarization image is obtained. 6. In image processing, small gray areas will produce a certain degree of interference. The LBP threshold is optimized according to the LBP histogram of the original seismic image. Values below the threshold are denoted as 0. Then, we refer to the LBP/VAR histogram to optimize the LBP/VAR threshold. In this paper, the upper and lower thresholds are set to 38 and 235, respectively. The gray value of pixels less than the lower threshold value of 38 is set to 0 (black), while the gray value of pixels greater than the upper threshold value of 235 is set to 255 (white). 7. According to Eqs. (4)-(8), the correlation features of LBP/VAR are further extracted. 8. The variable parameters are adjusted throughout the process by trial calculation, so that the results have better adaptability, stability, and effectiveness. The calculation flowchart is shown in Fig. 2.

Results and discussion
To verify the effectiveness of the LBP/VAR operator based on guided filtering for seismic fault detection, it was compared with different image edge-based fault detection operators, i.e., the coherent operator [2], the Canny operator based on guided filtering [5], the LBP/VAR operator [24], and the LBP/VAR operator based on median filtering. Four different tests were conducted to determine the effectiveness of the proposed model, i.e., fault detection using the general Marmousi model, a seismic section, a seismic slice, and actual fault images.

Theoretical fault model test
The Marmousi seismic model (shown in Fig. 3a) is a classic seismic forward geological model for evaluating seismic methods because it contains the characteristics of anisotropy, attenuation, and complex strata. It has widely been used in the field of geophysics and was selected to test the effectiveness of the LBP/VAR operator seismic fault detection method based on guided filtering.

Model test
We used split step Fourier forward modeling and migration to obtain the seismic profile of the Marmousi model, which had the CDP interval of 12.5 m, a sampling rate of 1 ms, and a 30-Hz ricker wavelet [28]. In this paper, we only used small areas with developed faults, and obvious fault responses were selected for processing. In this area, the number of seismic wave in-phase axes suddenly increases or disappears, and the seismic wave in-phase axes are staggered, bifurcated, and merged. There are phenomena such as the turning of in-phase axes under the theoretical model, which constitute a complex block fault system. The LBP/VAR, the median filtering-LBP/VAR, the guided filtering-Canny, and the guided filtering-LBP/VAR algorithm were compared, with the results shown in Fig. 3. The results show that the above operators can effectively detect faults, and the results are consistent with the theoretical model. Comparing Fig. 3b and c shows that the LBP/VAR operator exhibits better fault detection than the LBP/VAR operator based on median filtering. Although it improves the SNR, the LBP/VAR operator based on median filtering also blurs transverse discontinuities in some places, resulting in inferior detection results. Moreover, the detection results of the common Canny operator based on guided filtering (Fig. 3d) reveal that not all strata are depicted and blank areas remain, which is not conducive to the delineation of faults. Conversely, the LBP/VAR operator based on guided filtering detects faults slightly better than other methods (Fig. 3e). The detection results show fine stratigraphic characterization and easy identification of breakpoints.

Model added noise test
In this test, 20% and 30% Gaussian noise of the maximum effective amplitude was added to the Marmousi model (Fig. 3a). The test results of the LBP/VAR operator based on guided filtering were then analyzed (Fig. 4).
The results show that faults can be identified with different amounts of noise. Even when the noise reaches 20%, fault contours and formation boundaries can still be identified. For 20% noise, the LBP/VAR operator based on guided filtering exhibits slightly better fault detection than the LBP/VAR operator, as well as certain antiinterference ability.

Results of evaluation
We introduce structural similarity (SSIM) and peak signalto-noise ratio (PSNR) parameters to evaluate the image quality. SSIM is a good parameter for measuring image quality [29]. The simplest form of the equation is: Here, x and y are the averages of x and y, respectively, and σ x x and σ y y represent the standard deviations of x and y, respectively σ xy xy represents the covariance of x and y, whereas C 1 and C 2 are constants. The range of SSIM values is [0, 1]. The higher the SSIM, the better the quality of image.
PSNR is an objective standard used to measure image distortion or noise level. PSNR can be defined as follows: Here, Max is the maximum amplitude value of image, O(m,n) and MSE is the mean square error for the two images, which is defined as follows: Here, we used the original images without noise ,  O(m,n). The processed image, as the contrast images, C(m,n). Here, m is the horizontal pixel of the image and n is the vertical pixel of the image. Generally, PSNR is considered to be better if its value is greater than 30. The larger the PSNR between the two images, the more similar the two images are. We use Fig. 3b as the input image and Fig. 3c, Fig. 3d, Fig. 3e and Fig. 3f as the contrast image to calculate SSIM and PSNR, as shown in Table 1. Table 1 shows that the guided filtering-LBP/VAR method has the largest SSIM and the highest PSNR, and it is superior to the other methods regarding the quantitative evaluation of parameters PSNR and SSM.
Therefore, we chose this method when testing the other data.

Seismic profile detection
In the study area, faults are developed, including continuous strata, which constitute a complex block fault system. The detection of this section has high requirements for fault detection technology, which has a certain research value. The CDP interval is 12.5 m. Figure 5a shows the original seismic profile of inline1224, which has relatively welldeveloped faults and relatively continuous strata, forming a complex block fault system. Thus, this profile places higher requirements on fault detection technology. The coherence section obtained using landmarks is shown in Fig. 5b, which contains certain noise and other disorderly reflection phenomena. In Fig. 5b, the more continuous coaxial strong reflection region is depicted in blue in the coherent section (highly coherent region), and the strata typically reflect more continuous strata. The less coherent or incoherent regions (red and light blue) in the coherent section predominantly reflect the faults or cracks in the seismic section; however, some of the fault features are not obvious [25,30]. Figure 5c is obtained using the LBP/ VAR operator proposed in [24] to detect the fault. The overall section exhibits good quality and can depict the basic  Figure 5d is obtained by detecting faults with the LBP/VAR operator based on median filtering, which also filters out some effective information. The detection result of the canny operator with guided filtering is again not conducive to the delineation of faults (Fig. 5e). In contrast, the LBP/VAR operator based on guided filtering results in fine stratigraphic characterization and clearer small faults (Fig. 5f ).

Seismic slice detection
A planar graph of the study area after seismic interpretation is selected to detect and analyze the results. The results of each operator are shown in Fig. 6. The fault information can be better identified from the coherence attribute [25] in Fig. 4-11 (b). Again, the LBP/VAR operator based on the median filter preserves the structure, but blurs effective information (Fig. 6d), and the Canny operator based on guided filtering cannot effectively detect the study area plane (Fig. 6e), whereas the LBP/VAR operator based on guided filtering not only detects the faults detected by conventional operators, but also has certain advantages in the identification of small faults, which is helpful in analyzing the cutting relationship, extension range, and fault trend between faults in a plane (Fig. 6f ).

Actual geological section test
A geological section with clear fault characteristics (taken from https:// www. veer. com) was selected for analysis. This image is the Wolf Hill sandstone formation in Arizona, USA. Petrified Forest Member is a package of siliceous conglomeratic sandstones and interbedded mudstone called the Sonsela sandstone bed (Fig. 7a). The type area of the Sonsela sandstone bed is along the east flank of the Defiance Uplift north of Petrified Forest, where the unit is 120-200 feet thick and consists of two conglomeratic sandstone beds separated by siltstone. Since the fault features of the map are obvious, we used this image to verify the proposed method. The detection results are good, stratification is clear, and the faults are easily observed, and the results are shown in Fig. 7. Compared with other operators, the LBP/VAR operator based on guided filtering exhibits a better performance both in stratification definition and in coaxial discontinuity characterization.

Comparison of fault detection methods
The coherent body technique is a simple way to detect faults; however, it is only suitable for strata with obvious reflection characteristics and has certain limitations for strata with good continuity [25]. The LBP/VAR operator based on median filtering detects faults improves the SNR and blurs transverse discontinuities. The Canny operator based on guided filtering does not delineate all strata and leaves blank areas, which is not conducive to fault delineation. Overall, the LBP/VAR operator based on guided filtering results in a clearer detection of faults and the breakpoint location direction axis fault break. It also results in better continuity of non-tectonic faults after processing, more continuous phase axes, and preservation of structures as well as the simultaneous suppression of the effects of noise. Moreover, the proposed method enhances the image more clearly, boasts a certain anti-interference ability, and can improve horizon and fault interpretation.

Conclusions
Fault seismic responses are typically characterized by lateral discontinuities such as in-phase axial displacement or interruption of seismic waves. Improper filtering destroys the edge of the structure; therefore, it is crucial to preserve the structure, while suppressing noise. This study showed that the proposed LBP/VAR operator seismic fault detection method based on guided filtering can control the balance between noise attenuation and effective signal preservation as well as effectively detect faults in seismic data. The results obtained by the proposed method are better than those obtained using the coherent operator, LBP/ VAR operator, LBP/VAR operator based on median filtering, or Canny operator based on guided filtering. Specifically, the proposed method results in a higher SNR, with more obvious fault features that are consistent with actual faults, and a certain ability to identify small faults. Moreover, the operator is robust to noise; thus, the fault detection results can provide guidance for fault interpretation. Some input parameters (window radius r, and regularization parameter ε) are tuned in this paper for specific seismic datasets. Providing a general parameter based optimization method as input parameters and additional complicated fault characteristics will be integrated into the proposed fault detection method in future work.