CN103399348B - Based on the Denoising of Seismic Data method of Shearlet conversion - Google Patents

Abstract

The invention discloses a kind of Denoising of Seismic Data method based on Shearlet conversion, it comprises the steps: 1. to read in two-dimension earthquake cross-sectional data; 2. by two-dimension earthquake cross-sectional data expand to the square formation that length and width are odd number ; 3. frequency domain direction bank of filters is constructed; 4. by each transformation matrix respectively with signal vector is multiplied, then does two-dimension fourier inverse transformation respectively, obtains the Shearlet coefficient of all directions and yardstick

Description

Based on the Denoising of Seismic Data method of Shearlet conversion

Technical field

The present invention relates to seismic data processing technology field, be specifically related to a kind of Denoising of Seismic Data method converted based on discrete shearing wave.

Background technology

Seismic data processing is the key link in geophysics, but due to sampling environment and the restriction of Sampling techniques, the seismic signal of acquisition has noise usually, so the denoising of seismic signal is an important step.So-called Denoising of Seismic Data slackens noise exactly, improves the signal to noise ratio (S/N ratio) of geological data, resolution and fidelity, is convenient to later use seismic data interpretation signal content.

Shearlet conversion is proposed by Guo etc. the earliest, and it is a kind of New Multi-scale geometric analysis instrument inheriting curve ripple and profile ripple advantage, by generating the shearing wave function with different characteristic to basic function convergent-divergent, shearing and the conversion of translation equiaffine.For 2D signal, shearing wave can detect all singular points, and can self-defined filtering direction, and Shearlet has more simple mathematic(al) structure, for multidimensional data geometric representation provides theoretical tool more flexibly.Meanwhile, Shearlet has the near optimum None-linear approximation performance identical with Curvelet and multiresolution analysis characteristic.These good characteristics that Shearlet has, for the difficulties solved in signal transacting provides new Research Thinking.Application major part at present for shearlet conversion is image domains.Glenn Easley (2006) has delivered the sparse orientation representation discrete Shearlet conversion being used for image.Discrete shearlet conversion is used for Postprocessing technique by R.Gomathi (2010), proposes a kind of new generalized expectation-maximization (GEM) model, effectively can reduce computing time, improve image restoration quality.The direction major part of the domestic Shearlet of utilization conversion concentrates in the process to SAR image, Zhang little Hua, Chen Jiawei etc. of Xian Electronics Science and Technology University propose the non-local mean SAR image Speckle reduction based on non-lower sampling Shearlet and directional weighting neighborhood window, directional weighting neighborhood window based on a large amount of texture edges etc. is combined for the Speckle reduction of SAR image with non-local mean computation model, while suppressing coherent spot, effectively maintains texture, marginal information that SAR image comprises.The people such as Liu Shuaiqi, Hu Shaohai of Beijing Jiaotong University reach the object of denoising in being represented by the coefficient of Shearlet conversion introducing SAR image, also obtain good denoising effect.Shearlet conversion is used for image watermark and adds algorithm by Liao Xin, considers the different perception of human eye to image border district and texture area, suitably selects, efficiently solve the contradiction between robustness and vision invisibility to embed watermark intensity.

The noise of seismic signal can be divided into organized noise and random noise from the classification of its noise characteristic.Organized noise has typical space-time feature due to it, and its minimizing technology often adopts time domain/space and the method for carrying out filtering.As the spatial filterings such as f-x territory filtering can obtain the object of the ground roll noise well removing rule, Radon converts filtering can multiple suppression reflection effectively.And for random noise due to its time domain, spatial distribution is random, frequency distribution feature is not obvious, and process is difficulty comparatively.At present conventional traditional median filter method is mainly included, frequency domain filtering method, F-K filtering method to the method that seismic signal carries out denoising, and F-K the Method of Deconvolution, based on small wave converting method etc.But all there is the problem of distortion and seismic signal edge fog etc. in various degree in these methods.Since the nineties in 20th century, a collection of scholar artificially represented with Donoho, Candes etc. proposes the super wavelet multi_resolution analysis technology of the thought of combined with wavelet transformed thought and multi-scale geometric analysis.This technology is that seismic data processing field brings new thought and new possibility: Herrmann (2007) takes the lead in adopting Curvelet conversion to establish complete seismic signal multiscale analysis flow process.Neelamani (2008) have employed Curvelet conversion and has carried out the preliminary experiment of seismic data noise attenuation.Hennenfent (2008) proposes and adopts nonuniformity to adopt Curvelet transfer pair seismic signal to carry out denoising, to some extent solves the ill-defined phenomenon brought in previous Curvelet denoising.

Shearlet conversion is a kind of multi-scale geometric analysis method just put forward recent years, is constructed by the expansion of a basic function, shearing and translation transformation.Relatively convert with Curvelet, it has the more excellent approximation capability to signal profile.Wavelet transformation base is being not fixing direction, but can freely determine, can the singular point of better trace signals, and the accuracy requirement of seismic signal to singular value point is also higher.In to the process of seismic signal, there is more sparse expression.And utilize the best rarefaction representation of each component of seismic signal, can restraint speckle signal under the prerequisite ensureing signal original waveform significantly.

Summary of the invention

For above-mentioned prior art, order of the present invention is how to provide a kind of Denoising of Seismic Data method converted based on discrete shearing wave, by the noise separation in original seismic signal, earthquake information is fully utilized, and improves the accuracy of oil and gas detection.

For reaching above-mentioned purpose, the present invention adopts following technical scheme:

Based on a Denoising of Seismic Data method for Shearlet conversion, it is characterized in that, comprise the steps:

Step one: read in one group of two-dimension earthquake cross-sectional data, data are saved as matrix S, read the size of this matrix for [m n], decomposed class k and Shearlet simultaneously setting Shearlet decomposes the Meyer small echo adopted;

Step 2: two-dimension earthquake cross-sectional data S is expanded to according to following formula the square formation S that length and width are odd number ₁:

S ₁＝MS

Wherein, M is a transformation matrix, and size is:

A=max (m, n) is if max is odd number

A=max (m, n) is if+1 max is even number;

Step 3: structure frequency domain direction bank of filters;

Step 4: by each transformation matrix respectively with S ₁signal vector is multiplied, then does two-dimension fourier inverse transformation respectively, obtains the Shearlet coefficient c of all directions and yardstick _i,j.

Step 5: threshold process; Determine the threshold value of this subband according to the ratio of the noise variance under the variance of the coefficient of all directions subband under each yardstick and this directional subband of this yardstick, formula is as follows:

${\mathrm{\τ}}_{i,j}={\mathrm{\σ}}_{{\mathrm{\ϵ}}_{i,j}}^2/{\mathrm{\σ}}_{i,j,n}^2$

Wherein represent the variance of the n-th coefficient of i-th directional subband under jth yardstick, under expression j yardstick, the noise variance in i direction, is determined the threshold value of each subband, converts the shearlet coefficient obtained and do threshold process, obtain the conversion coefficient under all yardsticks to signal by ratio;

Step 6: obtain signal after denoising by doing shearlet inverse transformation through the shearlet conversion coefficient of threshold process.

In the present invention, described step 3 comprises the following steps:

1. the direction grid that following formula construction is horizontal and longitudinal:

Y＝X ^T

2. dimension calculation anisotropic filter number (this sentences 4 yardsticks is example) is passed through according to following formula

$n_s=\underset{j=0}{\overset{k-1}{\Σ}}{2}^{j+1}$

m _s＝n _s+1

3. to each yardstick, according to following formula, the new grid matrix after converting is obtained to the yardstick of structure X and Y and shear transformation:

$Y=(i\×\sqrt{j})X+\sqrt{j}Y$

X＝jY

Wherein, X, Y are respectively the grid of the both direction in 3.a, and i, j are respectively the numbering in the direction of directional filter banks and the numbering of yardstick;

4. under current scale, following Meyer small echo auxiliary function is constructed:

$v\left(x\right)=\left\{\begin{array}{cc}0& x<0\\ 35x^4-84x^5+70x^6-20x^7& 0\&le;x\&le;1\\ 1& x>1\end{array}\right.$

$b\left(\mathrm{\&omega;}\right)=\left\{\begin{array}{cc}sin\left(\frac{\mathrm{\&pi;}}{2}v\right(\left|2^{-j}\mathrm{\&omega;}\right|-1))& 1\&le;\left|\mathrm{\&omega;}\right|\&le;2\\ \cos \left(\frac{\mathrm{\&pi;}}{2}v\right(\frac{1}{2}\left|2^{-j}\mathrm{\&omega;}\right|-1))& 2\&le;\left|\mathrm{\&omega;}\right|\&le;4\\ 0& otherwise\end{array}\right.$

5. based on the frequency domain base in both direction on above-mentioned Auxiliary frequency domain and scaling function:

${\mathrm{\&psi;}}_1\left(2^{-j}\mathrm{\&omega;}\right)=\sqrt{b^2\left(2^{-j+1}\mathrm{\&omega;}\right)+b^2\left(2^{-j}\mathrm{\&omega;}\right)}$

${\mathrm{\&psi;}}_2\left(2^{-j}\mathrm{\&omega;}\right)=\left\{\begin{array}{cc}\sqrt{v(1+2^{-j}\mathrm{\&omega;})}& \mathrm{\&omega;}\&le;0\\ \sqrt{v(1-2^{-j}\mathrm{\&omega;})}& \mathrm{\&omega;}>0\end{array}\right.$

6. on described grid, construct the diagonal line hourglass shape window of two different directions, whole territory be divided into some parts:

$C^h=\{({\mathrm{\&omega;}}_1,{\mathrm{\&omega;}}_2)\&Element;R^2:|{\mathrm{\&omega;}}_1|\&GreaterEqual;\frac{1}{2},|{\mathrm{\&omega;}}_2|<|{\mathrm{\&omega;}}_1\left|\right\}$

$C^v=\{({\mathrm{\&omega;}}_1,{\mathrm{\&omega;}}_2)\&Element;R^2:|{\mathrm{\&omega;}}_1|\&GreaterEqual;\frac{1}{2},|{\mathrm{\&omega;}}_2|>|{\mathrm{\&omega;}}_1\left|\right\}$

7., on frequency window, each yardstick corresponding, obtains frequency domain filter according to following formula:

$SH({\mathrm{\&omega;}}_1,{\mathrm{\&omega;}}_2)={\mathrm{\&psi;}}_1\left({\mathrm{\&omega;}}_1\right){\mathrm{\&psi;}}_2\left(\frac{{\mathrm{\&omega;}}_2}{{\mathrm{\&omega;}}_1}\right)\mathrm{\&chi;}$

8. step is repeated 3. 4. 5. 6., until obtain all yardsticks and the transformation matrix in direction.

Compared with prior art, the present invention has following beneficial effect:

The present invention will do laplacian decomposition containing noisy seismic signal, then shearing wave function is utilized to carry out filtering process, obtain corresponding shearing wave coefficient, by threshold process filtering noise signal, the signal of denoising is being recovered by inverse non-lower sampling shearing wave conversion, good denoising effect can be obtained, there is good practical value.

Embodiment

Below in conjunction with embodiment, the invention will be further described.

Step one: step one: read in one group of two-dimension earthquake cross-sectional data, data are saved as matrix S, read the size of this matrix for [m n], decomposed class k and Shearlet simultaneously setting Shearlet decomposes the Meyer small echo adopted; Adopt Meyer small echo as mother wavelet in this patent.

Step 2: two-dimension earthquake data S is expanded the matrix S for length and width and odd number ₁;

Step 3: structure transverse direction and longitudinal grid, utilizes Meyer small echo to set up wave filter in each yardstick all directions, the frequency domain subdivision of shape pair signals;

Step 4: by each transformation matrix respectively with S ₁signal vector is multiplied, then does two-dimension fourier inverse transformation respectively, obtains the Shearlet coefficient c of all directions and yardstick _i,j.

Step 5: threshold process determines the threshold value of this subband according to the ratio of the noise variance under the variance of the coefficient of all directions subband under each yardstick and this directional subband of this yardstick

${\mathrm{\&tau;}}_{i,j}={\mathrm{\&sigma;}}_{{\mathrm{\&epsiv;}}_{i,j}}^2/{\mathrm{\&sigma;}}_{i,j,n}^2$

Wherein represent the variance of the n-th coefficient of i-th directional subband under jth yardstick, the noise variance in i direction under expression j yardstick.Determined the threshold value of each subband by ratio, the shearlet coefficient obtained is converted to signal and does threshold process, obtain the conversion coefficient under all yardsticks.

Step 6: obtain signal after denoising by doing shearlet inverse transformation through the shearlet conversion coefficient of threshold process

Step 7: the feasibility and the validity that adopt theoretical model test-based examination context of methods; After denoising, signal significantly remains most of peak point of original signal, effectively can reflect the geological information representated by original signal.And the noise major part of signal is all removed after denoising, also remain the effective constituent of original signal preferably simultaneously, confirm the reliability of scheme herein.

Step 8: the frequency domain subdivision of noisy geological data section is shown, the effective constituent of what low-frequency component was more remain signal, and the noise performance of signal is obvious in high frequency, effectively the effective constituent of signal is separated with noise contribution through shearlet conversion.For follow-up process lays a solid foundation.

Step 9: true pre-stack seismic road collection denoising experiment; The effective information of the noise removed not containing seismic signal, confirms the validity of context of methods.

Claims (2)

1., based on the Denoising of Seismic Data method of Shearlet conversion, it is characterized in that, comprise the steps:

Step one: read in one group of two-dimension earthquake cross-sectional data, data are saved as matrix S, read the size of this matrix for [m n], decomposed class k and Shearlet simultaneously setting Shearlet decomposes the Meyer small echo adopted;

Step 2: two-dimension earthquake cross-sectional data S is expanded to according to following formula the square formation S that length and width are odd number ₁:

S ₁＝MS

Wherein, M is a transformation matrix, and size is:

A=max (m, n) is if max is odd number

A=max (m, n) is if+1 max is even number;

Step 3: structure frequency domain direction bank of filters;

Step 4: by each transformation matrix respectively with S ₁signal vector is multiplied, then does two-dimension fourier inverse transformation respectively, obtains the Shearlet coefficient c of all directions and yardstick _i,j;

Step 5: threshold process; Determine the threshold value of this subband according to the ratio of the noise variance under the variance of the coefficient of all directions subband under each yardstick and this directional subband of this yardstick, formula is as follows:

Wherein represent the variance of the n-th coefficient of i-th directional subband under jth yardstick, under expression j yardstick, the noise variance in i direction, is determined the threshold value of each subband, converts the shearlet coefficient obtained and do threshold process, obtain the conversion coefficient under all yardsticks to signal by ratio;

Step 6: obtain signal after denoising by doing shearlet inverse transformation through the shearlet conversion coefficient of threshold process.

2. the Denoising of Seismic Data method based on Shearlet conversion according to claim 1, it is characterized in that, described step 3 comprises the following steps:

1. the direction grid that following formula construction is horizontal and longitudinal:

Y＝X ^T

2. dimension calculation anisotropic filter number is passed through according to following formula

m _s＝n _s+1

3. to each yardstick, according to following formula, the new grid matrix after converting is obtained to the yardstick of structure X and Y and shear transformation:

X＝jY

Wherein, X, Y are respectively the grid of both direction, and i, j are respectively the numbering in the direction of directional filter banks and the numbering of yardstick;

4. under current scale, following Meyer small echo auxiliary function is constructed:

5. based on the frequency domain base in both direction on above-mentioned Auxiliary frequency domain and scaling function:

6. on described grid, construct the diagonal line hourglass shape window of two different directions, whole territory be divided into some parts:

7., on frequency window, each yardstick corresponding, obtains frequency domain filter according to following formula:

8. step is repeated 3. 4. 5. 6., until obtain all yardsticks and the transformation matrix in direction.