CN105372707A - Method for attenuating multi-scale seismic data random noise - Google Patents

Abstract

The invention provides a method for attenuating the multi-scale seismic data random noise and belongs to the seismic signal processing field. The method comprises the steps of (1) conducting the Contourlet transformation on seismic data containing random noises so as to convert the seismic data into the Contourlet domain and obtain one low frequency sub-band (namely low-frequency coefficient) of the seismic data in the Contourlet domain and multiple directional sub-bands (namely high-frequency coefficient c<j, k>) of the seismic data in the Contourlet domain; (2) setting a basic threshold for the Contourlet high-frequency coefficient of the directional sub-bands at the current size of m*n; (3) calculating an adjustment factor for the Contourlet high-frequency coefficient of the directional sub-bands at the current size of m*n; (4) calculating an adaptive threshold for the Contourlet high-frequency coefficient of the directional sub-bands at the current size of m*n, wherein k = [mu]j, k lambda; (5) conducting the calculation of the Contourlet high-frequency coefficient on the high-frequency coefficient obtained in the step (1) to obtain a result c<j, k>.

Description

A kind of multi-scale seismic data random noise damped system

Technical field

The invention belongs to seismic data processing field, be specifically related to a kind of multi-scale seismic data random noise damped system, for random noise decay and the weak signal enhancement process of low signal-to-noise ratio seismic data.

Background technology

Through development for many years, Some Comments On Geophysical Work person proposes the treatment technology of a variety of compacting seismic data random noise, is broadly divided into time-space domain and transform domain two class methods.

Time-space domain random noise damped system mainly comprises overlap-add procedure, convolution filtering, Polynomial Fitting Technique, medium filtering, svd, chaotic oscillator detection technique etc.

Transform domain random noise attenuation method mainly comprises LK-conversion, the filtering of XF-territory, τ-p conversion, time-frequency method, Noise Elimination from Wavelet Transform technology, inverse Q filtering spectrum matrix technology of stablizing, the dimension method of becoming, model constrained method, complex trace analysis, zero phase spectra enhancing technology, non-physical can the technology such as implementation space predictive filtering.The wavelet analysis method being wherein particularly referred to as " school microscop " receives much concern, for many years, small wave converting method is that the research of signal transacting, image procossing and other nonlinear science brings revolutionary impact, and the analysis aspect that small wave converting method is applied to seismic data also achieves great successes.For the one-dimensional signal of " putting unusual ", little wave energy reaches the None-linear approximation rank of " optimum ", and when processing the 2D signal containing " line is unusual ", because two-dimentional separable wavelets is made up of the tensor product of one dimension small echo, although the two-dimensional wavelet transformation of such formation easily detects the point of discontinuity be positioned on edge, the continuity between marginal point but cannot be represented accurately.That is lacking direction property of higher-dimension wavelet basis, can not curve well in expression signal or face unusual, wavelet transformation the most effectively, the most sparsely can not express high dimensional signal.

Because various random noise damped system has himself relative merits and applicability, particularly when complex structure, signal to noise ratio (S/N ratio) is on the low side, traditional method is difficult to obtain satisfied result.Therefore research can suppress random noise, and injury-free other new method of signal of remaining valid again is extremely important.

Summary of the invention

The object of the invention is to solve the difficult problem existed in above-mentioned prior art, a kind of multi-scale seismic data random noise damped system is provided, according to significant wave and the random noise strong and weak difference in Contourlet territory correlativity, carry out random noise decay, improve seismic data signal to noise ratio (S/N ratio), strengthen weak useful signal.

The present invention is achieved by the following technical solutions: a kind of multi-scale seismic data random noise damped system, the steps include:

(1) contourlet transformation is carried out to the seismic data comprising random noise, geological data is transformed to Contourlet territory, obtain low frequency sub-band, the i.e. low frequency coefficient of this geological data in Contourlet territory, and multiple directions subband, i.e. high frequency coefficient c _{j, k};

(2) the basic threshold value that current size is the directional subband Contourlet high frequency coefficient of M × N is set $\mathrm{\&lambda;}=\mathrm{\&sigma;}\sqrt{2\log \mathrm{MN}};$

(3) Dynamic gene that current size is the directional subband Contourlet high frequency coefficient of M × N is calculated ${\mathrm{\&mu;}}_{j,k}=\sqrt{c_{j,k}^2/S_{j,k}^2};$

(4) the adaptive threshold λ that current size is the directional subband Contourlet high frequency coefficient of M × N is calculated _{j, k}=μ _{j, k}λ;

(5) to the high frequency coefficient c that step (1) obtains _{j, k}carry out the computing of Contourlet high frequency coefficient, obtain

(6) neighborhood window is obtained coefficient according to step (5) move that (coefficient moves one by one, each mobile one, will process each coefficient one by one.), often move once, repeat (2)-(5) step, until the process of all Contourlet high frequency coefficients completes.

(7) the Contourlet high frequency coefficient after the Contourlet low frequency coefficient utilizing step (1) to obtain and computing carry out Contourlet inverse transformation, obtain the seismic data after random noise attenuation;

(8) seismic data after random noise attenuation is exported.

In described step (2), the value of λ is determined by σ, and σ is noise criteria variance, the size of its empirical value to be 0-0.1, M and N be current directional subband.

In described step (3), c _{j, k}for Contourlet coefficient, automatically generate after contourlet transformation, for Contourlet coefficient c in this neighborhood window _{j, k}quadratic sum, specific as follows:

For the Contourlet coefficient c of jth yardstick, kth directional subband _{j, k}, A _{j, k}with c _{j, k}centered by the neighborhood window of a n × n, order for the quadratic sum of Contourlet coefficient in this neighborhood window, that is:

$S_{j,k}^2=\underset{(i,l)\&Element;A_{j,k}}{\mathrm{\&Sigma;}}{c}_{i,l}^2---\left(4\right).$

Described step (5) utilizes following formula to carry out the computing of Contourlet high frequency coefficient:

${\hat{c}}_{j,k}=\left\{\begin{array}{c}{c}_{j,k}\sqrt{1-{\mathrm{\&lambda;}}_{j,k}^2/c_{j,k}^2},\left|c_{j,k}\right|\&GreaterEqual;{\mathrm{\&lambda;}}_{j,k}\\ 0,\left|c_{j,k}\right|<{\mathrm{\&lambda;}}_{j,k}\end{array}\right.---\left(8\right).$

Compared with prior art, the invention has the beneficial effects as follows:

The image denoising territory very wide contourlet transformation of application is incorporated in seismic data denoising process by the present invention, propose the multi-scale seismic data random noise damped system based on contourlet transformation, the method can remove random noise in seismic data well, does not damage useful signal simultaneously.

Accompanying drawing explanation

Fig. 1 laplacian-pyramid filter decomposition process figure

Fig. 2 laplacian-pyramid filter synthetic schemes

Fig. 3 Contourlet neighbour coefficient relation schematic diagram

Fig. 4 adds the theogram of random noise

Fig. 5 the present invention removes the theogram after random noise

The random noise that Fig. 6 the present invention removes

Fig. 7 real seismic record

Fig. 8 the present invention removes the real seismic record after random noise

Random noise in the real seismic record that Fig. 9 the present invention removes

The step block diagram of Figure 10 the inventive method.

Embodiment

Below in conjunction with accompanying drawing, the present invention is described in further detail:

Contourlet transformation is widely applied in image denoising field, and seismic data denoising is more similar with image denoising, and contourlet transformation is applied to seism processing field by the present invention, decays for random noise.Contourlet transformation can carry out sparse expression to having the non-stationary seismic signal enriching texture features well, and have multiple dimensioned, multi-direction characteristic flexibly, significant wave shows as strong correlation in Contourlet territory, and the correlativity of random noise is extremely weak.

Contourlet transformation may also be referred to as low profile conversion, and it adopts piecewise smooth basis function to approach original image.Contourlet transformation is decomposed can be divided into two independently steps: (1) is according to the feature of image, adopt laplacian-pyramid filter (LP, LaplacianPyramid) multi-resolution decomposition is carried out to original two dimensional image, to catch the edge singular point in different scale image; (2) directional filter banks (DFB is used to the high-frequency image on each yardstick after laplacian decomposition, DirectionalFilterBank) travel direction decomposes, singular point on equidirectional is linked to be line, merges into same coefficient, i.e. Contourlet coefficient, synthesis conversion is the inverse process of decomposition transform.

Laplacian-pyramid filter decomposable process is:

(1) original signal is decomposed filtering H by low pass ₀produce the low frequency sub-band signal a of this signal;

(2) original signal is decomposed filtering H by high pass ₁produce the high frequency subband signals b of this signal.

After this low frequency sub-band signal produced is decomposed to previous step and carry out laplacian-pyramid filter (LP) decomposition, generate a low frequency signal and a high frequency subband signals.After iteration several times, laplacian-pyramid filter is decomposed original signal is resolved into a low frequency sub-band signal and a series of high frequency subband signals.Synthesis is the inverse process decomposed, and decomposition process figure is shown in Fig. 1, and synthesis flow is shown in Fig. 2, and wherein resolution filter and composite filter meet:

H ₀(z)G ₀(z)+H ₁(z)G ₁(z)＝1。(1)

It is 2 that high frequency subband signals b after laplacian pyramid filtering is decomposed by directional filter banks (DFB) resolves into ^lthe individual directional subband carving type, this directional subband is the matrix of coefficients of contourlet transformation, kth (k=0,1,2 ..., 2 ^l-1) signal of individual subband is:

$c_k\left(m\right)=(b\left(n\right)*h_k\left(n\right))(\&DownArrow;S_k)=\underset{n}{\mathrm{\&Sigma;}}b\left(n\right)h_k(S_{k}m-n)---\left(2\right)$

C _km matrix of coefficients that () is contourlet transformation, b (n) is the high frequency subband signals after laplacian pyramid filtering decomposition, S _kfor sampling matrix, h _kfor the wedge shape electric-wave filter matrix that anisotropic filter decomposes.

Directional filter banks synthesis expression formula is:

$\hat{b}\left(n\right)=\underset{k=0}{\overset{2^l-1}{\mathrm{\&Sigma;}}}\underset{m\&Element;Z^2}{\mathrm{\&Sigma;}}{c}_k\left(m\right)g_k(n-S_{k}m)---\left(3\right)$

Collection of functions two-dimensional discrete function space l ²(Z ²) one group of base, and collection of functions be called its reciprocal basis, and g _kand h _kmeet biorthogonal condition.

Multi-scale seismic data random noise damped system based on contourlet transformation is as follows:

According to the mutual relationship between the coefficient of signal after contourlet transformation, Contourlet coefficient can be divided into three kinds, i.e. neighbour coefficient, brother and sister's coefficient and paternal number, wherein neighbour coefficient is defined as: when studying a certain Contourlet coefficient, be called current coefficient, then neighbour coefficient refers to, is arranged in the Contourlet coefficient on the same directional subband adjacent position of same yardstick with current coefficient.Suitably can choose research range to determine the number of the neighbour coefficient of current coefficient, if research range chooses the rectangle frame of 3 × 3, then have 8 neighbour coefficients, edge less than 8, as shown in Figure 3.

Seismic data useful signal is after contourlet transformation, and coefficient in transform domain is mutually related and non-fully independence, if namely current C ontourlet coefficient is comparatively large, its neighbour coefficient may be also larger; And the ununified rule of random noise, its Contourlet coefficient does not have correlativity, based on this thought, the present invention proposes a kind of multi-scale seismic data random noise damped system based on contourlet transformation.The method, according to the different characteristics of Contourlet neighbour coefficient information, adjusts threshold value adaptively, and compacting random noise, strengthens useful signal.

For the Contourlet coefficient c of jth yardstick, kth directional subband _{j, k}(this coefficient automatically generates after conversion), A _{j, k}with c _{j, k}centered by the neighborhood window of a n × n, order for the quadratic sum of Contourlet coefficient in this neighborhood window, that is:

$S_{j,k}^2=\underset{(i,l)\&Element;A_{j,k}}{\mathrm{\&Sigma;}}{c}_{i,l}^2---\left(4\right)$

Wherein: $\begin{array}{c}i\&Element;[j-(n-1)/2,j+(n-1)/2]\\ l\&Element;[k-(n-1)/2,k+(n-1)/2]\end{array}---\left(5\right)$

Research neighbour coefficient correlativity take directional subband as processing unit, for a certain particular neighborhood window, if its central point c _{j, k}be positioned at the edge of this directional subband, then its neighbor coefficient distribution range can exceed this directional subband, and the fringing coefficient that during process, in subband, this coefficient of distance is nearest or angle factor replace the neighbor coefficient be beyond the boundary.

For each Contourlet coefficient c _{j, k}, M × N is the size of current coefficient place directional subband.Its adaptive threshold is set as:

λ _j，k＝μ _j，kλ(6)

Wherein for basic threshold value, σ is noise criteria variance, and Dynamic gene is:

${\mathrm{\&mu;}}_{j,k}=\sqrt{c_{j,k}^2/S_{j,k}^2}---\left(7\right)$

Corresponding adaptive thresholding value function is:

${\hat{c}}_{j,k}=\left\{\begin{array}{c}{c}_{j,k}\sqrt{1-{\mathrm{\&lambda;}}_{j,k}^2/c_{j,k}^2},\left|c_{j,k}\right|\&GreaterEqual;{\mathrm{\&lambda;}}_{j,k}\\ 0,\left|c_{j,k}\right|<{\mathrm{\&lambda;}}_{j,k}\end{array}\right.---\left(8\right)$

As shown in Figure 10, the multi-scale seismic data random noise damped system performing step based on contourlet transformation is:

1) contourlet transformation is carried out to the seismic data comprising random noise, geological data is transformed to Contourlet territory, obtain the low frequency sub-band (i.e. low frequency coefficient) of this geological data in Contourlet territory and multiple directions subband (i.e. high frequency coefficient).

2) the basic threshold value that current size is the directional subband Contourlet high frequency coefficient of M × N is set the value of λ is determined primarily of σ, and σ is noise criteria variance, the size of its empirical value to be 0-0.1, M and N be current directional subband.

For four direction subband, 0-45 degree and 180-225 degree are first directional subband, 45-90 degree and 225-270 degree are second directional subband, and 90-135 degree and 270-315 degree are the 3rd directional subband, and 135-180 degree and 315-360 degree are four direction subband.

3) Dynamic gene that current size is the directional subband Contourlet high frequency coefficient of M × N is calculated wherein c _{j, k}for Contourlet coefficient, automatically generate after contourlet transformation, for the quadratic sum of Contourlet coefficient in this neighborhood window, see formula (4).

4) the adaptive threshold λ that current size is the directional subband Contourlet high frequency coefficient of M × N is calculated _{j, k}=μ _{j, k}λ, λ are the 2nd) calculate in step, μ _{j, k}the 3rd) calculate in step).

5) computing of Contourlet high frequency coefficient is carried out according to formula (8).

6) neighborhood window is moved one by one that (coefficient moves one by one, each mobile one, will process each coefficient.), repeat 2-5 step, until the process of all Contourlet high frequency coefficients completes.

7) the Contourlet high frequency coefficient (namely step (6) obtains) after adopting Contourlet low frequency coefficient and computing carries out Contourlet inverse transformation, obtains the seismic data after random noise attenuation.

8) terminate.

Below by Benq theogram and real seismic record in the effect of the multiple dimensioned earthquake random noise damped system of contourlet transformation.In this example, the theogram (see Fig. 4) selecting to add random noise and a real seismic record (see Fig. 7) test the denoising effect of the present invention in seism processing.

As can be seen from Fig. 4, there is very strong random noise in theogram, the signal to noise ratio (S/N ratio) of data is very low.After the multiple dimensioned earthquake random noise damped system denoising based on contourlet transformation, the signal to noise ratio (S/N ratio) of theogram is improved significantly, and sees Fig. 5.Fig. 6 is the random noise that the present invention removes, and has no useful signal in random noise, and as can be seen here, the present invention not only has good suppression to random noise in seismic data, and has very high fidelity to signal.

Fig. 7 is a real seismic record in certain work area, and as seen from the figure, the existence of random noise has had a strong impact on the part in the quality of data, particularly figure in square frame, and the signal to noise ratio (S/N ratio) of data is very low, almost can't see useful signal.Fig. 8 be in Fig. 7 data through removing the result after random noise based on the multiple dimensioned random noise damped system of contourlet transformation, the overall looks now recorded are greatly improved, signal to noise ratio (S/N ratio) is also greatly improved, and the useful signal in figure in square frame is high-visible.Fig. 9 is the random noise removed, in noise, do not see useful signal.Can find by analyzing, based on the multiple dimensioned random noise damped system of contourlet transformation while effectively removing random noise, not damaging useful signal.

The denoising result of theogram and real seismic record shows, the present invention can remove the random noise in seismic data effectively, does not damage useful signal simultaneously.

Contourlet transformation is widely applied in image denoising field, seismic data denoising is more similar with image denoising, contourlet transformation is applied to seism processing field by the present invention, propose a kind of multiple dimensioned earthquake random noise damped system based on contourlet transformation, the method first carries out contourlet transformation to seismic data, according to the difference of the Contourlet coefficient of useful signal and random noise, adaptively the Contourlet coefficient of random noise is processed, Contourlet inverse transformation is done to the Contourlet coefficient after process and namely obtains the seismologic record after removing random noise.The denoising result of theogram and real seismic record shows, the present invention can remove the random noise in seismic data effectively, does not damage useful signal simultaneously.

Technique scheme is one embodiment of the present invention, for those skilled in the art, on the basis that the invention discloses application process and principle, be easy to make various types of improvement or distortion, and the method be not limited only to described by the above-mentioned embodiment of the present invention, therefore previously described mode is just preferred, and does not have restrictive meaning.

Claims (4)

1. a multi-scale seismic data random noise damped system, is characterized in that: described method comprises:

(1) contourlet transformation is carried out to the seismic data comprising random noise, geological data is transformed to Contourlet territory, obtain low frequency sub-band, the i.e. low frequency coefficient of this geological data in Contourlet territory, and multiple directions subband, i.e. high frequency coefficient c _{j, k};

(2) the basic threshold value that current size is the directional subband Contourlet high frequency coefficient of M × N is set $\mathrm{\&lambda;}=\mathrm{\&sigma;}\sqrt{2\log \mathrm{MN}};$

(3) Dynamic gene that current size is the directional subband Contourlet high frequency coefficient of M × N is calculated ${\mathrm{\&mu;}}_{j,k}=\sqrt{c_{j,k}^2/S_{j,k}^2};$

(4) the adaptive threshold λ that current size is the directional subband Contourlet high frequency coefficient of M × N is calculated _{j, k}=μ _{j, k}λ;

(5) high frequency coefficient obtained step (1) carries out the computing of Contourlet high frequency coefficient, obtains

(6) neighborhood window is obtained coefficient according to step (5) move one by one, often move once, repeat (2)-(5) step, until the process of all Contourlet high frequency coefficients completes;

(7) the Contourlet high frequency coefficient after the Contourlet low frequency coefficient utilizing step (1) to obtain and computing carry out Contourlet inverse transformation, obtain the seismic data after random noise attenuation;

(8) seismic data after random noise attenuation is exported.

2. multi-scale seismic data random noise damped system according to claim 1, it is characterized in that: in described step (2), the value of λ is determined by σ, σ is noise criteria variance, the size of its empirical value to be 0-0.1, M and N be current directional subband.

3. multi-scale seismic data random noise damped system according to claim 2, is characterized in that: in described step (3), c _{j, k}for Contourlet coefficient, automatically generate after contourlet transformation, for Contourlet coefficient c in this neighborhood window _{j, k}quadratic sum, specific as follows:

For the Contourlet coefficient c of jth yardstick, kth directional subband _{j, k}, A _{j, k}with c _{j, k}centered by the neighborhood window of a n × n, order for the quadratic sum of Contourlet coefficient in this neighborhood window, that is:

$S_{j,k}^2=\underset{(i,l)\&Element;A_{j,k}}{\mathrm{\&Sigma;}}{c}_{i,l}^2---\left(4\right).$

4. multi-scale seismic data random noise damped system according to claim 3, is characterized in that: described step (5) utilizes following formula to carry out the computing of Contourlet high frequency coefficient:

${\hat{c}}_{j,k}=\left\{\begin{array}{c}{c}_{j,k}\sqrt{1-{\mathrm{\&lambda;}}_{j,k}^2/c_{j,k}^2},\left|c_{j,k}\right|\&GreaterEqual;{\mathrm{\&lambda;}}_{j,k}\\ 0,\left|c_{j,k}\right|<{\mathrm{\&lambda;}}_{j,k}\end{array}\right.---\left(8\right).$