CN115980854A - Three-dimensional reflection coefficient inversion method with jitter artifact suppression function - Google Patents
Three-dimensional reflection coefficient inversion method with jitter artifact suppression function Download PDFInfo
- Publication number
- CN115980854A CN115980854A CN202211393012.6A CN202211393012A CN115980854A CN 115980854 A CN115980854 A CN 115980854A CN 202211393012 A CN202211393012 A CN 202211393012A CN 115980854 A CN115980854 A CN 115980854A
- Authority
- CN
- China
- Prior art keywords
- dimensional
- reflection coefficient
- seismic
- inversion
- wavelet
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Pending
Links
Images
Landscapes
- Geophysics And Detection Of Objects (AREA)
Abstract
The invention provides a three-dimensional reflection coefficient inversion method with a jitter artifact suppression function, which comprises the steps of firstly, combining logging data and well-side seismic channels to obtain seismic wavelets by using a damped least square inversion method; then, adding a three-dimensional reflection coefficient inversion target function into a construction constraint term by utilizing a transverse second-order difference operator; and finally splitting the inversion target function into a plurality of sub-problems of a frequency domain according to a split Bregman optimization framework, and performing alternate optimization solution on the sub-problems to obtain a three-dimensional reflection coefficient inversion result. According to the method, the correlation between adjacent seismic channels in the three-dimensional space is added into the inversion objective function through second-order difference, so that the shaking artifacts in the directions of the main survey line and the cross survey line can be suppressed simultaneously, the trackability of the top-bottom interface and the small structural features of the thin stratum in the strong noise-containing seismic data on the section and the plane is obviously enhanced, and the precision of geological interpretation by using the three-dimensional reflection coefficient is improved.
Description
Technical Field
The invention belongs to the field of oil-gas seismic exploration data processing, and particularly relates to a three-dimensional reflection coefficient inversion method with an all-dimensional jitter artifact suppression function.
Background
The reflection coefficient data volume has wide application in construction boundary interpretation and thin-layer top-bottom interface identification, so the three-dimensional reflection coefficient is regarded as a bridge between seismic data and an underground geological structure. At present, methods for obtaining a three-dimensional reflection coefficient data volume include two types, namely one-dimensional inversion by channel and two-dimensional inversion by line. The channel-by-channel one-dimensional inversion method is characterized in that an inversion target function is established only in the longitudinal time axis direction by adopting sparse constraint, and because the high similarity between adjacent seismic channels is not utilized in the two directions of a main survey line and an interconnection survey line, when the signal-to-noise ratio is low, the same interface has the characteristic of micro time dislocation on each adjacent seismic channel, so that the reflection coefficient has obvious jitter instability on the integral form. The actual underground interface generally presents a steady fluctuation state, and the jitter artifacts can make the boundary of the underground structure and the top-bottom interface of the thin stratum become quite fuzzy, thereby causing great trouble to geological interpretation.
In order to suppress the jitter artifact in the three-dimensional reflection coefficient data volume, various line-by-line inversion methods are proposed at present, and generally, the methods introduce correlation information between adjacent seismic channels into an inversion target function in the main survey line direction to realize a jitter suppression function, and then perform two-dimensional reflection coefficient inversion on line-by-line to obtain a final three-dimensional reflection coefficient volume. In the line-by-line inversion methods, the correlation between adjacent seismic traces is simulated by using a Markov random field, and the local linear characteristic of the same phase axis in a time-space domain isf-kSparsity in the domain embodies this correlation, some use covariance matrices to describe the lateral correlation between seismic traces, and others use event dip information extracted from seismic data to establish the correlation between adjacent traces. However, these methods of describing the correlation between seismic traces involve only the inline direction, and the correlation between seismic traces in the crossline direction is not contained in the inverted objective function, but ratherAnd the methods are difficult to directly expand to a three-dimensional space, so that when the signal-to-noise ratio is low, the problem of the jitter of the reflection coefficient is still very prominent, and the accuracy of geological interpretation is reduced.
Disclosure of Invention
In order to solve the above problems in the prior art, embodiments of the present invention provide a three-dimensional reflection coefficient inversion method with a jitter artifact suppression function, where correlation between adjacent seismic traces in a three-dimensional space is added as a lateral constraint to an inversion objective function by using a lateral second-order difference operator, so that a jitter suppression function is obtained in a line measurement direction and a line tie measurement direction at the same time, thereby solving the problems in the prior art. The method specifically comprises the following steps:
acquiring a seismic data volume and a longitudinal wave velocity and density logging curve;
calculating a logging reflection coefficient by using the longitudinal wave velocity and density logging curves, constructing a wavelet inversion target function by using a well side channel and the logging reflection coefficient in the seismic data volume, and calculating a seismic wavelet sequence by using a damped least square method;
constructing a wavelet three-dimensional convolution operator and a transverse second-order difference three-dimensional convolution operator according to the dimensionality of the seismic data volume and the seismic wavelet sequence obtained through calculation, and calculating three-dimensional Fourier transform of the three-dimensional convolution operator;
constructing a three-dimensional reflection coefficient inversion target function with a transverse constraint term by using the transverse second-order difference three-dimensional convolution operator;
and combining a splitting Bregman algorithm frame to perform frequency domain splitting on the three-dimensional reflection coefficient inversion target function, and performing alternate optimization on each subproblem to obtain a three-dimensional reflection coefficient inversion result.
Drawings
The accompanying drawings, which are incorporated herein and form a part of the specification, illustrate one embodiment of the method of the present invention and, together with the description, further serve to explain the principles and advantages of the method and to enable a person skilled in the pertinent art to make and use the invention. In the drawings:
FIG. 1 is a block flow diagram of a three-dimensional reflection coefficient inversion method with a jitter artifact suppression function according to the present invention;
FIG. 2 is a schematic diagram of a wavelet three-dimensional convolution operator constructed using a seismic wavelet sequence according to an embodiment of the present invention;
FIG. 3 is a schematic diagram of a method for constructing a three-dimensional transverse second-order difference convolution operator according to an embodiment of the present invention;
FIG. 4 is a schematic representation of a three-dimensional seismic data volume provided by an embodiment of the invention;
FIG. 5 is a schematic diagram of a three-dimensional reflection coefficient inversion result according to an embodiment of the present invention;
FIG. 6 is a cross-sectional view of a reflection coefficient inversion result provided by an embodiment of the present invention;
fig. 7 is a time slice diagram of a reflection coefficient inversion result according to an embodiment of the present invention.
Detailed Description
The technical solution of the present invention will be described in further detail below with reference to the accompanying drawings and examples. It should be understood that the particular embodiments chosen herein are illustrative only and are not limiting of the invention, i.e., the embodiments provided are only one embodiment of the invention and not all embodiments.
The present invention is described in detail below with reference to fig. 1 to 7.
Fig. 1 is a flowchart of a three-dimensional reflection coefficient inversion method with a jitter artifact suppression function according to an embodiment of the present invention, as shown in fig. 1, including:
step S101: acquiring a seismic data volume and a longitudinal wave velocity and density logging curve;
step S102: calculating a logging reflection coefficient by using the logging curve, constructing a wavelet inversion target function by using a well side channel and the logging reflection coefficient in the seismic data volume, and calculating a seismic wavelet sequence by using a damped least square method;
step S103: constructing a frequency domain wavelet three-dimensional convolution operator and a transverse second-order difference three-dimensional convolution operator according to the dimensionality of the seismic data volume and the seismic wavelets obtained through calculation;
step S104: constructing a three-dimensional reflection coefficient inversion target function with a transverse constraint term by using the transverse second-order difference three-dimensional convolution operator;
step S105: and performing frequency domain splitting on the three-dimensional reflection coefficient inversion target function by combining a splitting Bregman algorithm, and alternately optimizing each subproblem to obtain a three-dimensional reflection coefficient inversion result.
Specifically, the acquiring the seismic data volume and the longitudinal wave velocity and density log comprises the following steps:
step 1.1: seismic data volumeCorresponding data file and longitudinal wave speed->A density->And respectively reading the data files corresponding to the two logging curves into the memory of the computer.
Further, the calculating a logging reflection coefficient by using the logging curve, constructing a wavelet inversion target function by using a well side channel and the logging reflection coefficient in the seismic data volume, and calculating a seismic wavelet sequence by using a damped least square method includes:
step 2.1: calculating a sequence of well logging reflection coefficients using well logging wave impedanceAccording to the formula:
wherein the content of the first and second substances,denotes the firstiThe log wave impedance value at each point in time,Nthe number of time sampling points of the logging sequence is counted.
Step 2.2: using the sequence of well-logging reflection coefficientsTopritz matrixes are respectively constructed by seismic channels beside wellsSum column vector:
Wherein the content of the first and second substances,Mthe number of time sampling points of a seismic trace,for the first in the seismic channel beside the welliThe amplitude values of the individual sampling points.
Step 2.3: using the Toprlitz matrixVector formed by well-side seismic tracesConstructing wavelet vectorsThe inverse objective function of (2):
wherein the content of the first and second substances,represents a two-norm vector>Indicates that the formula in the brace is minimized->,/>Is a regularization parameter.
Step 2.4: the wavelet sequence is calculated by using a damped least square methodThe formula is as follows:
wherein, the upper label "T"denotes the matrix transpose," -1 "denotes the matrix inversion, and I is the identity matrix.
Further, the constructing of the frequency domain wavelet three-dimensional convolution operator and the transverse second-order difference three-dimensional convolution operator according to the dimension of the seismic data volume and the seismic wavelet obtained through calculation comprises the following steps:
step 3.1: structural and original seismic data volumeFour zero-valued data volumes of the same scale: />、/>、/>And &>Wherein the first dimension X is the direction of the survey line, adoptThe sampling point number is 1 toGThe second dimension Y is the direction of the contact survey line, and the number of the sampling points is 1 toFThe third dimension Z is the time direction, and the number of sampling points is 1 toM;
Step 3.2: changes are made toValue of the specific position element: />、WhereinpThe sampling point number corresponding to the wavelet zero time; change>、/>And &>Value of the specific position element: />、/>、/>、/>、/>;/>、/>、/>、/>;/>、/>、/>、/>。
Dimension of three-dimensional seismic data volumeG=F=MWhen =10, the sampling point number of wavelet zero time is given in FIG. 2pAn embodiment where the wavelet three-dimensional convolution operator W is constructed from the wavelet vector W at = 5; FIG. 3 shows the construction of the transverse second order difference three-dimensional convolution operator、/>And &>In other dimensions>、/>、/>And &>Can be derived from this embodiment without additional inventive effort.
Step 3.3: computing、/>、/>And &>Three-dimensional Fourier transform of>、/>、/>And &>And three-dimensional seismic data volume>Three-dimensional Fourier transform of>。
Further, constructing a three-dimensional reflection coefficient inversion target function with a transverse constraint term by using the transverse second-order difference three-dimensional convolution operator comprises:
wherein the content of the first and second substances,is a three-dimensional reflection factor>Represents a Frobenius norm,. Sup.>Represents a norm, < > is>Represents a three-dimensional time-domain convolution, <' > or>Indicates that the formula in the brace is minimized->,μFor the time-wise sparse regularization parameter,λthe parameters are regularized for lateral constraint.
Further, the frequency domain splitting is performed on the three-dimensional reflection coefficient inversion target function by combining the splitting Bregman algorithm, and the three-dimensional reflection coefficient inversion result obtained by alternately optimizing each subproblem comprises the following steps:
step 5.1: given a maximum number of iterationsKGiven an index numberk=0, given eight three-dimensional zero-valued data volumes on the same scale as the original seismic data:given a three-dimensional fast Fourier forward-backward transform operatorAnd &>。
Step 5.2: splitting the three-dimensional reflection coefficient inversion target function into the following frequency domain optimization sub-problems:
Step 5.3: according to the following formula and、/>、/>、/>and &>For reflection coefficient in three-dimensional Fourier domainAnd intermediate variables>、/>、/>And &>Alternate updating is carried out: />
Wherein the content of the first and second substances,and &>Are respectively>And &>Conjugation of (1);
step 5.4: order toIf is greater or greater>Repeating the alternating update of step 5.3, otherwise ending the update process and calculating the three-dimensional reflection coefficient->。
FIG. 4 is a schematic diagram of a three-dimensional seismic data volume in which relatively significant random noise is present, according to an embodiment of the present invention; FIG. 5 is a three-dimensional reflection coefficient inversion result obtained by the method of the present invention, wherein the top-bottom interfaces of each stratum in the inversion result are smooth and natural and do not contain a jitter artifact; FIG. 6 is a cross-section corresponding to line 100 in the inversion result shown in FIG. 5, wherein the extension of the top and bottom interfaces of the thin layer indicated by the arrows conforms to geological rules and does not contain jitter artifacts; FIG. 7 is a time slice at 530 ms in the inversion result shown in FIG. 5, and the thin layer indicated by the arrow can also be clearly tracked and explained on the plane, further illustrating that the method of the present invention can effectively suppress the reflection coefficient jitter artifact caused by noise.
Claims (5)
1. A three-dimensional reflection coefficient inversion method with a jitter artifact suppression function, comprising:
using the well-logging reflection coefficientWell-side seismic channelConstructing wavelet inversion target function and calculating seismic wavelet vector by using damped least square method;
Three-dimensional convolution operator for constructing frequency domain wavelets with same dimension as seismic data volumeAnd a transverse second order difference three-dimensional convolution operator->、/>And &>;
Constructing a three-dimensional reflection coefficient inversion target function with a transverse constraint term by using the transverse second-order difference three-dimensional convolution operator;
2. The method of claim 1, wherein the inversion of the objective function using the well-log reflection coefficients and the well-side seismic trace formation wavelets and the computation of the seismic wavelet series using the damped least squares method comprises:
using the reflection coefficient of the logTopritz matrixes are respectively constructed by seismic channels beside wellsSum column vector:
Wherein the content of the first and second substances,Mthe number of time sampling points of a seismic trace,for the first in the seismic channel beside the welliAmplitude values of the sampling points;
wherein the content of the first and second substances,represents a two-norm vector>Indicates that the formula in the brace is minimized->,/>Is a regularization parameter;
the method for solving the formula (2) by the damped least squares method comprises the following steps:
wherein, the upper label "T"denotes a matrix transpose," -1 "denotes a matrix inversion, and I is an identity matrix.
3. The method as claimed in claim 1, wherein the frequency domain wavelet three-dimensional convolution operator and the transverse second order difference three-dimensional convolution operator with the same dimension as the seismic data volume are constructed, and the method comprises the following steps:
structural and primitive seismic data volumesFour zero-value data blocks with the same scale>、/>、/>And &>The first dimension X of these data volumes is the direction of the survey line, and the number of sample points is 1, 2.GThe second dimension Y is the direction of the contact survey line, the number of the sampling point is 1, 2.,Fthe third dimension Z is the time direction, and the number of the sampling points is 1, 2.,M;
changing the element value:、/>whereinpThe sampling point number corresponding to the wavelet zero time; changing the element value: />、/>、/>、、/>;/>、/>、/>、/>;、/>、/>、/>;
4. the method of claim 1, wherein constructing a three-dimensional reflection coefficient inversion objective function with a transverse constraint term using the transverse second order difference three-dimensional convolution operator comprises:
establishing an objective function:
wherein the content of the first and second substances,is a three-dimensional reflection factor>Represents a Frobenius norm,. Sup.>Represents a norm, < > is>Represents a three-dimensional time-domain convolution, <' > or>Means for deciding when the formula in the parenthesis is minimized>,μFor the time-wise sparse regularization parameter,λregularize the parameters for lateral constraints.
5. The method according to claim 1, wherein the objective function is subjected to frequency domain splitting by combining with a splitting Bregman algorithm, and each sub-problem after splitting is alternately optimized to obtain a three-dimensional reflection coefficient inversion resultThe method is characterized by comprising the following steps:
given a maximum number of iterationsKGiven an index numberk=0, given eight three-dimensional zero-valued data volumes on the same scale as the original seismic data:given a three-dimensional fast Fourier forward-inverse transform operator->And &>;
Splitting the three-dimensional reflection coefficient inversion target function into the following frequency domain optimization sub-problems:
reflection coefficient in three-dimensional Fourier domain according to the following formulaAnd each intermediate variable->、/>、/>Andalternate updating is carried out: />
Wherein the content of the first and second substances,and &>Are respectively>And &>Conjugation of (1);
wherein, the first and the second end of the pipe are connected with each other,,/>is a sign function;
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202211393012.6A CN115980854A (en) | 2022-11-08 | 2022-11-08 | Three-dimensional reflection coefficient inversion method with jitter artifact suppression function |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202211393012.6A CN115980854A (en) | 2022-11-08 | 2022-11-08 | Three-dimensional reflection coefficient inversion method with jitter artifact suppression function |
Publications (1)
Publication Number | Publication Date |
---|---|
CN115980854A true CN115980854A (en) | 2023-04-18 |
Family
ID=85961397
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202211393012.6A Pending CN115980854A (en) | 2022-11-08 | 2022-11-08 | Three-dimensional reflection coefficient inversion method with jitter artifact suppression function |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN115980854A (en) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116952850A (en) * | 2023-09-20 | 2023-10-27 | 安徽农业大学 | Dynamic jitter reduction inversion method based on non-dispersive infrared detection |
-
2022
- 2022-11-08 CN CN202211393012.6A patent/CN115980854A/en active Pending
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN116952850A (en) * | 2023-09-20 | 2023-10-27 | 安徽农业大学 | Dynamic jitter reduction inversion method based on non-dispersive infrared detection |
CN116952850B (en) * | 2023-09-20 | 2024-01-05 | 安徽农业大学 | Dynamic jitter reduction inversion method based on non-dispersive infrared detection |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Chen | Fast dictionary learning for noise attenuation of multidimensional seismic data | |
Wang et al. | Hankel low-rank approximation for seismic noise attenuation | |
CN109001802B (en) | Seismic signal reconstructing method based on Hankel tensor resolution | |
US7392166B2 (en) | Method for more rapidly producing the representative stochastic model of a heterogeneous underground reservoir defined by uncertain static and dynamic data | |
US20060241920A1 (en) | Method of reconstructing a stochastic model, representative of a porous heterogeneous medium, to improve its calibration by production data | |
CN111856578A (en) | Wide-azimuth prestack seismic reflection mode analysis method of tensor depth self-coding network | |
Yang et al. | Denoising of distributed acoustic sensing data using supervised deep learning | |
CN109143331B (en) | Seismic wavelet extraction method | |
CN111368247A (en) | Sparse representation regularization prestack AVO inversion method based on fast orthogonal dictionary | |
CN115980854A (en) | Three-dimensional reflection coefficient inversion method with jitter artifact suppression function | |
CN104730576A (en) | Curvelet transform-based denoising method of seismic signals | |
CN111522063A (en) | Pre-stack high-resolution fluid factor inversion method based on frequency division joint inversion | |
CN110687597B (en) | Wave impedance inversion method based on joint dictionary | |
CN113077386A (en) | Seismic data high-resolution processing method based on dictionary learning and sparse representation | |
CN114444393A (en) | Logging curve construction method and device based on time convolution neural network | |
CN110146923A (en) | A kind of efficient high accuracy depth domain methods of seismic wavelet extraction | |
Chen et al. | Integration of cumulative-distribution-function mapping with principal-component analysis for the history matching of channelized reservoirs | |
CN116068644B (en) | Method for improving resolution and noise reduction of seismic data by using generation countermeasure network | |
CN113296150B (en) | High-dimensional closed-loop network seismic inversion method under logging constraint | |
CN110865409B (en) | Seismic wave impedance inversion method based on wave impedance low-rank regularization | |
CN113722893A (en) | Seismic record inversion method, device, equipment and storage medium | |
Zhou et al. | Fast Dictionary Learning Based on Data-Driven Tight Frame for 3D Seismic Data Denoising | |
Zhou et al. | Coherent noise attenuation by kurtosis-guided adaptive dictionary learning based on variational sparse representation | |
CN111352159B (en) | Nuclear norm and generalized total variation joint constraint seismic random noise suppression method | |
Chen et al. | Dictionary Learning for Single‐Channel Passive Seismic Denoising |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination |