Denoising of seismic data via multi-scale ridgelet transform
Denoising of seismic data via multi-scale ridgelet transform∗
Henglei Zhang 1 2
Tianyou Liu 1 2
Yuncui Zhang 0
0 The 5th Geophysical Prospecting Brigade, Zhongnan Petroleum Bureau , Xiangtan 411104 , China
1 Key Laboratory of Tectonics and Petroleum Resources (China University of Geosciences), Ministry of Education , Wuhan 430074 , China
2 Institute of Geophysics and Geomatics, China University of Geosciences , Wuhan 430074 , China
Noise has traditionally been suppressed or eliminated in seismic data sets by the use of Fourier filters and, to a lesser degree, nonlinear statistical filters. Although these methods are quite useful under specific conditions, they may produce undesirable effects for the low signal to noise ratio data. In this paper, a new method, multi-scale ridgelet transform, is used in the light of the theory of ridgelet transform. We employ wavelet transform to do sub-band decomposition for the signals and then use non-linear thresholding in ridgelet domain for every block. In other words, it is based on the idea of partition, at sufficiently fine scale, a curving singularity looks straight, and so ridgelet transform can work well in such cases. Applications on both synthetic data and actual seismic data from Sichuan basin, South China, show that the new method eliminates the noise portion of the signal more efficiently and retains a greater amount of geologic data than other methods, the quality and consecutiveness of seismic event are improved obviously as well as the quality of section is improved.
ridgelet transform; multi-scale; random noise; sub-band decomposition; complex Morlet wavelet CLC number; P315; 63 Document code; A
Random noise is a stubborn enemy in seismic
exploration and to suppress it is an important work
. Denoising is the premise and foundation for
high-resolution seismic processing, imaging, inversion
of lithology parameters and attributes analysis. Because
of complicated surface including water system in plain,
outcrop in carbonate area and in sand-shale area, the old
stratum age, poor signal reflection, as well as
complicated subsurface structures, these complicated geologic
conditions and interferences will lead to low signal to
noise ratio. So we can not pick up the first-breaks
precisely, can not obtain high-accuracy velocity analysis,
can not effectively imaging for pre-stack seismic data.
discussed the problems such as investigation
of surface structure, shot receiving factor, interference
wave and observation system design, and proposed the
corresponding countermeasure and technology
according to exploration practice. It can increase the profile
data’s quality significantly with the optimization of
seismic acquisition technology. However, the complex
surface and subsurface geological structure are still the
difficulties for seismic acquisition, processing and
Wavelet transform has been successfully used in
many scientific fields such as image compression, image
denoising, signal processing, computer graphics, and
(Leblanc and Morris, 2001; Liu et al,
2007b; Gao et al, 2006)
, which are largely attributed to
the capability of wavelet time-frequency analysis. When
2-D wavelet is used to process image data, it can only
detect information in horizontal, vertical, diagonal, and
several other limited directions, so the linear
characteristics of the signals can not be well described. To the
shortcomings of wavelet transform, on the foundation of
wavelet theory, Candès et al put forward ridgelet
(Candès and Donoho, 1999; Donoho, 2001)
could process signals that contain linear characteristics
and many changes in the direction. As an improvement
to wavelet analysis, ridgelet transform can get a better
result in seismic data processing
(Zhang et al, 2007, 2008;
Bao et al, 2007)
, but it also has its limitations. Candès
and Donoho (1999), Arivazhagan et al (2006) and Zhang
et al (2007) pointed out when used to process the
seismic data, the ridgelet transform could only deal with the
linear events, while the actual seismic form also exist
curving characteristic and many other characteristics.
In response to this limitation of ridgelet transform,
we firstly block the signals on the basis of multi-scale
decomposition, and then process every block signal in
different scales by ridgelet transform. The idea is that, in
sufficiently fine scale, a curving singularity can look
straight, so ridgelet analysis works well in it.
2 The method and principle
2.1 Ridgelet transform
To overcome the weakness of wavelet in higher
Candès and Donoho (1999)
alternative system of multi-resolution analysis, called
ridgelet, which can effectively deal with line-like
phenomena in two dimension.
The ridgelet transform in R2 can be defined as
Pick a smooth univariate function ψ : R→R with
sufficient decay and vanishing mean, ∫Rψ (x)dx = 0 .
For each a > 0, each b∈R and each θ∈[0, 2π), define a
bivariate function ψa,b,θ : R2→R2 by
ψ a,b,θ (x1, x2 ) =
1 ψ ⎜⎛ x1 cosθ + x2 sinθ − b ⎞
| a | ⎝ a ⎠
This function is constant along “ridges”, i.e., x1cos θ +
x2sinθ = const. Transverse to these ridges it is a wavelet;
hence the name ridgelet. Given an integrable bivariate
function f (x), define its ridgelet transform as
R f (a, b,θ ) =< f ψ, a,b,θ >= ∫R2ψ a,b,θ (x) f (x)dx.
Our hypotheses on ψ guarantee that Kψ=
∫R (|ψ (ω ) |2 / |ω |2 )dω < ∞ , and if suppose further that ψ
is normalized we have Kψ = ∫R (|ψ (ω ) |2 / |ω |2 )dω = 1 ,
so that the exact reconstruction formula is
f (x) = 4π1a3 ∫02π ∫−+∞∞ ∫0+∞
R f (a, b,θ )ψ a,b,θ (x)dadbdθ . (3)
The ridgelet transform is similar to the 2-D
continuous wavelet transform except that the point
parameters are replaced by line parameters. In 2D, points and
lines are related by Radon transform, thus the wavelet
and ridgelet transforms are linked via the Radon
transform. The Radon transform of an object f is the
collection of line integrals indexed by (θ, ρ) (θ∈ [0, 2π)),
which is given by
R(ρ ,θ ) =
∫R ∫R f (x1, x2 )δ (ρ − x2 sinθ − x1 cosθ )d x1d x2 ,
where δ is the dirac distribution.
2.2 Ridgelet transform based on multiscale analysis
(Candès and Donoho, 1999;
Arivazhagan et al, 2006; Zhang et al, 2007; Zhang and Liu,
had pointed out the ridgelet’s limitation when it is
used to process seismic data. For the complicated
seismic records, it may be applicable to use the calculus
based on the idea of separation. The idea is that, at
sufficiently fine scale, a curving singularity looks straight, so
ridgelet transform can work well in such cases. The
multi-scale ridgelet transform refers to the ridgelet
transform with certain length and width of the block at
certain scales; it can also be said tower structure by a
series of window-ridgelet transform, we can call it local
In this paper, we use wavelet analysis to do
sub-band decomposition. For the signal f, do the
(Mallat, 1989; Mallat and Huang,
f = AN + ∑ Dj ,
where AN is N-order approximation, and Dj is j-order
We now present a sketch of the algorithm: 1
apply the wavelet algorithm with N scales; 2 set n1=2M,
where n1 is the block number at 1 scale; 3 for j=1, 2, ⋅⋅⋅,
N, decompose the subband Dj with a block number nj
and apply the ridgelet transform to each block; 4 set
nj=2M−j+1, where nj is the block number at j scale.
Considering the random noise included in the
details, we note that the coarse description of the data AN is
not processed. We use the default value 2M = 2N in our
implementation. Finally, Figure 1 gives an overview of
the organization of the algorithm.
2.3 Non-linear thresholding
pioneered a denoising scheme by
using hard thresholding and soft thresholding. Hard
thresholding sets any coefficient less than the threshold
as zero, which is given by
while soft thresholding is subtracted from any
coefficient which is greater than the threshold, that is
⎪⎧R f (a,b,θ ) − t R f (a,b,θ ) ≥ t
Rˆ f (a,b,θ ) = ⎨ . (7)
⎩⎪0 R f (a,b,θ ) < t
From Figure 2, we can see the hard threshold
function is discontinuous at |Rf |=t, due to this discontinuity
at the threshold, the hard-thresholding function is known
to yield abrupt artifacts in the denoised signal, especially
when the noise level is significant, which will make
noise cannot be effectively reduced. The soft threshold
function has permanent bias and it loses part of edge
information in reconstruction due to neglecting
singularity detecting. To avoid the discontinuity caused by using
the hard-thresholding model and the biased estimation
caused by using the soft-thresholding model, the
non-linear thresholding is adopted as following:
⎪⎧R f (a, b,θ ) R f (a, b,θ ) ≥ t1
Rˆ f (a, b,θ ) = ⎨R f (a, b,θ ) ⋅ϕ (t, t1) t ≤ R f (a, b,θ ) < t1, (8)
⎩⎪0 R f (a, b,θ ) < t
where ϕ(t, t1) is the weight function damping from 1 to 0
at interval (|t|, |t1|).
3 Synthetic model data
First, we formulated the functions of the
signal-to-noise ratio (RSN) and the mean square error (rMSE).
Suppose x is original signal and y is noisy signal, they
are all 2-D signals with M rows and N columns. Then
∑ ∑ xi2j
RSN = M N
∑ ∑ ( yij − xij )2
∑ ∑ (xij − yij )2
rMSE = i=1 j=1
M ⋅ N
A pure Ricker wavelet distribution of synthetic
data is generated (Figure 3a), and the dominant
frequencies of the two waveforms are 35 Hz and 25 Hz. Adding
random noise to the pure data, the result (its
signal-to-noise ratio is 0.1) is shown in Figure 3b, which is
the data set that we used to test the denoising methods.
The first denoising procedure to be attempted is
wavelet transform and its result is in Figure 3c. The
random noise does remove to some extent, but the form of
the resulting time series has been corrupted. Specifically,
the waveforms are quite different. Some waveforms
have been reduced, which is more evident at the tail of
the form. Further, considerable noise still exists in the
de-noised data throughout the entire map. The objective
of the exercise is to recover the pure waveform from the
denoised data. The result of wavelet transform on the
synthetic data shows apparently that some noises still
exist and some reflection waves have been suppressed.
Applying the multi-scale ridgelet transform to the data
in Figure 3b, we get the waveform shown in Figure 3d.
This process has nearly eliminated all the random noise,
and the quality of section is improved. Furthermore, the
reconstruction of the signal is better than the result by
wavelet transform, and the waveforms are preserved
Figure 4 gives the results of spectrum analysis by
complex Morlet wavelet
(Sheng et al, 2007; Liu et al,
2007a; Lu et al, 2002)
for the synthetic data (50th trace).
The dominant frequencies of the two wavelets are 35 Hz
and 25 Hz, respectively (Figure 4a). The spectral
analysis of the noised data (Figure 4b) indicates that random
noise in the entire frequency band has obscured the
useful signal. Figure 4c is the spectral analysis for Figure 3c,
the random noise has been suppressed, but it loses some
useful signal (we can see the signal of 35 Hz has been
suppressed), so this method may lead to information
distortion. Figure 4d is the spectral analysis for Figure
3d, the multi-scale ridgelet transform can get better
noise suppression, while effective reflected wave is well
maintained; the result has better information fidelity.
The denoising ability of the two methods
mentioned above is efficiently assessed by calculating both
the mean-square error (MSE) and the signal-to-noise
ratio (SNR) of the data. A summary of the results is
listed in Table 1. As regards to the elimination of noise
and retention of significant signals, it is evident that the
multi-scale ridgelet transform method provides
significantly better result. In fact, the MSE of the multi-scale
ridgelet transform result is nearly one-half of the
wavelet transform’s. Furthermore, it is also clear that the SNR
of the multi-scale ridgelet transform result is nearly
about 1.64 times of the wavelet transform’s.
4 Field data processing
We processed seismic data from an oil field in
Sichuan basin, South China by the new method in this
Figure 5a shows the original seismic section with
much interference from the survey area and its
signal-to-noise ratio is low. Figure 5b shows the section
processed by wavelet transform. We can see that wavelet
thresholding method can reduce noise to some extent,
but it does not detect the seismic events effectually, so it
will lose some event information. Besides, the event’s
continuity of the profile is not good.
Figure 5c gives the processing result by multi-scale
ridgelet transform. It is clear that the noise is greatly
reduced, the quality and continuity of seismic events are
improved obviously, and the signal to noise ratio of
seismic section is improved perfectly. Figure 6 gives the
spectral analysis on a random column in Figures 5a, 5b
and 5c. In Figure 6c, the bandwidth of dominant
frequency and useful information have been reserved well
and its effect of random noise attenuation at relatively
low-frequency and relatively high-frequency is better
than wavelet transform method’s. The multi-scale
ridgelet transform can suppress the random noise drastically
with the precondition of maintenance for the effective
Origin seismic data (a) and denoised results by wavelet transform (b) and multi-scale ridgelet transform (c).
In every phase of the seismic data processing, such
as the estimate of residual static corrections, velocity
analysis, migration, etc, as a prerequisite, the seismic
data must be of a certain signal to noise ratio. On the
basis of multi-scale wavelet analysis, this paper
proposes a new wavelet transform. In sufficiently fine scale,
a curving singularity looks straight, so ridgelet works
well in it. In other words, this can overcome the
limitation of ridgelet transform that can only handle the
seismic signals containing linear events. Comparison of
multi-scale ridgelet transform denoising with wavelet
denoising shows the former can significantly eliminate
noise and increase retention of geologic signals,
therefore the quality and consecutiveness of seismic events
are improved obviously; meanwhile the signal to noise
ratio and resolution of seismic section are improved.
Acknowledgements This paper is supported by
China Petrochemical key project during the 11th
Five-year Plan as well as the Doctorate Fund of Ministry
of Education of China (No.20050491504). Thanks to
three anonymous reviewers for their constructive
comments on the manuscript, and also to Dr. S. Morris
Cooper for his valuable suggestions.
Arivazhagan S , Ganesan L and Subash Kumar T G ( 2006 ). Texture classification using ridgelet transform . Pattern Recognition Letters 27 : 1 875 − 1 883.
Bao Q Z , Gao J H and Chen W C ( 2007 ). Ridgelet domain method of ground-roll suppression . Chinese J Geophys 50 ( 4 ): 1 210 − 1 215 (in Chinese with English abstract) .
Candès E J and Donoho D L ( 1999 ). Ridgelets: a key to high-dimensional intermittency . Phil Trans Roy Soc London A 357: 2 495 − 2 509.
Donoho D L ( 1995 ). Denoising by soft-thresholding . IEEE Trans Inf 41 : 613 − 627 .
Donoho D L ( 2001 ). Ridge functions and orthonormal ridgelets . Journal of Approximation Theory 111 : 143 − 179 .
Gao J H , Mao J , Man W S , Chen W C and Zheng Q Q ( 2006 ). On the denoising method of prestack seismic data in wavelet domain . Chinese J Geophys 49 ( 4 ): 1 155 − 1 163 (in Chinese with English abstract) .
Leblanc G E and Morris W A ( 2001 ). Denoising of aeromagnetic data via the wavelet transform . Geophysics 66 ( 6 ): 1 793 − 1 804.
Li L X ( 2005 ). The problem and counter measure in seismic survey in southern marine carbonate area of China . Geophysical Prospecting for Petroleum 44 ( 5 ): 529 − 536 (in Chinese with English abstract).
Li Q Z ( 1993 ). The Path Toward Precise Exploration-Engineering Analysis of High-Resolution Seismic Exploration System . Oil Industry Press, Beijing, 112 − 115 (in Chinese).
Liu T Y , Sheng Q H , Yang Y S and Liu D W ( 2007a ). Application of frequency spectrum analysis of complex wavelet to seismic data processing . Oil Geophysical Prospecting 42 ( B08 ): 72 − 75 (in Chinese with English abstract).
Liu T Y , Yang Y S , Li Y Y , Feng J and Wu X Y ( 2007b ). The order-depression solution for large-scale integral equation and its application in the reduction of gravity data to a horizontal plane . Chinese J Geophys 50 ( 1 ): 275 − 281 (in Chinese with English abstract).
Lu X C , Gong S G , Zhou J and Sun M ( 2002 ). High resolution signal spectrum estimation based on Morlet wavelet . Journal of Wuhan University of Technology 26 ( 5 ): 622 − 635 (in Chinese with English abstract).
Mallat S ( 1989 ). A theory for multi-resolution signal decomposition: The wavelet representation . IEEE Transactions on Pattern Analysis and Machine Intelligence 11 ( 7 ): 674 − 693 .
Mallat S and Huang W L ( 1992 ). Singularity detection and processing with wavelets . IEEE Trans Information Theory 38 ( 2 ): 617 − 643 .
Sheng Q H , Liu T Y and Liu D W ( 2007 ). The spectrum decomposition with complex wavelet transforms and its applications in reservoir prediction . Geological Science and Technology Information 26 ( 6 ): 88 − 90 (in Chinese with English abstract).
Zhang H L and Liu T Y ( 2007 ). Study of ridgelet transform in improving S/N ratio of seismic data in South China . Journal of China University of Geosciences 18 : 507 − 509 .
Zhang H L , Song S and Liu T Y ( 2007 ). The ridgelet transform with non-linear threshold for seismic noise attenuation in South China . Applied Geophysics 4 ( 4 ): 271 − 275 .
Zhang H L , Song S and Liu T Y ( 2008 ). A method of processing for low SNR seismic data by ridgelet transform . Coal Geology & Exploration 36 ( 3 ): 63 − 66 (in Chinese with English abstract).