CN105319593A - Combined denoising method based on curvelet transform and singular value decomposition - Google Patents

Abstract

The invention provides a seismic data combined denoising method based on curvelet transform and singular value decomposition. The combined denoising method based on curvelet transform and singular value decomposition comprises carrying out direction control and noise variance reduction through curvelet transform; partially flattening lineups through improved singular value decomposition, sequentially processing each data point and carrying out effective suppression on the noise of the whole seismic section; and combining the curvelet transform and singular value decomposition through quadtree decomposition. The technique, combined with good denoising effect of singular value decomposition and anisotropy and multi-scale characteristics of curvelet transform, can improve the signal to noise ratio of the seismic data, protect the effective signals while suppressing the random noise, can effectively remove image artifacts, and keep dynamical features of the seismic wave, and is an effective noise suppression technique.

Description

Based on the associating denoising method of warp wavelet and svd

Technical field

The present invention relates to exploration geophysics field, particularly relate to a kind of associating denoising method based on warp wavelet and svd in seism processing.

Background technology

In seism processing, denoising is very the key link, and the quality of denoising effect determines the quality of skew, superposition, imaging, and then affects follow-up seismic data inverting and interpretation work.According to the difference of Noise Characteristic, noise can be divided into random noise and coherent noise two class.Seismic data denoising technology conventional at present comprises f-x predictive filtering, Karhunen-Loeve transformation, tau-p conversion, coherent enhancement, Polynomial Fitting Technique, svd, warp wavelet, wavelet transformation etc., various noise-removed technology has played very important effect in real data process, and achieves good effect.But, they more or less also exist, and some is not enough, comprise denoising not thoroughly, to bending or intersect that lineups effect is poor, to produce artifact, counting yield low etc., the limitation that often kind of technology has it to apply, therefore, it is very necessary for studying a kind of effective noise-removed technology.

Summary of the invention

The object of the invention is to solve above technology Problems existing, providing a kind of and can effectively remove noise, the Seismic Wave Dynamics Characteristics significantly improve the associating denoising method based on warp wavelet and svd of seismic data signal to noise ratio (S/N ratio) of remaining valid.

Associating denoising method of being somebody's turn to do based on warp wavelet and svd of the present invention comprise:

Step 1, utilizes warp wavelet threshold method to carry out denoising to seismic section;

Step 2, utilizes the singular value decomposition method improved to carry out denoising to seismic section;

Step 3, utilizes quadtree decomposition to be combined by the singular value decomposition method of warp wavelet and improvement.

The further optimisation technique scheme of such scheme:

A kind of associating denoising method based on warp wavelet and svd of optimization comprises:

In step 1, using two-dimentional real seismic record as input, to seismologic record march wave conversion, obtain corresponding bent wave system number; Subsequently at bent wave zone, threshold process is carried out to bent wave system number, the bent wave system number being less than threshold value is regarded as the bent wave system number that noise is corresponding, zero setting, the bent wave system number being greater than threshold value is regarded as and retains the threshold value that useful signal is corresponding; Finally to the bent wave system number march ripple inverse transformation after threshold process, obtain the seismologic record R1 after denoising;

In step 2, denoising is carried out to the singular value decomposition method that two-dimentional actual seismic data separate improves, lineups direction is followed the trail of in little forms, based on singular value curve undulatory property identification noise point, lineups SVD decomposition denoising is evened up in local subsequently, successively to each data point process, obtain the seismologic record R2 after processing;

In step 3, when utilizing quadtree decomposition method warp wavelet and svd to be combined, fundamental relation is represented by formula below:

R＝R1*c+R2*(1-c)

$c=\left\{\begin{array}{c}0.5+0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i>4\\ 0.5-0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i\≤4\end{array}\right.$

In formula: R represents the result finally obtained, R1 represents the result after warp wavelet process, and R1 represents the result after the svd process of improvement, S (v _i) represent that in quadtree decomposition matrix, size is v _ithe number of sub-block.

Wherein: in step 1, when utilizing warp wavelet threshold method to carry out denoising, the warp wavelet fundamental relation of frequency field is represented by formula below:

In formula: c (j, l, k) for bent wave system number, j be scale parameter, l is angle parameter, and k is direction parameter, for the input of warp wavelet, for yardstick 2 ^-j, direction θ _l, position is $x_k^{(j,l)}=R_{{\mathrm{\θ}}_l}^{-1}(k_1\·2^{-j},k_2\·2^{-j/2})$ Qu Bo;

In step 2, when utilizing singular value decomposition method to carry out denoising, fundamental relation is represented by formula below:

$A=\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{\mathrm{\σ}}_i{u}_i{v_i}^T$

In formula: seismologic record is expressed as the matrix A on M × N rank, M road altogether, often have N number of sampled point together, r is the order of matrix A, subscript T representing matrix transposition, u _imatrix A A ^ti-th proper vector, v _iit is matrix A ^ti-th proper vector of A.

The associating denoising method based on warp wavelet and svd of the second optimization comprises:

In step 1, using two-dimentional real seismic record as input, to seismologic record march wave conversion, obtain corresponding bent wave system number; Subsequently at bent wave zone, threshold process is carried out to bent wave system number, the bent wave system number being less than threshold value is regarded as the bent wave system number that noise is corresponding, zero setting, the bent wave system number being greater than threshold value is regarded as and retains the threshold value that useful signal is corresponding; Finally to the bent wave system number march ripple inverse transformation after threshold process, obtain the seismologic record R1 after denoising;

In step 2, utilize the singular value decomposition method improved to carry out denoising, now fundamental relation is represented by formula below:

$A={\mathrm{UDV}}^T=\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{\mathrm{\σ}}_i{u}_i{v_i}^T$

In formula: seismologic record is expressed as the matrix A on M × N rank, M road altogether, often have N number of sampled point together, U represents M*N rank orthogonal matrix, and V represents N*N rank orthogonal matrix, and D represents M*N rank diagonal matrix, and r is the order of matrix A, subscript T representing matrix transposition, u _imatrix A A ^ti-th proper vector, v _iit is matrix A ^ti-th proper vector of A;

The singular value decomposition method that utilization improves carries out denoising step to two-dimentional actual seismic data and is: first judge the direction calculating sampling point place place lineups; Then intercept small data body, and judge to calculate whether sampling point is noise by singular value slope variance curve, if noise, lineups if not noise, are then rotated to be level by zero setting in subrange, by svd denoising; Repeat above process successively to each data point on whole seismic section, the result obtained is just the seismologic record R2 after improved singular value decomposition method process;

In step 3, when utilizing quadtree decomposition method warp wavelet and svd to be combined, fundamental relation is represented by formula below:

R＝R1*c+R2*(1-c)

$c=\left\{\begin{array}{c}0.5+0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i>4\\ 0.5-0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i\≤4\end{array}\right.$

In formula: R represents the result finally obtained, R1 represents the result after warp wavelet process, and R1 represents the result after the svd process of improvement, S (v _i) represent that in quadtree decomposition matrix, size is v _ithe number of sub-block.

Wherein: in step 1, when utilizing warp wavelet threshold method to carry out denoising, the warp wavelet fundamental relation of frequency field is represented by formula below:

In formula: c (j, l, k) for bent wave system number, j be scale parameter, l is angle parameter, and k is direction parameter, for the input of warp wavelet, for yardstick 2 ^-j, direction θ _l, position is $x_k^{(j,l)}=R_{{\mathrm{\θ}}_l}^{-1}(k_1\·2^{-j},k_2\·2^{-j/2})$ Qu Bo

In above-mentioned two kinds of prioritization schemes in step 1, if find the ground unrest having high spud angle in the seismologic record after process, re-start warp wavelet to the data after denoising, setting direction coefficient is zero, carry out contrary flexure wave conversion again, the ground unrest of high spud angle is suppressed.

The associating denoising method based on warp wavelet and svd in the present invention, make use of warp wavelet and has good azimuth characteristic, anisotropy, multiple dimensioned feature, can protect useful signal while Attenuating Random Noise; This method is on the basis of warp wavelet simultaneously, it is combined with the singular value decomposition method improved, make use of the feature that the singular value decomposition method denoising effect of improvement is good, removes the artifact that may exist, further increases the quality of seismic section; Last this method, by the associating of warp wavelet and svd, obtains the seismologic record of high s/n ratio, is conducive to follow-up seism processing, explanation, make result more reliable.

Accompanying drawing explanation

Fig. 1 is the process flow diagram of a specific embodiment of the inventive method;

Fig. 2 is the principle schematic of the warp wavelet of the inventive method;

Fig. 3 is the geologic model figure set up;

Fig. 4 is the seismogram that Fig. 3 forward modeling goes out;

Fig. 5 is the seismogram after the inventive method is suppressed;

Fig. 6 is the noise residual plot (noise namely removed) removed through the inventive method;

Fig. 7 is a real seismic record figure;

Fig. 8 carries out the seismogram after noise compacting through the inventive method to real data;

Fig. 9 is the residual plot (noise namely removed) after the inventive method carries out noise compacting to real data.

Embodiment

For making above and other object of the present invention, feature and advantage can become apparent, cited below particularly go out preferred embodiment, and coordinate institute's accompanying drawings, be described in detail below.

Embodiment 1

As shown in Figure 1, Fig. 1 is the process flow diagram of the associating denoising method based on warp wavelet and svd of the present invention.In step 101, select and need the seismic section carrying out noise compacting.Flow process enters into step 102.

In step 102, warp wavelet threshold method is utilized to carry out denoising to seismic section.Fig. 2 is the principle schematic of the warp wavelet of the associating denoising method based on warp wavelet and svd of the present invention, and when utilizing warp wavelet threshold method to carry out denoising, the warp wavelet fundamental relation of frequency field is represented by formula below:

In formula: c (j, l, k) for bent wave system number, j be scale parameter, l is angle parameter, and k is direction parameter, for the input of warp wavelet, for yardstick 2 ^-j, direction θ _l, position is $x_k^{(j,l)}=R_{{\mathrm{\θ}}_l}^{-1}(k_1\·2^{-j},k_2\·2^{-j/2})$ Qu Bo.

When utilizing warp wavelet threshold method to carry out noise compacting, using two-dimentional real seismic record as input, to seismologic record march wave conversion, obtain corresponding bent wave system number; Subsequently at bent wave zone, threshold process is carried out to bent wave system number, the bent wave system number being less than threshold value is regarded as the bent wave system number that noise is corresponding, zero setting, the bent wave system number being greater than threshold value is regarded as and retains the threshold value that useful signal is corresponding; Finally to the bent wave system number march ripple inverse transformation after threshold process, obtain the seismologic record R1 after denoising; If find the ground unrest having high spud angle in the seismologic record after process, warp wavelet is re-started to the data after denoising, setting direction coefficient (coefficient at time-space domain correspondence high spud angle) is zero, then carries out contrary flexure wave conversion, and the ground unrest of high spud angle is suppressed.Flow process enters into step 103.

In step 103, utilize the singular value decomposition method improved, noise compression process is carried out to two-dimentional real seismic record.When utilizing singular value decomposition method to carry out denoising, fundamental relation is represented by formula below:

$A=\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{\mathrm{\σ}}_i{u}_i{v_i}^T$

In formula: seismologic record is expressed as the matrix A on M × N rank, M road altogether, often have N number of sampled point together, r is the order of matrix A, subscript T representing matrix transposition, u _imatrix A A ^ti-th proper vector, v _iit is matrix A ^ti-th proper vector of A.

When utilizing singular value decomposition method to carry out noise compacting, denoising is carried out to the singular value decomposition method that two-dimentional actual seismic data separate improves, lineups direction is followed the trail of in little forms, based on singular value curve undulatory property identification noise point, lineups SVD decomposition denoising is evened up in local subsequently, successively to each data point process, obtain the seismologic record R2 after processing.Flow process enters into step 104.

In step 104, utilize quadtree decomposition method, warp wavelet and svd are combined.When utilizing quadtree decomposition method warp wavelet and svd to be combined, fundamental relation is represented by formula below:

R＝R1*c+R2*(1-c)

$c=\left\{\begin{array}{c}0.5+0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i>4\\ 0.5-0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i\≤4\end{array}\right.$

In formula: R represents the result finally obtained, R1 represents the result after warp wavelet process, and R1 represents the result after the svd process of improvement, S (v _i) represent that in quadtree decomposition matrix, size is v _ithe number of sub-block.The result R finally obtained, is the result that the associating denoising method based on warp wavelet and svd obtains.

Fig. 3 is the geologic model figure set up, Fig. 4 is the seismogram just performed, Fig. 5 is the seismogram after suppressing based on the associating denoising method of warp wavelet and svd, and Fig. 6 is the noise residual plot through removing based on the associating denoising method of warp wavelet and svd.As can be seen from Fig. 3 to Fig. 6, the associating denoising method based on warp wavelet and svd is good for the effect of theoretical model, and noise obtains effective compacting, and does not damage the information of effective lineups.

Associating denoising method based on warp wavelet and svd is applied to actual seismic data.Fig. 7 is a real seismic record figure, Fig. 8 carries out the seismogram after noise compacting through the associating denoising method based on warp wavelet and svd to real data, and Fig. 9 carries out the residual plot after noise compacting through the associating denoising method based on warp wavelet and svd to real data.Comparative analysis Fig. 7 and Fig. 8, can find out, carries out noise pressing result well by the associating denoising method based on warp wavelet and svd to real seismic record, and medium and deep lineups become more continuous, and signal to noise ratio (S/N ratio) is significantly improved.

Embodiment 2

In step 1, utilize warp wavelet to carry out threshold method denoising, now the warp wavelet fundamental relation of frequency field is represented by formula below:

In formula: c (j, l, k) is bent wave system number, and j represents scale parameter, and l represents angle parameter, and k represents direction parameter, the input that f (x) is warp wavelet, for yardstick 2 ^-j, direction θ _l, position is qu Bo.

In step 1, Main is as follows: first to the direct transform of two-dimentional real seismic record march ripple, transformed to bent wave zone, obtains multiple bent wave system number corresponding with it; Subsequently threshold process is carried out to these bent wave system numbers, the bent wave system number being greater than threshold value is regarded as and retains the threshold value that useful signal is corresponding, the bent wave system number being less than threshold value is regarded as the bent wave system number that noise is corresponding, zero setting; Finally to the bent wave system number march ripple inverse transformation after process, reconstruct back two-dimentional real seismic record, the result obtained is just the seismologic record R1 after the process of warp wavelet denoising method.It should be noted that, if find the ground unrest having high spud angle in the seismologic record after process, warp wavelet is re-started to the data after denoising, setting direction coefficient (coefficient at time-space domain correspondence high spud angle) is zero, carry out contrary flexure wave conversion again, the ground unrest of high spud angle is suppressed.

In step 2, utilize the singular value decomposition method improved to carry out denoising, now fundamental relation is represented by formula below:

$A={\mathrm{UDV}}^T=\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{\mathrm{\σ}}_i{u}_i{v_i}^T$

In formula: seismologic record is expressed as the matrix A on M × N rank, M road altogether, often have N number of sampled point together, U represents M*N rank orthogonal matrix, and V represents N*N rank orthogonal matrix, and D represents M*N rank diagonal matrix, and r is the order of matrix A, subscript T representing matrix transposition, u _imatrix A A ^ti-th proper vector, v _iit is matrix A ^ti-th proper vector of A.

In step 2, the singular value decomposition method improved is utilized to carry out denoising to two-dimentional actual seismic data.Key step is: first judge the direction calculating sampling point place place lineups; Then intercept small data body, and judge to calculate whether sampling point is noise by singular value slope variance curve, if noise, lineups if not noise, are then rotated to be level by zero setting in subrange, by svd denoising; Repeat above process successively to each data point on whole seismic section, the result obtained is just the seismologic record R2 after improved singular value decomposition method process.

In step 3, when utilizing quadtree decomposition method warp wavelet and svd to be combined, fundamental relation is represented by formula below:

R＝R1*c+R2*(1-c)

$c=\left\{\begin{array}{c}0.5+0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i>4\\ 0.5-0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i\≤4\end{array}\right.$

In formula: R represents the result finally obtained, R1 represents the result after warp wavelet process, and R1 represents the result after the svd process of improvement, S (v _i) represent that in quadtree decomposition matrix, size is v _ithe number of sub-block.The result R finally obtained, is the result that the associating denoising method based on warp wavelet and svd obtains.

Comparative example 2 and embodiment 1, its difference is a part for step 2, but its implementation process and effect basically identical, do not elaborating.

Claims (6)

1., based on the associating denoising method of warp wavelet and svd, it is characterized in that comprising:

Step 1, utilizes warp wavelet threshold method to carry out denoising to seismic section;

Step 2, utilizes the singular value decomposition method improved to carry out denoising;

Step 3, utilizes quadtree decomposition to be combined by the singular value decomposition method of warp wavelet and improvement.

2. the associating denoising method based on warp wavelet and svd according to claim 1, is characterized in that:

In step 1, using two-dimentional real seismic record as input, to seismologic record march wave conversion, obtain corresponding bent wave system number; Subsequently at bent wave zone, threshold process is carried out to bent wave system number, the bent wave system number being less than threshold value is regarded as the bent wave system number that noise is corresponding, zero setting, the bent wave system number being greater than threshold value is regarded as and retains the threshold value that useful signal is corresponding; Finally to the bent wave system number march ripple inverse transformation after threshold process, obtain the seismologic record R1 after denoising;

In step 2, denoising is carried out to the singular value decomposition method that two-dimentional actual seismic data separate improves, lineups direction is followed the trail of in little forms, based on singular value curve undulatory property identification noise point, lineups SVD decomposition denoising is evened up in local subsequently, successively to each data point process, obtain the seismologic record R2 after processing;

In step 3, when utilizing quadtree decomposition method warp wavelet and svd to be combined, fundamental relation is represented by formula below:

R＝R1*c+R2*(1-c)

$c=\left\{\begin{array}{c}0.5+0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i>4\\ 0.5-0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i\≤4\end{array}\right.$

In formula: R represents the result finally obtained, R1 represents the result after warp wavelet process, and R1 represents the result after the svd process of improvement, S (v _i) represent that in quadtree decomposition matrix, size is v _ithe number of sub-block.

3. the associating denoising method based on warp wavelet and svd according to claim 2, is characterized in that:

In step 1, when utilizing warp wavelet threshold method to carry out denoising, the warp wavelet fundamental relation of frequency field is represented by formula below:

In formula: c (j, l, k) for bent wave system number, j be scale parameter, l is angle parameter, and k is direction parameter, for the input of warp wavelet, for yardstick 2 ^-j, direction θ _l, position is $x_k^{(j,l)}=R_{{\mathrm{\θ}}_l}^{-1}(k_1\·2^{-j},k_2\·2^{-j/2})$ Qu Bo;

In step 2, when utilizing singular value decomposition method to carry out denoising, fundamental relation is represented by formula below:

$A=\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{\mathrm{\σ}}_i{u}_i{v_i}^T$

In formula: seismologic record is expressed as the matrix A on M × N rank, M road altogether, often have N number of sampled point together, r is the order of matrix A, subscript T representing matrix transposition, u _imatrix A A ^ti-th proper vector, v _iit is matrix A ^ti-th proper vector of A.

4. the associating denoising method based on warp wavelet and svd according to claim 1, is characterized in that:

In step 1, using two-dimentional real seismic record as input, to seismologic record march wave conversion, obtain corresponding bent wave system number; Subsequently at bent wave zone, threshold process is carried out to bent wave system number, the bent wave system number being less than threshold value is regarded as the bent wave system number that noise is corresponding, zero setting, the bent wave system number being greater than threshold value is regarded as and retains the threshold value that useful signal is corresponding; Finally to the bent wave system number march ripple inverse transformation after threshold process, obtain the seismologic record R1 after denoising;

In step 2, utilize the singular value decomposition method improved to carry out denoising, now fundamental relation is represented by formula below:

$A={\mathrm{UDV}}^T=\underset{i=1}{\overset{r}{\mathrm{\Σ}}}{\mathrm{\σ}}_i{u}_i{v_i}^T$

In formula: seismologic record is expressed as the matrix A on M × N rank, M road altogether, often have N number of sampled point together, U represents M*N rank orthogonal matrix, and V represents N*N rank orthogonal matrix, and D represents M*N rank diagonal matrix, and r is the order of matrix A, subscript T representing matrix transposition, u _imatrix A A ^ti-th proper vector, v _iit is matrix A ^ti-th proper vector of A;

The singular value decomposition method that utilization improves carries out denoising step to two-dimentional actual seismic data and is: first judge the direction calculating sampling point place place lineups; Then intercept small data body, and judge to calculate whether sampling point is noise by singular value slope variance curve, if noise, lineups if not noise, are then rotated to be level by zero setting in subrange, by svd denoising; Repeat above process successively to each data point on whole seismic section, the result obtained is just the seismologic record R2 after improved singular value decomposition method process;

In step 3, when utilizing quadtree decomposition method warp wavelet and svd to be combined, fundamental relation is represented by formula below:

R＝R1*c+R2*(1-c)

$c=\left\{\begin{array}{c}0.5+0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i>4\\ 0.5-0.5\×\frac{1}{n\×n}\mathrm{\ΣS}\left(v_i\right)v_i\≤4\end{array}\right.$

In formula: R represents the result finally obtained, R1 represents the result after warp wavelet process, and R1 represents the result after the svd process of improvement, S (v _i) represent that in quadtree decomposition matrix, size is v _ithe number of sub-block.

5. the associating denoising method based on warp wavelet and svd according to claim 4, is characterized in that:

In step 1, when utilizing warp wavelet threshold method to carry out denoising, the warp wavelet fundamental relation of frequency field is represented by formula below:

In formula: c (j, l, k) for bent wave system number, j be scale parameter, l is angle parameter, and k is direction parameter, for the input of warp wavelet, for yardstick 2 ^-j, direction θ _l, position is $x_k^{(j,l)}=R_{{\mathrm{\θ}}_l}^{-1}(k_1\·2^{-j},k_2\·2^{-j/2})$ Qu Bo.

6. the associating denoising method based on warp wavelet and svd according to claim 2,3 or 4,5, it is characterized in that: in step 1, if find the ground unrest having high spud angle in the seismologic record after process, warp wavelet is re-started to the data after denoising, setting direction coefficient is zero, carry out contrary flexure wave conversion again, the ground unrest of high spud angle is suppressed.