阅读说明：本技术 一种基于波阻抗低秩正则化的地震波阻抗反演方法 (Seismic wave impedance inversion method based on wave impedance low-rank regularization ) 是由 陈雷平 李曙 贺达江 丁黎明 于 2019-12-02 设计创作，主要内容包括：本发明公开了一种基于波阻抗低秩正则化的地震波阻抗反演方法。本发明直接针对地震波阻抗建立低秩正则化项,并在此基础上建立具有波阻抗低秩约束的目标函数,最后巧妙地利用奇异值分解进行求解,利用波阻抗的低秩性进行地震反演,充分利用了地震波阻抗数据的低秩性先验信息且提高了地震波阻抗反演的准确性。(The invention discloses a seismic wave impedance inversion method based on wave impedance low-rank regularization. The method directly establishes a low-rank regularization item aiming at the seismic wave impedance, establishes a target function with low-rank wave impedance constraint on the basis, finally solves the problem by skillfully utilizing singular value decomposition, performs seismic inversion by utilizing the low-rank property of the wave impedance, fully utilizes the low-rank prior information of seismic wave impedance data and improves the accuracy of seismic wave impedance inversion.)

1. A seismic wave impedance inversion method based on wave impedance low-rank regularization is characterized by comprising the following steps:

step one, inputting seismic data and an initial wave impedance model, and constructing an initial low-rank approximate matrix of the wave impedance model: constructing an initial low-rank approximate matrix by adopting a similarity search-based method;

step two, establishing a wave impedance low-rank regularization target function

Let the entire wave impedance profile be divided into K sub-data blocks with overlap therebetween, wherein the ith sub-data blockData block noted d_iThen, the wave impedance optimal low rank approximation problem is expressed as:

Φ_iis a similarity measure and selection operator with d_iIs a reference block for selecting the sum d from the seismic wave impedance m_iQ sub-data blocks with the smallest euclidean distance of (d) are arranged into a matrix of size P × Q, P representing sub-data block d_iLength of (d); o_iObtaining an optimal low-rank approximate matrix to be solved; | | non-woven hair_FExpressing the Frobenius norm of the matrix, and rank () is the operation of solving the rank of the matrix; λ represents a low-rank constrained regularization parameter;

on the basis of the formula (1), constructing an objective function of a wave impedance low-rank regularization seismic wave impedance inversion problem:

wherein G represents a forward operator in seismic wave impedance inversion, S represents seismic data obtained by observation, and η is a regularization coefficient of a wave impedance low-rank regularization term used for balancing a data fidelity term

step three, solving the objective function

For the objective function of the formula (2), the solution is carried out in two steps, firstly, the optimal low rank approximation is solved, namely, the optimal value of the formula (1) is solved, and the optimal low rank approximation of the formula (1)Writing into:

solving equation (3) using singular value decomposition and hard threshold shrinkage, let Φ_iThe singular value of m is decomposed into:

Φ_im＝Αdiag(ω)Β^H(4)

a and B respectively represent two unitary matrixes generated after singular value decomposition; omega is a singular value vector obtained by singular value decomposition; diag (ω) represents a diagonal matrix spanned by singular value vectors ω obtained by singular value decomposition, and the element on the diagonal is ω; h represents a conjugate operator; optimal low rank approximation

where λ represents the threshold used in hard threshold contraction, HT_λ() Represents a hard threshold function defined by the formula:

where Γ represents the input HT_λ() Data in the function;

in finding O_iAfter the optimal solution of (3), O is fixed in the objective function of equation (2)_iThe update formula of the wave impedance inversion result is obtained as follows:

equation (7) is a least squares problem, and equation (7) is formulated with respect to

wherein G is^*And phi_i ^*Are G and phi respectively_iUsing conjugate gradient method to quickly solve to obtain inversion resultStep four, inversion results are obtained

2. The seismic wave impedance inversion method based on wave impedance low rank regularization of claim 1 wherein η is determined by cross validation, i.e., adjusted based on the difference between the inversion result and the true value, and a value η is selected to minimize the inversion error.

3. The seismic wave impedance inversion method based on wave impedance low-rank regularization as claimed in claim 1, wherein in the first step, the step of constructing an initial low-rank approximation matrix of wave impedance comprises the following steps: arranging initial wave impedance into a long column vector m according to columns; then determining a sub-data block d as a reference, and assuming that the length of the sub-data block d is P; sequentially selecting Q sub-data blocks with the minimum Euclidean distance from d from m; and arranging the Q sub-data blocks into a matrix, and finishing the construction of the initial low-rank approximate matrix.

Technical Field

The invention belongs to the field of seismic exploration, and particularly relates to a seismic wave impedance inversion method based on wave impedance low-rank regularization.

Background

Seismic wave impedance is an important earth medium parameter that is numerically equal to the product of the density of the subsurface rock and the velocity of the seismic waves as they propagate in the rock. It is closely related to the lithology, porosity, etc. of the rock. Therefore, how to effectively and accurately obtain the seismic wave impedance is always an important research content in seismic exploration. The seismic wave impedance inversion technology reversely deduces the wave impedance of the underground medium by using prior information according to data such as earthquake, well logging and the like obtained by observation, and is the most main way for obtaining the seismic wave impedance at present. With the development of seismic acquisition technology and signal processing technology, the quality of seismic data is higher and higher. The quality of the seismic inversion results is more and more obviously limited by the utilization of prior information. How to mine and effectively use prior information becomes a hot problem in the current seismic inversion field.

The low rank nature of the data is an important a priori information that has been widely focused and used in recent years. Such as: a large number of literature reports based on data low-rank prior exist in the fields of voice and image signal denoising, seismic data denoising/reconstruction and the like. In the field of seismic exploration, the low-rank characteristic is mainly used for the aspects of seismic data denoising, missing seismic channel completion, seismic data reconstruction and the like. In these studies, the low rank property of seismic data was utilized. Namely: the matrix formed by combining the seismic data of different sub-blocks with similar waveforms in the seismic data has the characteristic of low rank. When the low rank property is used, the denoising/reconstruction effect of the data is obviously improved compared with the traditional method. Thus, these methods experimentally demonstrate the advantage of low rank as a priori information for use in seismic data processing. It should be noted that these seismic data processing methods only utilize the low rank property of seismic data in the data domain. The output results after they are processed are still only seismic data.

Since seismic data has low rank properties between different sub-blocks, seismic data is closely related to the nature and structure of the subsurface rock formations. The former is a signalized representation of the latter. If there is similarity between different sub-blocks of the former, the corresponding position of the latter should also have such similarity. In this case, the matrix formed by the wave resistances of the sub-blocks, such as wave impedance, density, poisson ratio, and the like, has low rank. Low rank is a non-localized a priori information. In electromagnetic parametric inversion, the idea of using low rank approximation is mentioned in the article "A low-rank approximation for large-scale 3D controlled-source electromagnetic Gauss-Newton inversion" published in Geophysics by ManuelAmaya et al (2016). However, their idea is to use its low rank approximation to replace the Hessian matrix in the inversion process to reduce the computational complexity, and not to directly utilize the low rank of the model parameters to be inverted. At present, no literature report is available on the utilization of the low rank property of wave impedance, and the wave impedance is obtained by using a seismic inversion method according to seismic data.

The existing seismic wave impedance inversion method mainly utilizes the seismic data or the localized prior information of the wave impedance. For example, the sparsity prior information utilized by sparse inversion is a localized prior information. It is a utilization of sparsity of local regions of seismic data or wave impedance under the L0, L1, or Lp (0< p <1) norm, lacking the utilization of a priori information between different sub-regions. In recent years, the seismic wave impedance inversion method based on total variation regularization which is concerned with uses gradient sparsity prior information of wave impedance, but the prior information is also localized, and non-localized prior information is not used.

In sum, the low rank characteristic has outstanding advantages in the fields of signal/data denoising, reconstruction and the like as important non-local prior information, so that the low rank characteristic plays an important role in seismic data processing. The seismic wave impedance should also have low rank. However, none of the existing seismic wave impedance inversion methods utilize this a priori information.

Disclosure of Invention

In order to solve the problems, the invention discloses a seismic wave impedance inversion method based on wave impedance low-rank regularization. The method directly establishes the low-rank regularization item aiming at the seismic wave impedance, establishes the target function with the low-rank constraint of the wave impedance on the basis, finally solves the problem by skillfully utilizing singular value decomposition, performs seismic inversion by utilizing the low-rank property of the wave impedance, fully utilizes the low-rank property of seismic wave impedance data and improves the accuracy of seismic wave impedance inversion.

In order to achieve the purpose, the technical scheme of the invention is as follows:

a seismic wave impedance inversion method based on wave impedance low-rank regularization comprises the following steps:

step one, inputting seismic data and a wave impedance initial model, and constructing an initial low-rank approximate matrix of wave impedance: constructing an initial low-rank approximate matrix by adopting a similarity search-based method;

step two, establishing a wave impedance low-rank regularization target function

Let the whole wave impedance profile be divided into K sub-data blocks with overlap, in which the ith sub-data block is denoted as d_iThen, the wave impedance optimal low rank approximation problem is expressed as:

Φ_iis a similarity measure and selection operator with d_iIs a reference block for selecting the sum d from the seismic wave impedance m_iQ sub-data blocks with the smallest euclidean distance of (d) are arranged into a matrix of size P × Q, P representing sub-data block d_iLength of (d); o_iObtaining an optimal low-rank approximate matrix to be solved; | | non-woven hair_FExpressing the Frobenius norm of the matrix, and rank () is the operation of solving the rank of the matrix; λ represents a low-rank constrained regularization parameter;

on the basis of the formula (1), constructing an objective function of a wave impedance low-rank regularization seismic wave impedance inversion problem:

wherein G represents a forward operator in seismic wave impedance inversion, S represents seismic data obtained by observation, and η is a regularization coefficient of a wave impedance low-rank regularization term used for balancing a data fidelity term

And the magnitude of the role played by the low-rank regularization term in inversion, wherein the larger the value is, the larger the role played by the low-rank regularization term is;

step three, solving the objective function

For the objective function of the formula (2), the solution is carried out in two steps, firstly, the optimal low rank approximation is solved, namely, the optimal value of the formula (1) is solved, and the optimal low rank approximation of the formula (1)Writing into:

solving equation (3) using singular value decomposition and hard threshold shrinkage, let Φ_iThe singular value of m is decomposed into:

Φ_im＝Αdiag(ω)Β^H(4)

a and B respectively represent two unitary matrixes generated after singular value decomposition; omega is a singular value vector obtained by singular value decomposition; diag (ω) represents a diagonal matrix spanned by singular value vectors ω obtained by singular value decomposition, and the element on the diagonal is ω; h represents a conjugate operator; optimal low rank approximationThe solution of (a) is written as:

where λ represents the threshold used in hard threshold contraction, HT_λ() Represents a hard threshold function defined by the formula:

where Γ represents the input HT_λ() Data in the function;

in finding O_iAfter the optimal solution of (3), O is fixed in the objective function of equation (2)_iThe update formula of the wave impedance inversion result is obtained as follows:

equation (7) is a least squares problem, and equation (7) is formulated with respect to

The equation of (c):

wherein G is^*And phi_i ^*Are G and phi respectively_iUsing conjugate gradient method to quickly solve to obtain inversion result

Step four, inversion results are obtained

And restoring the original wave impedance model into a 2-dimensional section or a 3-dimensional data volume according to a process opposite to the process of constructing the long column vector m, namely completing the whole inversion process.

In a further improvement, η is determined by cross-validation by adjusting the inversion result to the difference between the true value and the true value, and selecting the value η that minimizes the inversion error, wherein the seismic data is the seismic signal recorded by the geophone.

In a further improvement, in the first step, the step of constructing an initial low-rank approximation matrix of wave impedance comprises the following steps: arranging an initial wave impedance model into a long column vector m according to columns; then determining a sub-data block d as a reference, and assuming that the length of the sub-data block d is P; sequentially selecting Q sub-data blocks with the minimum Euclidean distance from d from m; arranging the Q sub-data blocks into a matrix, and completing construction of an initial low-rank approximate matrix; the initial wave impedance is a 2-dimensional matrix or a 3-dimensional data volume.

Drawings

FIG. 1 is a flow chart of the construction of an initial low rank approximation matrix; in the figure, (a) represents a wave impedance profile, (b) represents similar blocks with similar euclidean distances, and (c) an initial low rank approximation matrix;

FIG. 2 is a flow chart of the present invention;

FIG. 3 is a plot of a portion of data taken from the Marmousi2 model as a true model of wave impedance;

FIG. 4 is an initial model of wave impedance;

FIG. 5 is seismic data for inversion;

FIG. 6 is a Rake wavelet with a dominant frequency of 40 Hz;

fig. 7 is a diagram of inversion results.

Detailed Description

The invention is further explained with reference to the drawings and the embodiments.