CN115113265A - Seismic data resolution improving method and device based on compressed sensing - Google Patents
Seismic data resolution improving method and device based on compressed sensing Download PDFInfo
- Publication number
- CN115113265A CN115113265A CN202110285637.XA CN202110285637A CN115113265A CN 115113265 A CN115113265 A CN 115113265A CN 202110285637 A CN202110285637 A CN 202110285637A CN 115113265 A CN115113265 A CN 115113265A
- Authority
- CN
- China
- Prior art keywords
- seismic
- spectrum
- amplitude spectrum
- record
- resolution
- 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
- 238000000034 method Methods 0.000 title claims abstract description 91
- 238000001228 spectrum Methods 0.000 claims abstract description 196
- 238000004422 calculation algorithm Methods 0.000 claims description 16
- 230000008569 process Effects 0.000 claims description 12
- 238000000354 decomposition reaction Methods 0.000 claims description 10
- 230000001131 transforming effect Effects 0.000 claims description 8
- 230000010355 oscillation Effects 0.000 claims description 5
- 230000010354 integration Effects 0.000 claims description 4
- 239000011159 matrix material Substances 0.000 claims description 4
- 230000008859 change Effects 0.000 claims description 3
- 230000002708 enhancing effect Effects 0.000 claims 1
- 230000006870 function Effects 0.000 description 28
- 230000003595 spectral effect Effects 0.000 description 19
- 238000010586 diagram Methods 0.000 description 9
- 238000012545 processing Methods 0.000 description 9
- 238000011160 research Methods 0.000 description 9
- 230000000694 effects Effects 0.000 description 7
- 238000005070 sampling Methods 0.000 description 6
- 238000004364 calculation method Methods 0.000 description 4
- 238000005516 engineering process Methods 0.000 description 4
- 238000005457 optimization Methods 0.000 description 4
- 238000004458 analytical method Methods 0.000 description 3
- 230000015572 biosynthetic process Effects 0.000 description 3
- 239000007789 gas Substances 0.000 description 3
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 2
- 241000446313 Lamella Species 0.000 description 2
- 238000001514 detection method Methods 0.000 description 2
- 238000000605 extraction Methods 0.000 description 2
- 238000001914 filtration Methods 0.000 description 2
- 238000005755 formation reaction Methods 0.000 description 2
- 230000036039 immunity Effects 0.000 description 2
- 230000006872 improvement Effects 0.000 description 2
- 238000012360 testing method Methods 0.000 description 2
- 230000009466 transformation Effects 0.000 description 2
- 101100148545 Caenorhabditis elegans snr-5 gene Proteins 0.000 description 1
- 206010011971 Decreased interest Diseases 0.000 description 1
- YZCKVEUIGOORGS-IGMARMGPSA-N Protium Chemical compound [1H] YZCKVEUIGOORGS-IGMARMGPSA-N 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 238000012512 characterization method Methods 0.000 description 1
- 230000009194 climbing Effects 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000009795 derivation Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 238000005553 drilling Methods 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
- 238000009499 grossing Methods 0.000 description 1
- 230000007774 longterm Effects 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000002156 mixing Methods 0.000 description 1
- 238000010606 normalization Methods 0.000 description 1
- 239000013307 optical fiber Substances 0.000 description 1
- 238000005215 recombination Methods 0.000 description 1
- 230000006798 recombination Effects 0.000 description 1
- 230000000717 retained effect Effects 0.000 description 1
- 238000002922 simulated annealing Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000003786 synthesis reaction Methods 0.000 description 1
- 238000000844 transformation Methods 0.000 description 1
- 230000001052 transient effect Effects 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/32—Transforming one recording into another or one representation into another
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
- G01V1/364—Seismic filtering
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
The invention provides a seismic data resolution improving method based on compressed sensing, which comprises the following steps: converting the seismic record into a frequency domain to obtain an amplitude spectrum of the seismic record, and decomposing the amplitude spectrum of the seismic record to obtain odd-even components of the amplitude spectrum of the seismic record; obtaining a secondary spectrum through the amplitude spectrum of the seismic record, and simulating the amplitude spectrum of the seismic wavelet according to the secondary spectrum; performing spectrum inversion based on the odd-even component of the seismic record amplitude spectrum and the amplitude spectrum of the seismic wavelet, and solving the stability of an inversion problem by adopting a compressive sensing method to obtain an accurate reflection coefficient sequence in the seismic record; and (3) performing convolution on the seismic wavelets of the broadband to obtain the seismic record with high resolution aiming at the reflection coefficient sequence in the seismic record. The method improves the resolution of seismic data, expands high-frequency components in the seismic data, improves the signal-to-noise ratio of the seismic data after the resolution is improved, has high thin reservoir identification capacity, and has good transverse continuity of the inverted seismic profile.
Description
Technical Field
The invention relates to the technical field of seismic data processing, in particular to a method and a device for improving resolution of seismic data based on compressed sensing.
Background
The resolution improvement processing of seismic data is an important research direction in geophysical exploration, and a plurality of scholars develop research aiming at the problem to form a series of seismic data frequency extension methods. Among them, the spectrum inversion has become a new method for effectively identifying the thin layer after years of development. The spectrum inversion is implemented by converting seismic data into a frequency domain, selecting a dominant frequency band of the data for inversion, and performing inversion in the frequency domain to solve a reflection coefficient. The reflection coefficient obtained by spectrum inversion has high resolution, and thin interbed of oil and gas reservoirs and the like can be effectively identified.
The study of spectral inversion starts earlier, and is mainly divided into several stages: partyka (1999) et al discovered that applying the spectral decomposition method to reservoir characterization, transforming seismic records to the frequency domain by discrete Fourier transform, describing transient thin layer thickness variations with amplitude spectra and discontinuity of lateral geological formations with phase spectra, can improve the resolution of seismic data to some extent. Marfurt and Kirlin (2001) add a time window function to perform seismic signal spectral decomposition to explain the tuning effect of the lamella and develop a new set of attributes to rapidly map the tuning effect of the lamella in three dimensions. Chorra et al (2006) perform spectral decomposition on seismic data, intercept local spectral data and perform inversion to obtain a reflection coefficient sequence, and the resolution of an inversion result is improved compared with that of an inversion result obtained by a conventional method. Puryear and Castagna (2008) carry out odd-even decomposition on the reflection coefficient, different weights are distributed to the reflection coefficient, an objective function of spectrum inversion is constructed, and the effectiveness of the method is verified through model simulation test and actual seismic data inversion. Chopira et al (2009) explain the application of relative acoustic impedance in thin layer reflection coefficient inversion by using examples, and acquire more stratum detail information through inversion, thereby more accurately explaining seismic data. Zhang and Castagna (2011) realize reflection coefficient inversion by establishing a dictionary representing a reflection coefficient model and constructing optimized seismic records by means of dictionary base recombination.
In the domestic field, some units and scholars are also conducting research in this field. Representative studies include: thin layer prediction and inversion method researches based on spectral inversion are earlier carried out in Qindvi, Han dynasty and the like (2009), a target function and a global inversion method are researched, the target function and the global inversion method comprise a random hill climbing inversion method, a Monte Carlo inversion method based on random search and the like, the spectral inversion of actual seismic data is carried out, and an inversion result with improved resolution is obtained. The research of the department of technology (2010) adopts a simulated annealing method to solve an inversion objective function, a reflection coefficient sequence is obtained, and the feasibility of the algorithm is verified through a series of algorithms. The Li Linxia, Tang Xiangrong Rong and the like (2012) research various spectral decomposition algorithm principles and simultaneously research the principle of a high-precision spectral inversion method, and perform actual seismic data spectral inversion. Liuwanjin and the like (2013) invert the three-dimensional seismic data volume of Daqing oil fields on the basis of the introduction spectrum inversion principle, and apply the seismic data volume with improved resolution obtained by inversion and reconstruction to seismic attribute analysis, thereby improving the attribute interpretation precision. Relevant research is carried out on a method for spectrum resolving and target function inversion by Zhou Jia Xiong and the like (2013), an inversion spectrum decomposition method based on L-BFGS is mainly researched, and the effectiveness of the method is proved through inversion processing of actual seismic data. Li Shachun, Li dong, etc. (2013) according to the characteristics of seismic record make analysis and research on adaptive optimization parameters, and adopt the weighted least square method to calculate deconvolution to improve seismic data resolution. Yuan-san-first (2013) tries a reflection coefficient inversion method based on sparse Bayesian learning, wherein reflection coefficient spectrum inversion and sparse Bayesian learning are combined, the number of non-zero elements of a reflection coefficient sequence is not preset, and the sparse reflection coefficient sequence is reconstructed by sequentially adding or deleting or re-estimating parameters. Sun Raymilk et al (2014) have studied reflection coefficient inversion, proposed some constraint methods of linear spectrum inversion, and inverted to sparse reflection coefficients under the Cauchy condition constraint.
At present, the information of the seismic wavelets needs to be known in a spectrum inversion implementation method, and the frequency, the phase and the like of the seismic wavelets have great influence on spectrum inversion. Also, in conventional spectral inversion processes, lateral variations in the seismic wavelets can cause the lateral continuity of the spectral inversion to deteriorate.
Therefore, the invention provides a method and a device for improving the resolution of seismic data based on compressed sensing.
Disclosure of Invention
In order to solve the above problems, the present invention provides a method for improving resolution of seismic data based on compressed sensing, the method comprising the following steps:
the method comprises the following steps: transforming the seismic record to a frequency domain to obtain an amplitude spectrum of the seismic record, and decomposing the amplitude spectrum of the seismic record to obtain an odd-even component of the amplitude spectrum of the seismic record;
step two: obtaining a secondary spectrum through the amplitude spectrum of the seismic record, and simulating the amplitude spectrum of the seismic wavelet according to the secondary spectrum;
step three: performing spectrum inversion based on the odd-even component of the seismic record amplitude spectrum and the amplitude spectrum of the seismic wavelet, and solving the stability of an inversion problem by adopting a compressive sensing method to obtain an accurate reflection coefficient sequence in the seismic record;
step four: and (3) performing convolution on the seismic wavelets of the broadband to obtain the seismic record with high resolution aiming at the reflection coefficient sequence in the seismic record, and finally obtaining seismic data with improved resolution.
According to an embodiment of the present invention, the step one specifically includes the following steps:
and transforming the seismic record to a frequency domain by utilizing fast Fourier transform to obtain an amplitude spectrum of the seismic record, and decomposing the amplitude spectrum of the seismic record according to the odd-even decomposition theory of signals to obtain the odd-even component of the amplitude spectrum of the seismic record.
According to an embodiment of the present invention, the second step specifically includes the following steps:
separating a low-frequency part and a high-frequency part of the secondary spectrum by a preset low-pass filter, eliminating the high-frequency part which reacts violent oscillation of the reflection coefficient amplitude spectrum, reserving the low-frequency part which reacts smooth wavelet amplitude energy, and taking the low-frequency part as the amplitude spectrum of the seismic wavelet.
According to an embodiment of the invention, the preset low-pass filter is represented by the following formula:
wherein f is frequency, f c J is an imaginary unit, n is a constant, and A is the amplitude value and is a function of f.
According to an embodiment of the present invention, in the third step, a gradient projection sparse reconstruction algorithm in compressed sensing is selected to solve an objective function, and the objective function is represented by the following formula:
setting:
r＝u-v(u≥0,v>0)
and r is L 1 Norm:
the objective function is updated as:
further:
obtaining the target function:
wherein R is a reflection coefficient, R is a seismic record, and G is a seismic wavelet matrix.
According to one embodiment of the invention, the objective function is subjected to a compressed sensing-based inversion process by:
step a, initializing and setting parameters in the objective function;
b, calculating the proportion alpha of odd components in the reflection coefficient sequence 0 And mid (. alpha.) is used in combination min ,α 0 ,α max ) Substitute for alpha 0 ；
Step c, in the sequence { alpha 0 ,βα 0 ,β 2 α 0 … } and assigning a first value that holds the predetermined formula to alpha (k) ；
And d, performing convergence judgment, terminating iteration when the result meets the precision requirement, and if the result does not meet the precision requirement, changing k to k +1 and returning to the step b.
According to an embodiment of the present invention, the preset formula is expressed as:
in addition:
wherein α is a gradient change factor.
According to one embodiment of the invention, the broadband seismic wavelets are represented by the formula:
further obtaining:
where g represents a parameter of the Rake wavelet frequency, and p and q represent the integration limits of the parameter g.
According to another aspect of the invention, there is also provided a storage medium containing a series of instructions for carrying out the steps of the method as described in any one of the above.
According to another aspect of the present invention, there is also provided a compressed sensing-based seismic data resolution enhancement apparatus for performing seismic data resolution enhancement processing by a method as described in any one of the above, the apparatus comprising:
the device comprises a first module, a second module and a third module, wherein the first module is used for transforming the seismic record to a frequency domain to obtain an amplitude spectrum of the seismic record, and decomposing the amplitude spectrum of the seismic record to obtain an odd-even component of the amplitude spectrum of the seismic record;
the second module is used for solving a secondary spectrum through the amplitude spectrum of the seismic record and simulating the amplitude spectrum of the seismic wavelet according to the secondary spectrum;
the third module is used for carrying out spectrum inversion based on the odd-even component of the seismic record amplitude spectrum and the amplitude spectrum of the seismic wavelet, and solving the stability of an inversion problem by adopting a compressed sensing method to obtain an accurate reflection coefficient sequence in the seismic record;
and the fourth module is used for performing convolution on the seismic wavelets of the broadband to obtain the seismic record with high resolution aiming at the reflection coefficient sequence in the seismic record, and finally obtaining the seismic data with improved resolution.
The method and the device for improving the resolution of the seismic data based on the compressive sensing can well improve the resolution of seismic data, expand high-frequency components in the seismic data, improve the signal-to-noise ratio of the seismic data after the resolution is improved, have high thin reservoir identification capacity, and have good transverse continuity of an inverted seismic profile.
Additional features and advantages of the invention will be set forth in the description which follows, and in part will be obvious from the description, or may be learned by practice of the invention. The objectives and other advantages of the invention will be realized and attained by the structure particularly pointed out in the written description and claims hereof as well as the appended drawings.
Drawings
The accompanying drawings, which are included to provide a further understanding of the invention and are incorporated in and constitute a part of this specification, illustrate embodiments of the invention and together with the description serve to explain the principles of the invention and not to limit the invention. In the drawings:
FIG. 1 shows a flow diagram of a compressed sensing-based seismic data resolution enhancement method according to one embodiment of the invention;
FIG. 2 illustrates a flow diagram of a compressed sensing-based seismic data resolution enhancement method according to another embodiment of the invention;
FIG. 3 is a graph showing a quadratic spectrum comparison of a synthetic seismic record with a known wavelet and an amplitude spectrum comparison of the derived seismic wavelet with the known wavelet in the absence of noise according to one embodiment of the present invention;
FIG. 4 is a graph showing a quadratic spectrum comparison of a synthetic seismic record with a known wavelet and a comparison of the amplitude spectrum of the derived seismic wavelet with the known wavelet for a SNR of 10 in accordance with one embodiment of the present invention;
FIG. 5 is a graph showing a secondary spectrum comparison of a synthetic seismic record with a known wavelet and an amplitude spectrum comparison of the derived seismic wavelet with the known wavelet for an SNR of 5 in accordance with one embodiment of the present invention;
FIG. 6 shows a schematic diagram of a pre-defined low pass filter according to an embodiment of the invention;
FIG. 7 shows a graph comparing a reflection coefficient model with inversion results according to an embodiment of the invention;
FIG. 8 shows a synthetic seismic record versus inversion results according to an embodiment of the invention;
FIG. 9 shows a schematic of raw seismic data according to one embodiment of the invention;
FIG. 10 shows a reflection coefficient profile from spectral inversion according to one embodiment of the invention;
FIG. 11 shows a schematic of a broadband seismic wavelet according to one embodiment of the present invention;
FIG. 12 shows an amplitude spectrum of a broadband seismic sub-wave in accordance with an embodiment of the invention;
FIG. 13 shows a schematic diagram of the improved resolution results of spectral inversion according to one embodiment of the invention;
FIG. 14 shows an amplitude spectrum contrast plot of a raw seismic record and a spectrum inversion seismic record according to one embodiment of the invention; and
FIG. 15 is a block diagram of a compressed sensing-based seismic data enhanced resolution device according to an embodiment of the invention.
Detailed Description
In order to make the objects, technical solutions and advantages of the present invention more apparent, the following embodiments of the present invention are further described in detail with reference to the accompanying drawings.
In the spectrum inversion in the prior art, the interference condition of odd and even reflection coefficients is mainly utilized to effectively improve the resolution of seismic data and achieve the purpose of identifying thin layers. In the time domain, the seismic reflection coefficient may be expressed as:
in the formula: r is a radical of hydrogen 1 ，r 2 Expressed as the formation reflection coefficients at different times; δ (t) is a pulse function. According to the property of the pulse function, Fourier transform (Fourier transform) is performed on the above formula (1) to obtain:
it is clear that:
Re[R(f)]＝2r e cos(πfT) (3)
Im[R(f)]＝2r o sin(πfT) (4)
in the formula: r is e ,r o Respectively, a reflection coefficient pair (r) 1 ,r 2 ) Even and odd components.
According to the seismic record convolution model, under the condition of no noise, the seismic record s (t), the seismic wavelet omega (t) and the reflection coefficient r (t) satisfy the following conditions:
in the formula: "+" indicates convolution operation; "×" indicates product operation; FT is Fourier positive transformation; s (f) is a seismic recording spectrum; therefore, the spectrum inversion removes the influence of the seismic wavelets from the seismic records according to the seismic records and the frequency spectrum information of the seismic wavelets to obtain a reflection coefficient sequence, and the spectrum inversion target function is as follows:
to better embody the weight of the parity inversion, the above equation is decomposed into a matrix in parity form:
in the formula: re [. C]And Im [ ·]Respectively representing the real part and the imaginary part of the frequency spectrum; f. of i Represents the sampling frequency; alpha is alpha e And alpha o Representing the ratio of even and odd components of the reflection coefficient. The frequency spectrums of the seismic records and the seismic wavelets are obtained by directly carrying out Fourier transform on the frequency spectrums, and the real part and the imaginary part of the reflection coefficient frequency spectrum can be obtained by calculation through a derivation formula.
Equations (7) and (8) are abbreviated as a minimized objective function, and a standard form can be obtained by adding a reflection coefficient sparse constraint:
finally, the above-mentioned constrained optimization problem can be rewritten into a final objective function (unconstrained optimization problem) by the lagrange multiplier method:
in the formula: λ represents the regularization parameter used to balance the weight of the result to the quadratic and sparse terms.
As described above, in the prior art, the information of the seismic wavelet needs to be known in the spectrum inversion implementation method, and the frequency, the phase, and the like of the seismic wavelet have a large influence on the spectrum inversion. Also, in conventional spectral inversion processes, the lateral variation of the seismic wavelets can cause the lateral continuity of the spectral inversion to deteriorate.
Based on the problems in the spectrum inversion of the prior art, the invention provides a spectrum inversion algorithm based on a compressed sensing method, and adopts secondary spectrum based on seismic data to invert so as to improve the resolution of seismic data. The seismic data processed by improving the resolution ratio provides more reliable information for identification and detection of thin reservoirs and small target geologic bodies, and provides technical support for prediction and description of high-precision reservoirs.
FIG. 1 shows a flow diagram of a compressed sensing-based seismic data resolution enhancement method according to one embodiment of the invention.
As shown in fig. 1, in step S101, the seismic record is transformed into the frequency domain to obtain an amplitude spectrum of the seismic record, and the amplitude spectrum of the seismic record is decomposed to obtain parity components of the amplitude spectrum of the seismic record.
In one embodiment, in step S101, the seismic records are transformed to the frequency domain using a fast Fourier transform to obtain the amplitude spectrum of the seismic records, and the amplitude spectrum of the seismic records is decomposed to obtain the parity components of the amplitude spectrum of the seismic records according to the parity decomposition theory of the signals.
As shown in fig. 1, in step S102, a secondary spectrum is obtained from the amplitude spectrum of the seismic record, and the amplitude spectrum of the seismic wavelet is simulated based on the secondary spectrum.
In particular, the spectrum of the seismic wavelets is required in the spectrum inversion process. In step S102, the invention provides a method for simulating a seismic sub-wave spectrum by using a secondary spectrum of seismic channel data, and solves the defect that seismic wavelets need to be estimated in conventional spectrum inversion, so that the transverse continuity of the seismic data can be well improved.
Further, the quadratic spectrum method is based on the assumption of spectral modeling: the amplitude spectrum of the seismic wavelet is smooth or much smoother than that of the reflection coefficient, so that the amplitude spectrum of the seismic wavelet and the reflection coefficient sequence is necessarily separable in the frequency domain. In the aspect of frequency analysis, the low-frequency part of the seismic record secondary spectrum reflects the energy of a smooth wavelet amplitude spectrum, and the high-frequency component reflects the energy of a violent oscillation reflection coefficient sequence amplitude spectrum, so that the low-frequency part and the high-frequency component are separated by adopting a preset low-pass filter. The low-pass filter is preset to reject short-term fluctuations in the signal (sharp oscillations in the reflection coefficient amplitude spectrum) and preserving the long-term trend (smooth sub-spectrum) is a means of smoothing the signal.
In one embodiment, in step S102, a low frequency part and a high frequency part of the quadratic spectrum are separated by a predetermined low pass filter, the high frequency part reflecting the violent oscillation of the reflection coefficient amplitude spectrum is eliminated, the low frequency part reflecting smooth wavelet amplitude energy is retained, and the low frequency part is used as the amplitude spectrum of the seismic wavelet.
Further, the preset low-pass filter is expressed by the following formula:
wherein f is frequency, f c Is the cut-off frequency of the low-pass filtering, j is an imaginary unit, n is a constant, A is the amplitude value,as a function of f.
As shown in fig. 1, in step S103, spectrum inversion is performed based on the parity component of the seismic record amplitude spectrum and the amplitude spectrum of the seismic wavelet, and a compressive sensing method is used to solve the stability of the inversion problem, so as to obtain an accurate reflection coefficient sequence in the seismic record.
In one embodiment, in step S103, a gradient projection sparse reconstruction algorithm in compressed sensing is selected to solve the objective function, and the objective function is represented by the following formula:
setting:
r＝u-v(u≥0,v>0)
and r is L 1 Norm:
the objective function is updated as:
further:
obtaining an objective function:
wherein R is a reflection coefficient, R is a seismic record, and G is a seismic wavelet matrix.
Specifically, step S103 implements a compressed sensing-based inversion algorithm to perform spectrum inversion. Solving sparse reflectances in spectral inversionNumber series, i.e. with respect to L 1 And (3) a norm minimization problem solving algorithm, namely solving the formula (10) by selecting a gradient projection sparse reconstruction algorithm in compressed sensing:
the formula (10) is modified as follows:
the reflection coefficient sequence r to be inverted is decomposed into a positive number sequence u and a negative number sequence-v, and then r can be written as r-u-v (u is more than or equal to 0, v is more than or equal to>0) Thus L of r 1 The norm can be written as:
in the formula: l n ＝[1,1…,1] 1×n If equation (12) is introduced into equation (11), equation (11) is solved by solving the following objective function:
the following transformations are further performed:
obtaining a new target function expression form:
the solution for the spectral inversion is obtained by solving the new objective function (15).
Further, an inversion process based on compressed sensing is performed on the target function formula (15) through the following steps, and a solution of formula (15) is finally obtained through iteration of the following steps:
step a, initializing and setting parameters in the objective function. In particular toTo say, given z (0) λ is 0.01, and a parameter β ∈ (0,1) is selected, μ ∈ (0,0.5), and k is 0.
Step b, calculating the proportion alpha of the odd component in the reflection coefficient sequence 0 And mid (. alpha.) is used in combination min ,α 0 ,α max ) Substituted for alpha 0 . Specifically, α is calculated by using equations (16) and (17) 0 。
Step c, in the sequence { alpha 0 ,βα 0 ,β 2 α 0 … } and assigning a first value that holds the predetermined formula to alpha (k) And make an order
Specifically, the preset formula is expressed as:
wherein α is a gradient change factor.
And d, performing convergence judgment, terminating iteration when the result meets the precision requirement, and if the result does not meet the precision requirement, changing k to k +1 and returning to the step b.
As shown in fig. 1, in step S104, a seismic record with high resolution is obtained by convolution of broadband seismic wavelets for the reflection coefficient sequence in the seismic record, and finally, seismic data with improved resolution is obtained.
In one embodiment, in step S104, a broadband seismic wavelet is represented by the following formula:
further obtaining:
where g represents a parameter of the Rake wavelet frequency, and p and q represent the integration limits of the parameter g.
FIG. 2 shows a flow diagram of a compressed sensing-based seismic data resolution enhancement method according to another embodiment of the invention.
As shown in fig. 2, the seismic record is first transformed into the frequency domain by FFT to obtain the amplitude spectrum of the seismic record, and the amplitude spectrum of the seismic record is decomposed according to the parity decomposition theory of the signal to obtain the parity component of the amplitude spectrum of the seismic record.
As shown in fig. 2, the amplitude spectrum of the seismic wavelet is then simulated by obtaining a quadratic spectrum from the amplitude spectrum obtained from the input seismic record. From a convolution model of seismic data, the quadratic spectrum of a seismic record can be well approximated as the amplitude spectrum of a seismic wavelet.
As shown in fig. 2, spectrum inversion is performed by using the odd-even component of the seismic record amplitude spectrum and the seismic wavelet amplitude spectrum, and a compressed sensing method is used for solving the stability of the inversion problem in the inversion process, so that the compressed sensing spectrum inversion is realized to obtain an accurate reflection coefficient sequence in the seismic record.
As shown in fig. 2, the seismic records with high resolution are obtained by convolution of a broadband seismic wavelet for the reflection coefficient sequence obtained by inversion, and finally the seismic data with improved resolution are obtained.
Fig. 3-5 show graphs of a quadratic spectrum comparison of a 5 SNR, a 10 SNR synthetic seismic record with a known wavelet and an amplitude spectrum comparison of the derived seismic wavelet with the known wavelet, respectively, in the absence of noise, according to one embodiment of the invention. The left hand side of FIGS. 3-5 are secondary spectra of synthetic seismic records compared to known wavelets. The right hand graphs of FIGS. 3-5 are amplitude spectra of the computed seismic wavelets compared to known wavelets.
FIG. 3 shows a process of obtaining wavelet amplitude spectrum of synthetic seismic record by using the quadratic spectrum method provided by the invention, the invention designs an unequal interval reflection coefficient sequence with sampling interval of 0.001s, number of sampling points of 474 and number of reflection coefficients of 17, and the unequal interval reflection coefficient sequence is subjected to convolution with a Ricker wavelet (Ricker) with main frequency of 30Hz to obtain the synthetic seismic record.
Fig. 3-5 respectively obtain the secondary spectrum for the synthetic seismic record without noise, with SNR 5 and SNR 10, select a proper position to perform low-pass filtering on the secondary spectrum of the synthetic seismic record, respectively perform normalization processing on the actual wavelet and the obtained result, and finally respectively compare the actual wavelet with the amplitude spectrum of the actual wavelet, so that the wavelet obtained by using the secondary spectrum method provided by the invention can be well matched with the amplitude spectrum of the actual wavelet, and the requirement of use precision is met.
FIG. 6 shows a schematic diagram of a default low pass filter according to an embodiment of the invention.
In one embodiment, the present invention employs a Butterworth (Butterworth) low pass filter having the formula:
wherein: n is an integer and is called as the order of the filter, the larger n is, the better the approximation of the pass band and the stop band is, and the steeper the transition band is; f. of c Is the cut-off frequency. After optimization calculation, the Butterworth low-pass filter with n-3 has the best effect.
FIG. 7 shows a reflection coefficient model compared to inversion results according to one embodiment of the invention. The spectrum inversion is carried out on the synthetic seismic records with different noise interferences by using a compressed sensing algorithm (gradient projection sparse reconstruction algorithm), and the obtained reflection coefficient is compared with the original reflection coefficient, so that the spectrum inversion has good noise immunity and good stability. FIG. 7 is a) -b) -c) -d) from top to bottom, wherein a) is a sparse reflectance model; b) obtaining the SNR (signal to noise ratio) 10 reflection coefficient inversion result; c) obtaining the SNR (signal to noise ratio) 5 reflection coefficient inversion result; d) the SNR is 2 reflection coefficient inversion result.
FIG. 8 shows synthetic seismic records compared to inversion results according to one embodiment of the invention. The reflection coefficient and Ricker wavelet (f) obtained by the spectrum inversion method provided by the invention s 30Hz) to obtain the calculated seismic records, and comparing the seismic records. The seismic record obtained by the spectrum inversion method provided by the invention has higher inosculation degree with the actual synthetic seismic record, and the inversion precision is higher. FIG. 8 shows a) -b) -c) -d) from top to bottom, wherein a) is a synthetic seismic record; b) the SNR is 10 inversion result; c) the SNR is equal to 5 inversion result; d) the SNR is 2 inversion result.
FIG. 9 shows a schematic of raw seismic data according to one embodiment of the invention. Spectrum inversion of actual seismic data, applying a spectrum inversion method based on a compressed sensing algorithm to the actual seismic data (fig. 9), performing convolution on a reflection coefficient profile (fig. 10) obtained by calculating the actual data and a broadband wavelet (fig. 11), wherein a band-pass wavelet can be regarded as a result of harmonic synthesis of each frequency in a pass band, Ricker wavelets with different frequency parameters are synthesized into a novel wavelet-Shu's wavelet by the professor of Shu Shong, and the calculation formula is as follows:
further obtaining:
in the formula: g is a parameter characterizing the frequency of the Ricker wavelet; p and q are the integration limits of the parameter g in Shu's wavelet (20). Fig. 11 and 12 are graphs of the time domain shu wavelet and its amplitude spectrum with parameters p being 12 and q being 150, respectively.
The invention uses Shu's son broadband wavelets (figure 11), and finally obtains the result of spectrum inversion for improving resolution (figure 12). As can be seen, the in-phase axis of the optical fiber is thin, the inverted fault is obvious, and the details are highlighted. While the lateral continuity is also good.
FIG. 14 shows a graph of amplitude spectra comparison of an original seismic record and a spectrum inverted seismic record, in accordance with one embodiment of the invention. After the spectrum inversion is carried out on the original seismic record, the amplitude spectrum of the resolution-improved result of the original seismic record and the spectrum inversion are respectively calculated, and it can be seen that after the spectrum inversion method provided by the invention is carried out, the frequency band is obviously widened, and the effectiveness and the superiority of the spectrum inversion for improving the resolution are proved.
In summary, the invention successfully applies the compressed sensing method to the seismic data resolution enhancement processing technology, and simultaneously combines the secondary spectrum estimation seismic wavelet method using the seismic channels in the implementation process to construct the seismic data resolution enhancement processing technology which can combine the compressed sensing method and the spectrum inversion; the method can well solve the problem of transverse discontinuity of the existing spectrum inversion in the process of improving the resolution of seismic data by utilizing the technology, and effectively invert to obtain a seismic reflection coefficient profile and a seismic profile after the resolution is improved. The seismic data processed by improving the resolution ratio has higher resolution ratio, more reliable information can be provided for the subsequent identification and detection of thin reservoirs and small target geologic bodies, technical support is provided for the prediction and description of high-precision reservoirs, and the risk of loss of interest is reduced in oil and gas exploration drilling.
In order to verify the inversion effect of the spectrum inversion problem to the two models by the method provided by the invention, a one-dimensional coefficient pulse coefficient model shown as a) in fig. 7 is established, the sampling interval is 0.002s, 17 groups of reflection coefficients are totally arranged in the model, 474 sampling points are totally arranged, the synthetic seismic record obtained by convolution of the reflection coefficient model and a Ricker wavelet (Ricker) with the main frequency of 30Hz is shown as a) in fig. 8, and the noise with different energy is added into the synthetic seismic record to analyze the influence degree of the noise on the spectrum inversion, so as to obtain the synthetic record of different noises.
The process of the spectrum inversion method is a 'wavelet removing' process, the influence of seismic wavelets on an inversion result is removed, and the seismic wavelets are extracted before the spectrum inversion calculation. Wavelet extraction was first performed on the synthetic recordings with different noise, resulting in the comparison of results shown in fig. 3-5. It can be seen that: the extracted wavelet amplitude spectrum has better goodness of fit with the actual wavelet amplitude spectrum, and the influence of noise interference on wavelet extraction is not large.
In order to extract wavelets more accurately, the present invention applies a butterworth low pass filter (fig. 6), and after many experiments, the butterworth low pass filter with n-3 hours has the best effect.
B) -d) in FIG. 7 are the results of the reflection coefficients under different intensity noise interferences, respectively, compared with the convolution data of the known seismic wavelets (b-d in FIG. 8). The inversion result verifies that the method has good anti-noise effect, the position and the amplitude of the reflection coefficient can be matched more accurately, and the inversion result is stable and accurate.
In order to further analyze the effect of the application of the actual data based on the spectrum inversion of the compressed sensing, the method is applied to the processing of the actual seismic data (figure 9), and the inverted reflection coefficient is shown in figure 10. The resolution-improved seismic section is obtained by convolution of the reflection coefficient section obtained by using the Shu's broadband wavelet and the spectrum inversion as shown in FIG. 11 (FIG. 13). Comparing fig. 9 with fig. 13, it can be seen that compared with the original seismic profile, the seismic profile obtained by inversion has higher resolution, thin layers and small-scale geological bodies which are indistinguishable and have a thickness smaller than the tuned thickness on the original seismic profile are identified, and high-resolution and richer geological information is provided.
In summary, the spectrum inversion method provided by the invention is characterized in that the spectrum is inverted in the frequency domain by transforming the seismic data to the frequency domain and selecting the dominant frequency band of the data for inversion, and the reflection coefficient is solved in the frequency domain by inversion. The reflection coefficient obtained by spectrum inversion has high resolution, and thin interbed of oil and gas reservoirs and the like can be effectively identified. The invention realizes the compressed sensing theory, adopts random sampling to obtain discrete samples of signals, and then reconstructs the spectrum inversion method of the signals through a nonlinear reconstruction algorithm. The practical data application result shows that the spectrum inversion is carried out by utilizing the compressed sensing, the resolution ratio of seismic data can be well improved, and basic data are provided for high-precision reservoir prediction. The invention has the following advantages:
(1) the reflection coefficient profile or reflection coefficient volume can be inverted by using only part of the seismic data, and the resolution is far higher than that of the original seismic data.
(2) No other constraint conditions such as logging information are needed, no initial model assumption is needed, no transverse mixing is needed, no assumption is needed on the position and amplitude of the reflection coefficient, and no reflection coefficient spectrum assumption is needed.
(3) The compressed sensing algorithm used by the invention has strong noise immunity, and the result is more stable and accurate through the inversion test of the seismic data with different intensity noise interferences. Meanwhile, the processing cost is low according to the characteristics of the compressed sensing algorithm.
(4) The seismic wavelet of the frequency domain is solved by utilizing the quadratic spectrum of the seismic data, the influence of the seismic wavelet on spectrum inversion is solved, and the spectrum inversion result is more accurate. Meanwhile, the method expands the length of the seismic data to 2 times of the length of the original seismic record, and solves the problem that the reflection coefficient must be even when the reflection coefficient is decomposed in odd-even mode.
According to another aspect of the invention, there is also provided a storage medium containing a series of instructions for performing the steps of the method for improving resolution of seismic data based on compressed sensing as described in any one of the above.
FIG. 15 is a block diagram of a compressed sensing-based seismic data enhanced resolution device according to an embodiment of the invention. As shown in fig. 15, the compressed sensing-based seismic data resolution enhancement apparatus 1500 includes a first module 1501, a second module 1502, a third module 1503, and a fourth module 1504.
The first module 1501 is configured to transform the seismic records to the frequency domain to obtain amplitude spectra of the seismic records, and decompose the amplitude spectra of the seismic records to obtain parity components of the amplitude spectra of the seismic records. The second module 1502 is configured to obtain a secondary spectrum from the amplitude spectrum of the seismic record, and simulate the amplitude spectrum of the seismic wavelet according to the secondary spectrum. The third module 1503 is configured to perform spectrum inversion based on the parity component of the seismic record amplitude spectrum and the amplitude spectrum of the seismic wavelet, and perform stability solution of an inversion problem by using a compressive sensing method to obtain an accurate reflection coefficient sequence in the seismic record. The fourth module 1504 is configured to perform convolution on the broadband seismic wavelets for the reflection coefficient sequence in the seismic record to obtain a seismic record with high resolution, and finally obtain seismic data with improved resolution.
In conclusion, the method and the device for improving the resolution of the seismic data based on the compressed sensing can well improve the resolution of the seismic data, expand high-frequency components in the seismic data, improve the signal-to-noise ratio of the seismic data after the resolution is improved, have high thin reservoir identification capacity, and have good transverse continuity of the inverted seismic profile.
It is to be understood that the disclosed embodiments of the invention are not limited to the particular structures, process steps, or materials disclosed herein but are extended to equivalents thereof as would be understood by those ordinarily skilled in the relevant arts. It is also to be understood that the terminology used herein is for the purpose of describing particular embodiments only, and is not intended to be limiting.
Reference in the specification to "one embodiment" or "an embodiment" means that a particular feature, structure, or characteristic described in connection with the embodiment is included in at least one embodiment of the invention. Thus, the appearances of the phrase "one embodiment" or "an embodiment" in various places throughout this specification are not necessarily all referring to the same embodiment.
Although the embodiments of the present invention have been described above, the above description is only for the convenience of understanding the present invention, and is not intended to limit the present invention. It will be understood by those skilled in the art that various changes in form and details may be made therein without departing from the spirit and scope of the invention as defined by the appended claims.
Claims (10)
1. A method for enhancing resolution of seismic data based on compressed sensing, the method comprising the steps of:
the method comprises the following steps: transforming the seismic record to a frequency domain to obtain an amplitude spectrum of the seismic record, and decomposing the amplitude spectrum of the seismic record to obtain odd-even components of the amplitude spectrum of the seismic record;
step two: obtaining a secondary spectrum through the amplitude spectrum of the seismic record, and simulating the amplitude spectrum of the seismic wavelet according to the secondary spectrum;
step three: carrying out spectrum inversion based on the odd-even component of the seismic record amplitude spectrum and the amplitude spectrum of the seismic wavelet, and solving the stability of an inversion problem by adopting a compressive sensing method to obtain an accurate reflection coefficient sequence in the seismic record;
step four: and (3) performing convolution on the seismic wavelets of the broadband to obtain the seismic record with high resolution aiming at the reflection coefficient sequence in the seismic record, and finally obtaining seismic data with improved resolution.
2. The method for improving the resolution of seismic data based on compressed sensing of claim 1, wherein the first step specifically comprises the steps of:
and transforming the seismic record to a frequency domain by utilizing fast Fourier transform to obtain an amplitude spectrum of the seismic record, and decomposing the amplitude spectrum of the seismic record according to the odd-even decomposition theory of signals to obtain the odd-even component of the amplitude spectrum of the seismic record.
3. The method for improving the resolution of seismic data based on compressed sensing according to claim 1, wherein the second step specifically comprises the following steps:
separating a low-frequency part and a high-frequency part of the secondary spectrum by a preset low-pass filter, eliminating the high-frequency part which reacts violent oscillation of a reflection coefficient amplitude spectrum, reserving the low-frequency part which reacts smooth wavelet amplitude energy, and taking the low-frequency part as the amplitude spectrum of the seismic wavelet.
4. The compressed sensing-based seismic data up-resolution method of claim 3, wherein the pre-set low-pass filter is represented by the formula:
wherein f is frequency, f c The cut-off frequency of the low-pass filter is j in an imaginary unit, n is a constant, and A is the amplitude value and is a function of f.
5. The method for improving the resolution of the seismic data based on the compressed sensing as claimed in claim 1, wherein in the third step, a gradient projection sparse reconstruction algorithm in the compressed sensing is selected to solve the objective function, and the objective function is represented by the following formula:
setting:
r＝u-v(u≥0,v＞0)
and r is L 1 Norm:
the objective function is updated as:
further:
obtaining the target function:
wherein R is a reflection coefficient, R is a seismic record, and G is a seismic wavelet matrix.
6. The compressed sensing-based seismic data up-resolution method of claim 5, wherein the objective function is subjected to a compressed sensing-based inversion process by:
step a, initializing and setting parameters in the objective function;
b, calculating the proportion alpha of odd components in the reflection coefficient sequence 0 And mid (. alpha.) is used in combination min ,α 0 ,α max ) Substitute for alpha 0 ；
Step c, in the sequence { alpha 0 ,βα 0 ,β 2 α 0 … } and assigning a first value that holds the predetermined formula to alpha (k) ；
And d, performing convergence judgment, terminating iteration when the result meets the precision requirement, and if the result does not meet the precision requirement, changing k to k +1 and returning to the step b.
7. The compressed sensing-based seismic data resolution enhancement method of claim 6, wherein the predetermined formula is expressed as:
in addition:
wherein α is a gradient change factor.
8. The compressed sensing-based seismic data resolution enhancement method of claim 1, wherein the broadband seismic wavelets are represented by the formula:
further obtaining:
where g represents a parameter of the Rake wavelet frequency, and p and q represent the integration limits of the parameter g.
9. A storage medium characterized in that it contains a series of instructions for carrying out the steps of the method according to any one of claims 1 to 8.
10. An apparatus for resolution enhancement of seismic data based on compressed sensing by the method of any of claims 1-8, the apparatus comprising:
the device comprises a first module, a second module and a third module, wherein the first module is used for transforming the seismic record to a frequency domain to obtain an amplitude spectrum of the seismic record, and decomposing the amplitude spectrum of the seismic record to obtain an odd-even component of the amplitude spectrum of the seismic record;
the second module is used for solving a secondary spectrum through the amplitude spectrum of the seismic record and simulating the amplitude spectrum of the seismic wavelet according to the secondary spectrum;
the third module is used for carrying out spectrum inversion based on the odd-even component of the seismic record amplitude spectrum and the amplitude spectrum of the seismic wavelet, and solving the stability of an inversion problem by adopting a compressed sensing method to obtain an accurate reflection coefficient sequence in the seismic record;
and the fourth module is used for performing convolution on the seismic wavelets of the broadband to obtain the seismic record with high resolution aiming at the reflection coefficient sequence in the seismic record, and finally obtaining the seismic data with improved resolution.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110285637.XA CN115113265A (en) | 2021-03-17 | 2021-03-17 | Seismic data resolution improving method and device based on compressed sensing |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110285637.XA CN115113265A (en) | 2021-03-17 | 2021-03-17 | Seismic data resolution improving method and device based on compressed sensing |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN115113265A true CN115113265A (en) | 2022-09-27 |
Family
ID=83323906
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110285637.XA Pending CN115113265A (en) | 2021-03-17 | 2021-03-17 | Seismic data resolution improving method and device based on compressed sensing |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN115113265A (en) |
-
2021
- 2021-03-17 CN CN202110285637.XA patent/CN115113265A/en active Pending
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Gómez et al. | A simple method inspired by empirical mode decomposition for denoising seismic data | |
| Cicone et al. | Numerical analysis for iterative filtering with new efficient implementations based on FFT | |
| Chen et al. | Double-sparsity dictionary for seismic noise attenuation | |
| Neelamani et al. | Coherent and random noise attenuation using the curvelet transform | |
| CN102221708B (en) | Fractional-Fourier-transform-based random noise suppression method | |
| Dong et al. | Signal-to-noise ratio enhancement for 3C downhole microseismic data based on the 3D shearlet transform and improved back-propagation neural networks | |
| Lorentzen et al. | History matching the full Norne field model using seismic and production data | |
| CN103954992B (en) | Deconvolution method and device | |
| CN112630835A (en) | High-resolution post-stack seismic wave impedance inversion method | |
| CN102928875B (en) | Wavelet extraction method based on fractional number order Fourier | |
| CN113077386A (en) | Seismic data high-resolution processing method based on dictionary learning and sparse representation | |
| Qian et al. | Unsupervised seismic footprint removal with physical prior augmented deep autoencoder | |
| CN110687597B (en) | Wave impedance inversion method based on joint dictionary | |
| Bai et al. | Seismic signal enhancement based on the low‐rank methods | |
| Zhong et al. | Multiscale residual pyramid network for seismic background noise attenuation | |
| CN112213782B (en) | Processing method and device for sub-phase seismic data and server | |
| Yang et al. | Deep learning with fully convolutional and dense connection framework for ground roll attenuation | |
| Bai et al. | A U-Net based deep learning approach for seismic random noise suppression | |
| CN115113265A (en) | Seismic data resolution improving method and device based on compressed sensing | |
| CN111766624B (en) | Seismic data frequency extension processing method and device, storage medium and electronic equipment | |
| CN107678065B (en) | The guarantor for improving seismic resolution constructs well control space the Method of Deconvolution and device | |
| CN115236733A (en) | DAS-VSP data background noise suppression method based on deep learning | |
| US12050295B2 (en) | Bandwith extension of geophysical data | |
| Liu et al. | Adaptive time-reassigned synchrosqueezing transform for seismic random noise suppression | |
| CN108693558A (en) | Seismic data processing technique and device |
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 |