CN102879824B - Quick sparse Radon transformation method based on iterative shrinkage - Google Patents

Abstract

The invention discloses a quick sparse Radon transformation method based on iterative shrinkage. The quick sparse Radon transformation method comprises the following steps: firstly, setting an initial variable value; secondly, constructing a transformation operator L and calculating generalized inverse (LTL)-1LT of the transformation operator L; thirdly, treating a seismic channel set d to be treated by utilizing the generalized inverse (LTL)-1LT of the transformation operator L; and lastly, judging if all channel sets in a seismic data cube are treated, if not, continuing to treat the seismic channel set d to be treated by utilizing the generalized inverse (LTL)-1LT of the transformation operator L, and if so, ending. According to the quick sparse Radon transformation method, for one seismic data cube collected by adopting the same collection parameters, the generalized inverse of the transformation operator L only needs to be calculated once, then the transformation operator L and the generalized inverse (LTL)-1LT of the transformation operator L are applied to all seismic channel sets, thereby greatly reducing calculated amount; and the iterative shrinkage algorithm only includes product operation of simple matrixes and vectors and threshold operation, greatly reduces the calculated amount relative to the conventional sparse Radon transformation, and better adapts to treatment of practical seismic data.

Description

A kind of rapid sparse Radon transform method based on iterative shrinkage

Technical field

The invention belongs to technical field of data processing, relate to the process of geological data in seismic exploration technique, particularly a kind of rapid sparse Radon transform method based on iterative shrinkage.

Background technology

Radon conversion (RT) is widely used in seismic data process, such as the treatment of multiple reflection and linear coherent noises compacting, wave field separation, the problems such as data normalization (Trad, D., T.Ulrych, andM.Sacchi, 2003, Latest views of the sparse Radon transform:Geophysics, 68,386 – 399).The mathematical definition of RT conversion is as follows:

d=Lm （1）

Wherein, d is known geological data, and L is known RT transformation operator, and m is RT transform domain model to be asked.RT transformation operator L converts parameter according to the acquisition parameter of geological data and RT to determine.A seismic data volume is made up of many seismic channel sets, utilizing RT transfer pair seismic data volume to carry out process is realize by processing respectively single seismic channel set, for the different earthquake road collection adopting identical acquisition parameter to gather, RT transformation operator L is identical.

For RT conversion, system of equations (1) is an over-determined systems, and conventional RT conversion is exactly the least square solution m=(LTL) asking for m ^-1l ^td.In order to improve the effect of RT conversion in seismic data process further, there has been proposed sparse Radon and converting, being exactly when solving equation group (1), supposing that solution of equation m is sparse.If utilize the openness of a L1 norm measure signal, Radon conversion is represented as following sparse indirect problem:

$\underset{m}{\arg \min }\mathrm{\λ}{\left|\right|m\left|\right|}_1+{\left|\right|d-\mathrm{Lm}\left|\right|}_2^2---\left(2\right)$

In formula, λ is a coefficient, is used between the openness of model m and the fitting precision of data d, obtain one and trades off.The 1st, the right of formula (2) utilizes L1 norm to measure the openness of m, and the 2nd is the fitting precision utilizing L2 norm to carry out metric data d.

Traditional sparse Radon conversion is by iteration again weighted least square algorithm (Scales, J., Gersztenkorn, A., and Treitel, S., 1988, Fast l _psolution of large, sparse, linear systems:Application to seismic travel time tomography:J.Comp.Phys., 75,314-333) realize, walk in iteration at t+1, need to solve following minimum L2 norm problem:

$\underset{m(t+1)}{\arg \min }{\left|\right|w_{m\left(t\right)}(d-\mathrm{Lm}(t+1))\left|\right|}_2^2---\left(3\right)$

Wherein, w _{m (t)}it is the weighting matrix that the model solution m (t) walking iteration acquisition by t constructs.Due in every single-step iteration, weighting matrix w _{m (t)}be change, so solution formula (3) all can relate to the problem that a large-scale matrix is inverted, cause traditional sparse Radon transform method calculated amount large, sparse Radon conversion rate is slower.

Summary of the invention

The present invention is intended at least solve the technical matters existed in prior art, especially innovatively proposes a kind of rapid sparse Radon transform method based on iterative shrinkage.

In order to realize above-mentioned purpose of the present invention, the invention provides a kind of rapid sparse Radon transform method based on iterative shrinkage, it comprises the steps:

S1: variable initial value is set;

S2: tectonic transition operator L and calculate the generalized inverse (L of transformation operator L ^tl) ^-1l ^t;

S3: the generalized inverse (L of the transformation operator L utilizing step S2 to obtain ^tl) ^-1l ^tto pending seismic channel set d process;

S4: judge that in seismic data volume, whether all roads collection is all disposed, and if not, returns step S3, if be all disposed, then terminates.

The seismic data volume adopting identical acquisition parameter to gather for one based on the rapid sparse Radon transform method of iterative shrinkage of the present invention, only needs to carry out once generalized inverse asking for transformation operator L, then by transformation operator L and generalized inverse (L thereof ^tl) ^-1l ^tto be saved in calculator memory and to be applied to all seismic channel sets, significantly reducing calculated amount.

In the preferred embodiment of the present invention, carry out process to seismic channel set d to comprise the steps:

S31: input a pending seismic channel set d, described seismic channel set d are n dimension, described n be greater than 1 positive integer;

S32: establish number of iterations t=0, the generalized inverse (L of the transformation operator L utilizing step S2 to obtain ^tl) ^-1l ^task for n-dimensional model m (t)=(L ^tl) ^-1l ^td;

S33: make t=t+1, obtains upgrading rear model m (t+1),

m(t+1)=T _α{m(t)+β(LTL) ^-1L ^T[d-Lm(t)]}

Wherein, T _α: R ⁿ→ R ⁿbe contraction operator, be defined as:

$T_{\mathrm{\α}}{\left\{m\right\}}_{\mathrm{ij}}={\left(\right|m_{\mathrm{ij}}|-\mathrm{\α}{\tilde{m}}_{\mathrm{ij}})}_{+}\mathrm{sgn}\left(m_{\mathrm{ij}}\right)$

Wherein, α is threshold value, and 0< α < 1, m={m _ij, a n being tieed up mean filter and is applied to | m|, obtains filter result , and

${\left(x\right)}_{+}=\left\{\begin{array}{cc}x& x\&GreaterEqual;0\\ 0& x<0\end{array}\right.;$

S34: judging whether iterations t reaches maximum iteration time N, as do not reached, returning step S33; If reached, export the result of current seismic road collection.

Iterative shrinkage algorithm of the present invention only comprises the product calculation of simple matrix and vector, threshold operation, and relative to traditional sparse Radon conversion, algorithm of the present invention significantly reduces calculated amount, more adapts to the process of actual seismic data.

Additional aspect of the present invention and advantage will part provide in the following description, and part will become obvious from the following description, or be recognized by practice of the present invention.

Accompanying drawing explanation

Above-mentioned and/or additional aspect of the present invention and advantage will become obvious and easy understand from accompanying drawing below combining to the description of embodiment, wherein:

Fig. 1 is the rapid sparse Radon transform method schematic flow sheet that the present invention is based on iterative shrinkage;

Fig. 2 is the Comparative result figure adopting prior art and the present invention respectively Prof. Du Yucang seismic channel set to be carried out to multiple reflection experiment of pressing.

Embodiment

Be described below in detail embodiments of the invention, the example of described embodiment is shown in the drawings, and wherein same or similar label represents same or similar element or has element that is identical or similar functions from start to finish.Being exemplary below by the embodiment be described with reference to the drawings, only for explaining the present invention, and can not limitation of the present invention being interpreted as.

The present invention proposes a kind of rapid sparse Radon transform method based on iterative shrinkage, be more applicable to the huge seismic data process of data volume, as shown in Figure 1, should comprise the steps: based on the rapid sparse Radon transform method of iterative shrinkage

S1: variable initial value is set, in the present embodiment, need the variable arranging initial value to comprise threshold coefficient α, iteration step length β, maximum iteration time N and RT convert parameter;

S2: convert parameter and earthquake data acquisition parametric configuration transformation operator L according to RT, calculates the generalized inverse (L of transformation operator L ^tl) ^-1l ^t;

S3: the generalized inverse (L of the transformation operator L utilizing step S2 to obtain ^tl) ^-1l ^tto pending seismic channel set d process;

S4: judge that in seismic data volume, whether all roads collection is all disposed, and if not, returns step S3, if be all disposed, then terminates.

In the preferred embodiment of the present invention, step S3 comprises the steps:

S31: input a pending seismic channel set d, this seismic channel set d is n dimension, and n be greater than 1 positive integer, namely rapid sparse Radon conversion of the present invention can be two-dimentional, also can be three-dimensional and above higher-dimension, when the dimension n that this rapid sparse Radon converts is greater than 2, m, d are the column vector form of n-dimensional model and raw data, L is that n ties up Radon operator, and n ties up Radon conversion and is used for converting n dimension data d, obtains n-dimensional model m;

S32: establish number of iterations t=0, the generalized inverse (L of the transformation operator L utilizing step S2 to obtain ^tl) ^-1l ^task for n-dimensional model m (t)=(L ^tl) ^-1l ^td;

S33: make t=t+1, obtains upgrading rear model m (t+1),

m(t+1)=T _α{m(t)+β(L ^TL) ^-1L ^T[d-Lm(t)]} （4）

Wherein, T _α: R ⁿ→ R ⁿbe contraction operator, be defined as:

$T_{\mathrm{\&alpha;}}{\left\{m\right\}}_{\mathrm{ij}}={\left(\right|m_{\mathrm{ij}}|-\mathrm{\&alpha;}{\tilde{m}}_{\mathrm{ij}})}_{+}\mathrm{sgn}\left(m_{\mathrm{ij}}\right)---\left(5\right)$

Wherein, α is threshold value, and 0< α < 1, m={m _ij, a n being tieed up mean filter and is applied to | m|, obtains filter result and

${\left(x\right)}_{+}=\left\{\begin{array}{cc}x& x\&GreaterEqual;0\\ 0& x<0\end{array}\right.---\left(6\right)$

S34: judging whether iterations t reaches maximum iteration time N, as do not reached, returning step S33; If reached, export the result of current seismic road collection.

In the present embodiment, should implement in frequency domain based on the rapid sparse Radon transform method of iterative shrinkage, in the embodiment that the present invention is other, should implement in time domain based on the rapid sparse Radon transform method of iterative shrinkage,

Iterative shrinkage algorithm of the present invention only comprises the product calculation of simple matrix and vector, threshold operation, relative to traditional sparse Radon conversion, reduces calculated amount, more adapts to the process of actual seismic data.

Fig. 2 is the Comparative result figure adopting prior art and the present invention respectively Prof. Du Yucang seismic channel set to be carried out to multiple reflection experiment of pressing, and in Fig. 2, a is the result utilizing least square Radon transfer pair Prof. Du Yucang seismic channel set to carry out multiple reflection experiment of pressing; B is the result that traditional frequency domain sparse Radon transfer pair Prof. Du Yucang seismic channel set carries out multiple reflection experiment of pressing; C is the result that sparse Radon transfer pair Prof. Du Yucang seismic channel set of the present invention carries out multiple reflection experiment of pressing.In figure a, b, c, from left to right, be input seismic channel set respectively, export RT model, primary reflection is estimated, primary reflection evaluated error.

When utilizing sparse Radon transform method of the present invention to carry out emulation experiment, convolution model Prof. Du Yucang is utilized to go out a two-dimension earthquake road collection, the i.e. Far Left subgraph of Fig. 2 a, concentrate in the two-dimension earthquake road of Prof. Du Yucang, comprise two multiple reflections and two primary reflections, multiple reflection increases with offset distance, the lineups that whilst on tour increases, primary reflection increases with offset distance, and whilst on tour reduces or constant lineups.In the present embodiment, utilize Radon transform method to carry out the flow process of multiple reflection compacting, comprising: 1) seismic channel set is transformed to RT territory, obtain RT model; 2) utilize multiple reflection and the separability of primary reflection in RT territory, designing filter, primary reflection is suppressed by RT model; 3) multiple reflection is estimated by anti-Radon conversion; 4) finally concentrate from seismic traces the multiple reflection deducting estimation, obtain the estimation of primary reflection.

As shown in Figure 2, least square Radon converts, the square error that traditional frequency domain sparse Radon conversion and sparse Radon proposed by the invention convert the primary reflection estimation obtained is respectively 0.5727,0.4886 and 0.2924, can find out, sparse Radon conversion proposed by the invention, in 3 kinds of methods, obtains best primary reflection and estimates.In addition, as can be seen from Figure 2, it is the most sparse that sparse Radon proposed by the invention converts the RT model obtained.

Composition graphs 2 is visible, the seismic data volume adopting identical acquisition parameter to gather for one based on the rapid sparse Radon transform method of iterative shrinkage of the present invention, only need to carry out once generalized inverse asking for transformation operator L, then by transformation operator L and generalized inverse (L thereof ^tl) ^-1l ^tto be saved in calculator memory and to be applied to all seismic channel sets, significantly reducing calculated amount, improve the accuracy that primary reflection is estimated.

In the description of this instructions, specific features, structure, material or feature that the description of reference term " embodiment ", " some embodiments ", " example ", " concrete example " or " some examples " etc. means to describe in conjunction with this embodiment or example are contained at least one embodiment of the present invention or example.In this manual, identical embodiment or example are not necessarily referred to the schematic representation of above-mentioned term.And the specific features of description, structure, material or feature can combine in an appropriate manner in any one or more embodiment or example.

Although illustrate and describe embodiments of the invention, those having ordinary skill in the art will appreciate that: can carry out multiple change, amendment, replacement and modification to these embodiments when not departing from principle of the present invention and aim, scope of the present invention is by claim and equivalents thereof.

Claims (2)

1., based on a rapid sparse Radon transform method for iterative shrinkage, it is characterized in that, comprise the steps:

S1: variable initial value is set, wherein, described variable comprises threshold coefficient α, and iteration step length β, maximum iteration time N and RT convert parameter;

S2: tectonic transition operator L and calculate the generalized inverse (L of transformation operator L ^tl) ^-1l ^t;

S3: the generalized inverse (L of the transformation operator L utilizing step S2 to obtain ^tl) ^-1l ^tto pending seismic channel set d process, specifically comprise:

S31: input a pending seismic channel set d, described seismic channel set d are n dimension, described n be greater than 1 positive integer,

S32: establish number of iterations t=0, the generalized inverse (L of the transformation operator L utilizing step S2 to obtain ^tl) ^-1l ^task for n-dimensional model m (t)=(L ^tl) ^-1l ^td,

S33: make t=t+1, obtains upgrading rear model m (t+1) and is:

m(t+1)＝T _α{m(t)+β(L ^TL) ^-1L ^T[d-Lm(t)]},

Wherein, T _α: R ⁿ→ R ⁿbe contraction operator, be defined as:

$T_{\mathrm{\&alpha;}}{\left\{m\right\}}_{\mathrm{ij}}={\left(\right|m_{\mathrm{ij}}|-\mathrm{\&alpha;}{\tilde{m}}_{\mathrm{ij}})}_{+}\mathrm{sgn}\left(m_{\mathrm{ij}}\right),$

Wherein, α is threshold value, and 0< α <1, m={m _ij, a n being tieed up mean filter and is applied to | m|, obtains filter result and

${\left(x\right)}_{+}=\left\{\begin{array}{cc}x& x\&GreaterEqual;0\\ 0& x<0\end{array}\right.,$

S34: judging whether iterations t reaches maximum iteration time N, as do not reached, returning step S33, if reached, exports the result of current seismic road collection;

S4: judge that in seismic data volume, whether all roads collection is all disposed, and if not, returns step S3, if be all disposed, then terminates.

2., as claimed in claim 1 based on the rapid sparse Radon transform method of iterative shrinkage, it is characterized in that, the described rapid sparse Radon transform method based on iterative shrinkage is implemented or implements in time domain in frequency domain.