Lu2006.pdf

GEOPHYSICS, VOL. 71, NO. 6 共NOVEMBER-DECEMBER 2006兲; P. V163–V169, 7 FIGS., 1 TABLE. 10.1190/1.2335873

Local linear coherent noise attenuation based on local polynomial approximation

Wenkai Lu1, Wenpo Zhang2, and Dongqi Liu3

the data, which is usually inherent in this process, producing a wormy appearance in the output data. Other useful operators and separation methodologies are presented in many papers. A localized 2D filter in the Fourier projection domain 共FPF兲 proposed by Lu 共2001兲 shows good local properties and less filtering distortion in linear coherent noise suppression. In the FPF 共 ␶,p,x兲, where x is the offset of the input data, the linear coherent noise is identified by a simple threshold detector. After that, a 1D notch filter is used to remove the components corresponding to the coherent noise. Filtering approaches in the ␶ -p-x domain are discussed also by Mitchell and Kelamis 共1990兲. Abma and Claerbout 共1995兲 use a least-squares inversion method to predict both signal and noise in seismic data. Nemeth et al. 共2000兲 present a migration filtering method to separate signal and coherent noise according to their propagation paths. Linville and Meek 共1995兲 introduce a procedure for optimally removing localized coherent noise when the noise is localized in time. Dossary et al. 共2001兲 describe an adaptive linear noise attenuation technique based on patchwise ␶ -p transforms that provide good ground-roll suppression while preserving the signal, especially in the high-frequency components of the data. Lu 共2006兲 proposes an adaptive filter based on singular value decomposition and texture direction detection and demonstrates good results. The radial trace domain 共Brown and Claerbout, 2000; Henley, 2003兲 proved to be an effective domain in which to attenuate the coherent noise on both prestack and poststack seismic data. Ulrych et al. 共1999兲 give a good tutorial about the separation of seismic signal and noise. All of the methods mentioned above assume that the noise is linearly coherent. However, in some field data, the coherent noise is locally linear only because of complex acquisition circumstances 共for example, mountain and desert environments兲. In addition, it may be useful to have a more general noise model, allowing for amplitude variations of the noise and for noise events that have different moveouts but interfere in the time/distance 共t-x兲 domain. Similar to the FPF method 共Lu, 2001兲, our method uses a threshold detector to scan

ABSTRACT We propose a new technique for the attenuation of locally coherent noise. We assume that the moveout of the noise is locally linear and approximate its amplitude variations with offset using piecewise 共local兲 polynomial models. Thus, our method consists of three steps: detection of the noise 共locally linear coherent noise, LLCN兲, amplitude estimation by a local polynomial approximation 共LPA兲, and subtraction of the estimated coherent noise from the original data. Applying the proposed method to synthetic data and to a field data set shows that the LPA filter has good ability to model LLCN and is insensitive to the filter parameters. Comparisons of the results obtained by our method with those from the traditional frequency-wavenumber filter and the localized 2D filter in the Fourier projection domain 共FPF兲 show that the new method outperforms both traditional methods in situations with complex coherent noise.

INTRODUCTION In seismic applications, undesired coherent energy 共airwaves, ground roll, etc.兲 may be difficult to suppress during seismic data processing. Several methods have been developed and used in practice for filtering coherent energy. In general, these methods can be classified into two groups: global filtering and local filtering. Two of the best-known examples of global filtering methods are Radon transform 共 ␶ -p兲 filtering and frequency-wavenumber 共f-k兲 filtering. In either case, the data are transformed first into the corresponding domain 共␶ -p or f-k兲. In these two domains, filtering can be applied, either to reject or pass coherent components according to their velocities. The most undesirable aspect of a global filter is the mixing of

Manuscript received by the Editor June 18, 2004; revised manuscript received April 11, 2006; published online October 18, 2006. 1 Tsinghua University, Department of Automation, Institute of Information, Beijing 100084, China. E-mail: [email protected]. 2 Peking University, Department of Geophysics, Beijing 100871, China, and Exploration and Development Research Institute, Liaohe Oil Field Branch Cooperation, CNPC, Panjin City, Liaoning Province 124010, China. E-mail: [email protected]. 3 PetroChina Research Institute of Petroleum Exploration and Development, Institute of Geophysics, Beijing 100083, China, and Exploration and Development Research Institute, Liaohe Oil Field Branch Cooperation, CNPC, Panjin City, Liaoning Province 124010, China. E-mail: [email protected]. © 2006 Society of Exploration Geophysicists. All rights reserved.

V163

V164

Lu et al.

for the noise in a predefined velocity range, a procedure that extracts the local noise trajectory coordinates in the time-velocity domain. Subsequently, a low-rank polynomial approximation, which is a 1D filter, is implemented point by point along the detected noise trajectory to estimate the noise component. In this paper, we describe first the theory of the proposed method. Then we show the performance of our method on synthetic and field data, with comparisons to the traditional f-k and FPF approaches.

THEORY Suppose the coherent noise can be separated from the desired signal by its moveout velocity. We denote the velocity range of the coherent noise as V. For one given velocity ␯, in the range of velocities V, and a particular zero-offset time t0, a 1D signal can be extracted from the seismic shot gather d共t,x兲 along the following trajectory:

tx = t0 + x/␯ ,

共1兲

where x is the offset, and tx is the time location of the extracted 1D signal at offset x. For a given time t0, we sample the velocity range V at equal increments. The trajectory is determined by equation 1 using the time t0 and the velocity ␯. The extracted 1D signal, which is denoted as y共x, ␯,t0兲, is expressed as

y共x, ␯,t0兲 = d共tx,x兲 = d共t0 + x/␯,x兲.

共2兲

As mentioned above, the seismic record consists of three components: signal, coherent noise, and random noise. Equation 2 can be rewritten as

y共x, ␯,t0兲 = s共x, ␯,t0兲 + nc共x, ␯,t0兲 + nr共x, ␯,t0兲,

共3兲

where s共x, ␯,t0兲 is the desired signal, nc共x, ␯,t0兲 is the coherent noise that is supposed to be modeled as locally linear coherent noise 共LLCN兲 in this paper, and nr共x, ␯,t0兲 is the random noise. Because the 1D signal y共x, ␯,t0兲 is extracted along a trajectory defined by a coherent noise velocity, it will consist mostly of the coherent noise nc共x, ␯,t0兲. To remove undesired coherent energy, ␶ -p and f-k filters use only differences in moveout velocity. When these two methods are used to estimate and then subtract the noise, this is equivalent to estimating the coherent noise corresponding to 共 ␯,t0兲 as the mean value of the 1D signal y共x, ␯,t0兲,

nˆc共x, ␯,t0兲 = m共␯,t0兲 =

兺x y共x, ␯,t0兲/N,

共4兲

where N is the number of the traces involved in the above estimation and m共 ␯,t0兲 is the mean value of the extracted 1D signal. Similar to ␶ -p and f-k filters, the FPF method 共Lu, 2001兲 also uses equation 4 for the coherent noise estimation. In the FPF method, a simple detector based on a given threshold and the ratio between 兩m共 ␯,t0兲兩 and 兺x 兩y共x, ␯,t0兲兩/N is applied to decide whether the 1D signal y共x, ␯,t0兲 is mostly coherent noise or not. To improve the noise estimate in equation 4, we introduce a local polynomial approximation 共LPA兲-based estimator. To estimate the coherent noise at the offset x, only a few traces 共denoted as X兲 around offset x are used to make a low-order polynomial approximation. Suppose M traces are included and the order of the polynomial is L共L ⬍ M ⱕ N兲; the LPA is expressed as L

nˆc共xi, ␯,t0兲 =

ak共xi − x兲k, 兺 k=0

i = 1, . . ., M ,

where ak, k = 0, . . ., L are the coefficients of the polynomial.

共5兲

We rewrite equation 5 in vector-matrix form:

共6兲

nˆ c = Ca, where nˆ c = 关nˆc共xi, ␯,t兲;i = 1, . . ., M兴T, a = 关a0,a1, . . ., aL兴T and

C=

冤

1 共x1 − x兲 ¯

共x1 − x兲L

1 共x2 − x兲 ¯

共x2 − x兲L

]

]

]

]

1 共x M − x兲 ¯ 共x M − x兲L

冥

.

共7兲

The polynomial coefficients are estimated by minimizing the following objective function:

E = 储nˆ c − y储22 ,

共8兲

where y = 关y共xi, ␯,t兲;i = 1, . . ., M兴 . A least-square solution for the coefficients of the polynomial is obtained as 共Bjorck,1996兲: T

a共␯兲 = 共CTC兲−1CTy.

共9兲

Using the polynomial obtained above, an estimation of the coherent noise at offset x, nˆc共x, ␯,t0兲, is obtained. The energy ratio between the original extracted 1D signal and its approximation is given by

冑 兺 关nˆ 共x, ␯,t 兲兴 e共␯兲 = . 关y共x, ␯ ,t 兲兴 兺 冑 c

0

2

x苸X

0

共10兲

2

x苸X

Because we assume that only the coherent noise at velocity ␯ in the extracted 1D signal is predictable, a simple detector based on a given threshold ␴ is applied to decide whether the signal at offset x is mostly coherent noise or not 共Lu, 2001兲. The estimated coherent noise is

nˆ共x, ␯,t0兲 =

再

nˆc共x, ␯,t0兲 if e共␯兲 ⱖ ␴; 0

if e共␯兲 ⬍ ␴ .

冎

共11兲

By processing seismic data trace by trace, an estimation of the coherent noise corresponding to 共 ␯,t0兲 is obtained. The filtered signal is expressed as

sˆ共x, ␯,t0兲 = y共x, ␯,t0兲 − nˆ共x, ␯,t0兲.

共12兲

In summary, the proposed method and the FPF filter 共Lu, 2001兲 detect the coherent noise in the 共t0, ␯,x兲 domain using a threshold, while ␶ -p and f-k filters exploit the moveout velocity only. Furthermore, the LPA filter achieves a better coherent noise estimate by using equation 5, in comparison to these methods that estimate the coherent noise using equation 4.

EXAMPLES To demonstrate the performance of the proposed method, we apply it first to a synthetic data set based on the convolutional model. We create a synthetic shot with six hyperbolic reflections 共Figure 1a兲. The sample interval in time is 4 ms, and the dominant frequency of the Ricker wavelet used to simulate these reflections is 30 Hz. In Figure 1b, three linear events 共N1, N2, and N3兲 are superimposed on the data shown in Figure 1a to simulate the coherent noise. The dominant frequency of the Ricker wavelet used to simulate the coherent noise is 10 Hz. The reflectivities used to create these linear events vary with offset 共Table 1兲. The filtered results obtained using the f-k,

Local linear coherent noise attenuation FPF, and LPA techniques are shown in Figure 1c–e, respectively. The noise removed by the f-k, FPF, and LPA filters are shown in Figure 2a–c, respectively. In this circumstance, where the amplitudes of the coherent noise decrease with offset, the result obtained by the proposed method contains less residual noise compared to those ob-

V165

tained by the other two methods. In this experiment, each method is performed in two steps: noise estimation followed by subtraction of the noise from the data. Figure 3 shows an original shot gather of a field data set with 360 traces acquired in a desert in northwest China. The length of each

Figure 1. Synthetic data example. 共a兲 Shot gather with only six reflections. 共b兲 Shot gather with reflections and coherent noise. Result obtained by 共c兲 the f-k filter, 共d兲 the FPF filter, and 共e兲 the LPA filter.

Figure 2. Synthetic data example. The noise removed by 共a兲 the f-k filter, 共b兲 the FPF filter, and 共c兲 the LPA filter.

Figure 3. An original shot gather with complex coherent noise.

Table 1. The reflectivity versus offset x for creation of the linear events in the synthetic data. Linear event Reflectivity versus offset

N2

N1 0.3 − 0.001*

共 兲

x x − 0.00001* 25 25

2

0.25 − 0.001*

N3 x 25

0.2 − 0.0001*

共 25 兲 x

2

V166

Lu et al.

trace is 6000 ms. The sample interval is 4 ms. It is seen that the coherent noise pattern is complicated in this gather. First we apply an f-k filter to suppress the coherent noise. The result obtained by this f-k filter and the noise removed are shown in Figure 4a and b, respectively. The coherent noise is attenuated, but there are large artifacts over the whole profile. To reduce the signal distortion, the f-k filter

cannot be designed to remove the coherent noise completely. In Figure 4a, the residual coherent noise is visible. The FPF method is also applied on the same shot gather; its result is shown in Figure 4c. The threshold is 0.2. It is seen that the coherent noise is attenuated while the signal is well preserved. However, there are still some visible residual coherent noise events in Figure 4c because the FPF method

Figure 4. Results obtained by linear coherent noise suppression techniques. 共a兲 The result after f-k coherent noise attenuation. 共b兲 The noise removed by the f-k filter. 共c兲 The result after FPF coherent noise attenuation. 共d兲 The noise removed by the FPF filter. 共e兲 The result after LPAcoherent noise attenuation. 共f兲 The noise removed by our LPA method.

Local linear coherent noise attenuation cannot provide a good estimation of the noise for this complicated noise pattern. Figure 4d shows the noise removed by the FPF method. Figure 4e and f shows the result obtained by the proposed LPA method and the noise removed, respectively. The threshold ␴ is 0.4, whereas M = 21 traces are used to design the polynomial filter. The order of the polynomial is L = 4. It appears that the coherent noise is removed well by our space and time-velocity localized filter without damaging the signals.

V167

We now examine the performance of the proposed method with different choices of the parameters L and M and threshold ␴. In each experiment, only one parameter is altered, whereas the other two parameters are fixed. Some results of our experiments are shown in Figure 5. When the parameters are set to L = 0, M = 21, and ␴ = 0.4, the filtered result and the removed noise are shown in Figure 5a and b, respectively. It is seen that a lower-rank polynomial does

Figure 5. Sensitivity of the proposed method to different choices of parameters. 共a兲 Filtered result and 共b兲 the noise removed when M = 21, ␴ = 0.4, and L = 0. 共c兲 Filtered result and 共d兲 the noise removed when L = 4, ␴ = 0.4, and M = 11. 共e兲 Filtered result and 共f兲 the noise removed when L = 4, M = 21, and ␴ = 0.6.

V168

Lu et al.

not give a good approximation of the LLCN in this case. There is still some visible residual coherent noise in Figure 5a. When the parameters are set to L = 4, M = 11, and ␴ = 0.4, the filtered result and the removed noise are shown in Figure 5c and d, respectively. In Figure 5d, there are some visible reflection signals removed as part of the noise. When the parameters are set to L = 4, M = 21, and ␴ = 0.6, the filtered result and the removed noise are shown in Figure 5e and f, respectively. There is visible residual coherent noise in Figure 5e because the estimated noise is too small. It is clear that the parameters used to obtain the result shown in Figure 4e 共L = 4, M = 21, and ␴ = 0.4兲 give the best filtered result. The results shown in Figure 5 demonstrate the sensitivity of the LPA method to the parameters. To evaluate the performance of the proposed method further, we show the stack sections before and after applications of these three

methods in Figure 6. The parameters used to obtain Figure 4 are adopted to process all shot gathers, and the same mute parameters and stacking velocity are used to obtain these stack profiles. Figure 6a shows the stack section of the original data. The stack sections of the filtered results obtained by the f-k, FPF, and LPA filters are shown in Figure 6b–d, respectively. In Figure 6a, the seismic events are obscured by the coherent noise over the whole profile. After coherent noise attenuation by our method, the seismic events shown in Figure 6d are enhanced with less residual noise compared to the results obtained by the other two methods. The stack sections of the noise removed by the f-k, FPF, and LPA filters are shown in Figure 7a–c, respectively. In comparison with the f-k filter, the FPF and LPA filters preserve the signal better, and the LPA filter removes more coherent noise than the FPF filter.

Figure 6. Stack sections. 共a兲 Before coherent noise attenuation and after application of 共b兲 the f-k filter, 共c兲 the FPF filter, and 共d兲 the proposed method.

Local linear coherent noise attenuation

V169

CONCLUSIONS We propose an LPA-based method for attenuation of LLCN. The proposed method is local in the 共t0, ␯,x兲 domain and detects and estimates LLCN effectively. The performance of the proposed method was evaluated with different choices of three key parameters of the LPA filter. From these experiments, it is seen that the LPA filter is robust and not very sensitive to these three parameters. One important condition for our method to be effective is that the noise is locally linear coherent, such that the noise trajectories can be fitted piecewise by lines and the amplitude variation along the trajectory can be approximated by polynomials. In general, this condition is satisfied in seismic data processing. Note that our method does not assume constant amplitudes for the noise along the locally linear trajectories. Comparison of the profiles resulting from the application of the f-k, FPF, and LPA filters shows that the best results are obtained using our proposed technique.

ACKNOWLEDGMENTS This work was sponsored by the National Natural Science Foundation of China 共40474040兲. The author is grateful to the Associate Editor and three anonymous reviewers for their constructive remarks on this manuscript.

REFERENCES Abma, R., and J. Claerbout, 1995, Lateral prediction for noise attenuation by T-X and F-X techniques: Geophysics, 60, 1887–1896. Al-Dossary, S., B. Maddison, A. A. Buali, Y. Luo, M. Al Faraj, and Q. Li, 2001, Linear adaptive noise attenuation: 71st Annual International Meeting, SEG, Expanded Abstracts, 1993–1996. Bjorck, A., 1996, Numerical methods for least squares problems: Society for Industrial and Applied Mathematics. Brown, M., and J. Claerbout, 2000, Ground roll and radial trace transform — Revisited: Stanford Exploration Project Report SEP-103, 205–231. Henley, D. C., 2003, Coherent noise attenuation in the radial trace domain: Geophysics, 68, 1408–1416. Linville, A. F., and R. A. Meek, 1995, A procedure for optimally removing localized coherent noise: Geophysics, 60, 191–203. Lu, W., 2001, Localized 2D filter based linear coherent noise attenuation: IEEE Transactions on Image Processing, 10, 1379–1382. Lu, W., 2006, Adaptive noise attenuation of seismic images based on singular value decomposition and texture direction detection: Journal of Geophysics and Engineering, 3, 28–34. Mitchell, A. R., and P. G. Kelamis, 1990, Efficient tau-p hyperbolic velocity filtering: Geophysics, 55, 619–625. Nemeth, T., H. Sun, and G. T. Schuster, 2000, Separation of signal and coherent noise by migration filtering: Geophysics, 65, 574–583. Ulrych, T. J., M. D. Sacchi, and J. M. Graul, 1999, Signal and noise separation: Art and science: Geophysics, 64, 1648–1656.

Figure 7. Stack sections. The noise removed by 共a兲 the f-k filter, 共b兲 the FPF filter, and 共c兲 the LPA filter.