CN114779335B - Elastic wave direct envelope inversion method based on anisotropic total variation constraint - Google Patents
Elastic wave direct envelope inversion method based on anisotropic total variation constraint Download PDFInfo
- Publication number
- CN114779335B CN114779335B CN202210336492.6A CN202210336492A CN114779335B CN 114779335 B CN114779335 B CN 114779335B CN 202210336492 A CN202210336492 A CN 202210336492A CN 114779335 B CN114779335 B CN 114779335B
- Authority
- CN
- China
- Prior art keywords
- wave
- longitudinal
- envelope
- model
- total variation
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 46
- 238000000354 decomposition reaction Methods 0.000 claims abstract description 8
- 238000005457 optimization Methods 0.000 claims description 6
- 238000004458 analytical method Methods 0.000 claims description 4
- 238000012937 correction Methods 0.000 claims description 3
- 230000010354 integration Effects 0.000 claims description 3
- 230000003068 static effect Effects 0.000 claims description 3
- 238000002945 steepest descent method Methods 0.000 claims description 3
- 238000001514 detection method Methods 0.000 claims description 2
- 150000003839 salts Chemical class 0.000 description 8
- 230000000694 effects Effects 0.000 description 7
- 230000006870 function Effects 0.000 description 6
- 238000004364 calculation method Methods 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 238000003384 imaging method Methods 0.000 description 2
- 238000012545 processing Methods 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 238000001228 spectrum Methods 0.000 description 2
- 239000004215 Carbon black (E152) Substances 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000013135 deep learning Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000006073 displacement reaction Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000001914 filtration Methods 0.000 description 1
- 229930195733 hydrocarbon Natural products 0.000 description 1
- 150000002430 hydrocarbons Chemical class 0.000 description 1
- 230000000877 morphologic effect Effects 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 230000009919 sequestration 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
- G01V1/30—Analysis
- G01V1/303—Analysis for determining velocity profiles or travel times
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
- G01V2210/622—Velocity, density or impedance
- G01V2210/6222—Velocity; travel time
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 relates to an elastic wave direct envelope inversion method based on anisotropic total variation constraint, which is to obtain high-precision longitudinal wave velocity structure and transverse wave velocity structure of an elastic strong scattering medium through the anisotropic total variation constraint and the elastic wave direct envelope inversion. Firstly, carrying out wave field mode decomposition on an elastic wave field to respectively obtain longitudinal wave fields and transverse wave fields, and calculating envelope fields to respectively obtain forward longitudinal wave envelope fields and transverse wave envelope fields; then, according to the gradient expression, longitudinal and transverse wave velocity gradients of elastic wave direct envelope inversion can be obtained, and a longitudinal and transverse wave velocity model can be updated; then, applying anisotropic total variation constraint to the current longitudinal and transverse wave velocity models to obtain constrained longitudinal and transverse wave velocity update models; and finally, taking the direct envelope inversion result of the anisotropic total variation constraint elastic wave as an initial model, and carrying out the full waveform inversion of the anisotropic total variation constraint elastic wave to obtain a high-precision longitudinal wave velocity model of the elastic strong scattering medium.
Description
Technical Field
The invention relates to a method for applying anisotropic total variation constraint in elastic wave direct envelope inversion to improve inversion effect of longitudinal and transverse wave velocity structures of a strong scattering medium.
Background
Strong scattering media are generally capable of forming good sequestration media for hydrocarbon resources, and therefore high-precision parametric modeling methods for strong scattering media have attracted great attention in the field of exploration geophysics. The full waveform inversion method is the method with highest parameter modeling precision in the current seismic exploration field, and can fully utilize the kinematics and dynamics information of the seismic wave field. In recent years, it is desirable to perform high-precision parametric modeling on strong scattering media (such as salt domes) encountered in actual exploration by using a full waveform inversion method, so as to improve the imaging quality of the strong scattering body and a shielding region thereof. However, the conventional full waveform inversion method is based on the weak scattering approximation of the born approximation, and the actual seismic data often lacks effective low-frequency information, so that the parameter modeling of the strong scattering medium is difficult to perform directly by using the conventional full waveform inversion method. In recent years, scholars at home and abroad improve the theoretical framework of the full waveform inversion method from different angles so as to adapt to the parameter inversion of a strong scattering medium, and most of the researches are based on the acoustic medium and mainly comprise a Laplace-Fourier domain waveform inversion method, a standard setting method, a full variation constraint method, a deep learning method, a direct envelope inversion method and the like. The multi-parameter modeling research on the elastic strong scattering medium is still in a starting exploration stage, and the existing method has no ideal modeling effect on the inside and the lower speed of the elastic strong scattering body. In summary, there is currently a lack of a method that can model high accuracy longitudinal and transverse wave velocities of strongly scattering media without prior information and in the absence of low frequency data.
The direct envelope inversion method can perform large-scale parameter modeling on the strong scattering medium under the condition that the seismic data lack low frequency and no priori information, and has relatively high calculation efficiency. The biggest difference between the direct envelope inversion method and the full waveform inversion method is that the direct envelope inversion method defines a direct envelope sensitive kernel function (different from a waveform sensitive kernel function of conventional envelope inversion), and low-frequency envelope data disturbance can be directly mapped into large-scale parameter disturbance of a strong scattering medium. However, when the current elastic wave direct envelope inversion method models the longitudinal wave speed and the transverse wave speed of the elastic strong scattering medium, ideal internal speed information of the strong scattering body is difficult to obtain, so that the lower boundary of the strong scattering body and the speed modeling effect below the lower boundary of the strong scattering body are influenced.
Disclosure of Invention
The invention aims to overcome the defects of the prior art, and provides a novel elastic wave multi-parameter inversion method of a strong scattering medium, so that the problem of high-precision modeling of longitudinal and transverse wave speeds of a strong scattering reservoir in oil-gas resource seismic exploration is solved, and a high-precision longitudinal and transverse wave speed model is provided for complex strong scattering reservoir imaging.
The idea of the invention is that: the advantage that the anisotropic total variation constraint can sharpen the model boundary and enhance the uniformity in the layer is fully utilized, and the model boundary and the uniformity in the enhancement layer are introduced into the elastic wave direct envelope inversion process. In each iteration of elastic wave direct envelope inversion, anisotropic total variation constraint is respectively applied to the longitudinal wave velocity model after iteration update, so that the inside of the strong scatterer tends to be uniform and the boundary is more obvious. Finally, the anisotropic total variation constraint elastic wave full waveform inversion technology is combined to obtain a high-precision longitudinal and transverse wave velocity model of the strong scattering medium, and the defects of the prior art are overcome.
The invention aims at realizing the following technical scheme:
firstly, preparing preprocessed elastic wave multi-component observation seismic data; estimating a source wavelet by using observed seismic data, and giving an initial model of longitudinal wave speed and transverse wave speed (no prior information is needed); forward modeling is carried out on the initial velocity model to obtain simulated seismic data and forward seismic wave fields; performing wave field mode decomposition on the forward seismic wave field to obtain longitudinal and transverse wave forward wave fields respectively, obtaining corresponding envelope fields to obtain longitudinal wave and transverse wave forward transmission envelope fields respectively; performing difference on the observed envelope data and the simulated envelope data to obtain an envelope companion source; carrying out accompanying envelope field calculation and mode decomposition on the initial longitudinal wave velocity model and the transverse wave velocity model to obtain a longitudinal wave accompanying envelope field and a transverse wave accompanying envelope field; respectively calculating longitudinal wave and transverse wave velocity gradients by using the forward envelope field and the accompanying envelope field; step length is obtained, and a longitudinal wave speed model and a transverse wave speed model are updated; respectively applying anisotropic total variation constraint to the updated longitudinal wave velocity model and the updated transverse wave velocity model to obtain constrained longitudinal wave velocity model and constrained transverse wave velocity model; performing iterative updating until a stopping condition is met, and obtaining a large-scale longitudinal wave speed structure of the strong scattering medium; and carrying out anisotropic total variation constraint elastic wave full waveform inversion by taking a large-scale longitudinal and transverse wave velocity model of the strong scattering medium as an initial model to obtain a high-precision longitudinal and transverse wave velocity modeling result of the strong scattering medium.
The elastic wave direct envelope inversion method based on anisotropic total variation constraint is realized through a MATLAB platform;
the invention discloses an elastic wave direct envelope inversion method based on anisotropic total variation constraint, which comprises the following steps:
a. installing a MATLAB software platform;
b. carrying out static correction and denoising pretreatment on the data to obtain high-quality elastic wave multi-component observation seismic data;
c. carrying out wavelet estimation on the seismic data, and extracting the source wavelet of each shot of data;
d. obtaining a rough background longitudinal wave speed range and a rough background transverse wave speed range through background speed analysis, generating a background speed model, wherein the background speed model does not contain any priori information of a strong scatterer and is respectively used as an inverted initial longitudinal wave speed model v p0 And an initial shear wave velocity model v s0 ；
e. Calculating elastic wave multicomponent simulated seismic data on an initial longitudinal wave velocity model, and enveloping the simulated data to obtain simulated envelope data
The upper corner mark i represents the component in the i direction, and for the two-dimensional case, refers to the x (horizontal) and z (vertical) directions; taking envelope of the observed seismic data to obtain observed envelope data +.>
Calculating an objective function sigma of elastic wave direct envelope inversion through a formula (1) EDEI ：
Wherein, the summation symbol subscript sr represents integration of all the seismic sources and the detection points, T represents time, and T represents total recording time length;
f. calculating simulated seismic wave field on the initial model, performing wave field mode decomposition on the simulated seismic wave field to obtain forward longitudinal wave field and forward transverse wave field, and respectively taking envelope to obtain longitudinal wave forward envelope field
And transverse wave forward envelope field->
g. Calculating the difference between the simulated envelope data and the observed envelope data to obtain an accompanying source, and reversely transmitting the accompanying source to obtain an accompanying envelope field; decomposing the wave field mode of the adjoint envelope field to obtain a longitudinal wave adjoint envelope field
And transverse wave concomitant envelope field->
h. Zero-delay cross-correlation is carried out on the longitudinal wave forward envelope field and the longitudinal wave accompanying envelope field to obtain the longitudinal wave velocity gradient of the elastic wave direct envelope inversion, and the longitudinal wave velocity gradient is shown in a formula (2):
in the formula, v p Represents longitudinal wave velocity and ρ represents density. The shear wave velocity gradient for elastic wave direct envelope inversion can be calculated from equation (3):
in the formula, v s Represents transverse wave velocity, & represents dot product, μ represents shear modulus;
i. and selecting a proper step length, and updating the longitudinal wave speed model by adopting a steepest descent method. Assuming that the current iteration number is m, using longitudinal wave velocity and transverse wave velocity models obtained by updating the current iteration respectively
And->
It is shown that the process of applying the anisotropic total variation constraint to both is equivalent to solving the optimization problem shown in equations (4) and (5):
wherein J is 1 And J 2 Respectively represent the objective function of applying anisotropic total variation constraint to longitudinal and transverse wave velocity model, alpha 1 And alpha 2 Respectively representThe update step length of the longitudinal wave speed and the transverse wave speed,
and->
Update amounts of longitudinal wave speed and transverse wave speed respectively of mth iteration lambda 1 And lambda (lambda) 2 Weight coefficients for applying anisotropic total variation constraint to longitudinal and transverse wave speeds respectively, the terms represent the two norms ATV Indicating an anisotropic total variation norm. The specific expression for calculating the anisotropic total variation norm for the velocity model v is shown in formula (6):
where nz and nx represent the number of grid points in the vertical and horizontal directions of the model, respectively. Solving the optimization problem shown in formulas (4) and (5) to obtain a longitudinal wave speed model and a transverse wave speed model after the anisotropic total variation constraint as a speed model after the current iteration update;
j. on the updated model, carrying out iteration stop condition judgment; if the stopping condition is not met, taking the updated longitudinal wave speed model as an initial model, and returning to the step e to continue iterative computation; if the stopping condition is met, outputting a result of a strong scattering medium large-scale longitudinal and transverse wave speed structure v after anisotropic total variation constraint pt And v st ；
k. V is set as pt And v st And (3) as an initial model, carrying out full waveform inversion of the anisotropic full-variation constraint elastic wave to obtain a final inversion result, namely a high-precision longitudinal wave speed structure and a high-precision transverse wave speed structure of the strong scattering medium.
Compared with the prior art, the invention has the beneficial effects that: according to the invention, the anisotropic total variation constraint is introduced into the elastic wave direct envelope inversion process, so that the internal speed of the strong scatterer is uniform and the boundary information is outstanding, thereby improving the speed modeling effect of the boundary depiction of the strong scatterer and the shielding region below, and further obtaining the high-precision longitudinal and transverse wave speed inversion result of the elastic strong scattering medium.
Has the following advantages: 1. according to the invention, by introducing the anisotropic total variation constraint, high-quality strong scatterer speed information is easier to obtain in the elastic wave direct envelope inversion process. 2. The anisotropic total variation constraint elastic wave direct envelope inversion can better characterize the boundary information of the strong scatterer, in particular the lower boundary of the strong scatterer and the speed information of a shielding area below the lower boundary of the strong scatterer. 3. The elastic wave direct envelope inversion and the elastic wave full waveform inversion method based on the anisotropic full variation constraint are connected in series, inversion of a large-scale macro structure and a small-scale detail structure of the medium can be considered, and finally, high-precision longitudinal wave speed structures and transverse wave speed structures of the strong scattering medium can be obtained. 4. The invention can carry out high-precision multi-parameter modeling of the elastic strong scattering medium under the condition that the low-frequency information of the seismic data is missing and no model priori information exists.
Drawings
FIG. 1 is a flow chart of an elastic wave direct envelope inversion method based on anisotropic total variation constraints;
FIG. 2 is a graph of a true velocity model;
(a) A true longitudinal wave velocity model map, (b) a true transverse wave velocity model map;
FIG. 3 is a view of a source wavelet and its spectrogram;
(a) A source wavelet map (b) of a source wavelet spectrum;
FIG. 4 is an initial velocity model diagram;
(a) An initial longitudinal wave velocity model map, (b) an initial transverse wave velocity model map;
FIG. 5 is a conventional elastic wave full waveform inversion result;
(a) Inverting a longitudinal wave speed result by using a conventional elastic wave full waveform; (b) inverting the transverse velocity results from the conventional elastic wave full waveform;
FIG. 6 is an unconstrained elastic wave direct envelope inversion result;
(a) Unconstrained elastic wave directly envelopes and inverts longitudinal wave speed results; (b) Unconstrained elastic wave directly envelopes and inverts transverse wave speed results; (c) final inversion of the longitudinal wave velocity; (d) transverse wave velocity final inversion results;
FIG. 7 is an elastic wave direct envelope inversion result based on anisotropic total variation constraints;
(a) Elastic wave direct envelope inversion longitudinal wave speed results based on anisotropic total variation constraint; (b) Elastic wave direct envelope inversion transverse wave speed results based on anisotropic total variation constraint; (c) final inversion of the longitudinal wave velocity; (d) transverse wave velocity final inversion results.
Detailed Description
The invention is described in further detail below with reference to the drawings and examples.
The invention discloses an elastic wave direct envelope inversion based on anisotropic total variation constraint, which comprises the following steps:
a. installing MATLAB software platform under win7 or Linux system requires that MATLAB R2016a and above versions be employed and parallel toolkits (Parallel Computing Toolbox) are provided.
b. Data preprocessing is carried out, static correction processing is carried out on the data, and the influence of the undulating surface on the reflection phase axis is corrected; denoising the data to remove micro-vibration, low-frequency and high-frequency background noise and other random noise; interference waves are removed, including sound waves, surface waves, industrial electric interference, ghost reflection, multiple reflection, side surface waves, bottom waves, reverberation, ringing and the like. And finally observing the seismic data by using high-quality elastic wave multi-component.
c. The wavelet estimation is carried out on the seismic data, and the estimation method can adopt a direct wave estimation method, an autocorrelation method and the like to extract the source wavelet of each shot of data.
d. Obtaining a rough background longitudinal wave speed range and a rough background transverse wave speed range through background speed analysis, generating a background speed model, wherein the background speed model does not contain any priori information of a strong scatterer and is respectively used as an inverted initial longitudinal wave speed model v p0 And an initial shear wave velocity model v s0 。
e. Calculating elastic wave multicomponent simulated seismic data on an initial longitudinal wave velocity model, and enveloping the simulated data to obtain simulated envelope data
The upper corner mark i represents the component in the i direction, and for the two-dimensional case, refers to the x (horizontal) and z (vertical) directions; taking envelope of the observed seismic data to obtain observed envelope data +.>
Calculating an objective function sigma of elastic wave direct envelope inversion through a formula (7) EDEI ：
Where the summation symbol subscript sr denotes integration of all sources and detectors, T denotes time, and T denotes the total recording time length.
f. Calculating simulated seismic wave field on the initial model, performing wave field mode decomposition on the simulated seismic wave field to obtain forward longitudinal wave field and forward transverse wave field, and taking envelope to obtain longitudinal wave forward envelope field
And transverse wave forward envelope field->
Wherein the wavefield pattern decomposition employs the methods of equations (8) and (9):
wherein u represents an elastic wave displacement vector, u p Representing a longitudinal wave field, u s Represents the transverse wave field, ρ represents density, λ and μ represent the lame constants, v represents the gradient operation, v represents the divergence operation, and v represents the rotation operation.
g. Calculating simulated envelope data and observed envelopesData difference to obtain the companion source f s I.e.
The companion source back-transmits to obtain a companion envelope field; decomposing the wave field mode of the adjoint envelope field to obtain a longitudinal wave adjoint envelope field
And transverse wave concomitant envelope field->
The wave field mode decomposition method is shown in formulas (8) and (9).
h. Zero-delay cross-correlation is carried out on the longitudinal wave forward envelope field and the longitudinal wave accompanying envelope field to obtain the longitudinal wave velocity gradient of the elastic wave direct envelope inversion, and the longitudinal wave velocity gradient is shown in a formula (11):
in the formula, v p Represents longitudinal wave velocity and ρ represents density. The shear wave velocity gradient for elastic wave direct envelope inversion can be calculated by equation (12):
in the formula, v s Represents transverse wave velocity, represents dot product, and μ represents shear modulus.
i. And selecting a proper step length, and updating the longitudinal wave speed model by adopting a steepest descent method. Assuming that the current iteration number is m, using longitudinal wave velocity and transverse wave velocity models obtained by updating the current iteration respectively
And->
It is shown that the process of applying the anisotropic total variation constraint to both is equivalent to solving the optimization problem shown in equations (13) and (14):
wherein J is 1 And J 2 Respectively represent the objective function of applying anisotropic total variation constraint to longitudinal and transverse wave velocity model, alpha 1 And alpha 2 Respectively represent the update step length of the longitudinal wave speed and the transverse wave speed,
and->
Update amounts of longitudinal wave speed and transverse wave speed respectively of mth iteration lambda 1 And lambda (lambda) 2 Weight coefficients for applying anisotropic total variation constraint to longitudinal and transverse wave speeds respectively, the terms represent the two norms ATV Indicating an anisotropic total variation norm. The specific expression for calculating the anisotropic total variation norm for the velocity model v is shown in formula (15):
where nz and nx represent the number of grid points in the vertical and horizontal directions of the model, respectively. Solving the optimization problem shown in formulas (13) and (14) to obtain longitudinal wave velocity models and transverse wave velocity models after the anisotropic total variation constraint as velocity models after the current iteration update.
j. On the updated model, carrying out iteration stop condition judgment; if the stopping condition is not met, taking the updated longitudinal wave speed model as an initial model, returning to the step e, and continuing iterationCalculating; if the stopping condition is met, outputting a result of a strong scattering medium large-scale longitudinal and transverse wave speed structure v after anisotropic total variation constraint pt And v st 。
k. V is set as pt And v st And (3) as an initial model, carrying out full waveform inversion of the anisotropic full-variation constraint elastic wave to obtain a final inversion result, namely a high-precision longitudinal wave speed structure and a high-precision transverse wave speed structure of the strong scattering medium.
Example 1
The overall flow of the present invention is shown in figure 1.
It is assumed that the subsurface true longitudinal and transverse wave velocity model is shown in fig. 2a and 2b, respectively. In the real speed model, the background speed is lower, and a high-speed strong scattering salt dome body is arranged in the middle. The forward modeling is performed on the real model to obtain an observation seismic record, the waveform and the frequency spectrum of the source wavelet are shown in figures 3a and 3b respectively, the high-pass filtering processing is performed on the Rake wavelet to simulate the condition of low-frequency information loss in the actual seismic acquisition, the low-frequency information below 3Hz is cut off, and the main frequency of the source wavelet is about 9Hz. And (4) carrying out background velocity analysis by using the observation record to obtain the approximate range of the background longitudinal wave velocity and the transverse wave velocity, and establishing a longitudinal wave and transverse wave initial velocity model as shown in fig. 4a and 4 b. The prior information of the strong scatterer is not contained in the initial longitudinal wave velocity model and the transverse wave velocity model.
In order to compare inversion effects of the method, conventional elastic wave full waveform inversion is firstly carried out on an initial longitudinal wave velocity model and a transverse wave velocity model, and inversion results of the longitudinal wave velocity model and the transverse wave velocity model are shown in figures 5a and 5b respectively. Therefore, due to the lack of low-frequency information, the conventional elastic wave full-wave inversion can only obtain partial information of the top interface of the strong scatterer, and can not recover the speed information inside the strong scatterer and the morphological information of the strong scatterer.
Then, an unconstrained elastic wave direct envelope inversion method is carried out, and longitudinal wave speed inversion results and transverse wave speed inversion results are shown in fig. 6a and 6b respectively. It can be seen that although the strong scatterer morphology is recovered to some extent, there is still a deviation in the internal velocity of the salt dome from the true value. And under the unconstrained condition, taking the results of fig. 6a and 6b as an initial model, performing conventional elastic wave full waveform inversion, and finally performing longitudinal and transverse wave velocity inversion, wherein the results are shown in fig. 6c and 6 d. It can be seen that the internal speed uniformity of the salt dome is poor and the lower boundary of the salt dome is not fully delineated.
Inversion is then performed using the method of the invention. Taking the model shown in fig. 4 as an initial model, according to steps (e) to (j), the longitudinal and transverse wave velocity models which can be obtained are shown in fig. 7a and 7b, and it can be seen that the large-scale longitudinal and transverse wave velocity construction information of the strong scatterer is successfully recovered, and the velocity information in the salt dome is relatively close to a true value. Using fig. 7a and 7b as an initial model, performing step (k), namely performing anisotropic total variation constraint elastic wave full waveform inversion, and finally performing longitudinal and transverse wave velocity inversion, wherein the results are shown in fig. 7c and 7 d. It can be seen that the boundary information and internal velocity of the strong scatterer are better inverted. The inversion effect of the inner part, the lower boundary and the region under the salt dome is better than that under the unconstrained condition. In summary, the final inversion result of the method provided by the invention is that the internal speed and the boundary of the salt dome body are accurate, and the overall effect is obviously better than that of the conventional method (fig. 5 and 6).
Claims (1)
1. An elastic wave direct envelope inversion method based on anisotropic total variation constraint is characterized in that under the condition that low-frequency information is absent in seismic data and model prior information is absent, the anisotropic total variation constraint is applied to an elastic wave direct envelope inversion process to obtain high-quality large-scale longitudinal wave velocity structures and transverse wave velocity structures of elastic strong scattering media; taking the obtained high-quality large-scale longitudinal and transverse wave velocity structure as an initial model of anisotropic total variation constraint elastic wave full waveform inversion to obtain a high-precision longitudinal and transverse wave velocity structure of an elastic strong scattering medium;
the method comprises the following steps:
a. installing a MATLAB software platform;
b. carrying out static correction and denoising pretreatment on the data to obtain high-quality elastic wave multi-component observation seismic data;
c. carrying out wavelet estimation on the seismic data, and extracting the source wavelet of each shot of data;
d. by background velocity analysis, obtainGenerating a background velocity model by using the approximate background longitudinal wave velocity range and the transverse wave velocity range, wherein the background velocity model does not contain any prior information of strong scatterers and is respectively used as an inverted initial longitudinal wave velocity model v p0 And an initial shear wave velocity model v s0 ；
e. Calculating elastic wave multicomponent simulated seismic data on an initial longitudinal wave velocity model, and enveloping the simulated seismic data to obtain simulated envelope data
The upper corner mark i represents the component in the i direction, and for the two-dimensional case, refers to the x horizontal and z vertical directions; taking envelope of the observed seismic data to obtain observed envelope data +.>
Calculating an objective function sigma of elastic wave direct envelope inversion through a formula (1) EDEI ：
Wherein, the summation symbol subscript sr represents integration of all the seismic sources and the detection points, T represents time, and T represents total recording time length;
f. calculating simulated seismic wave field on the initial model, performing wave field mode decomposition on the simulated seismic wave field to obtain forward longitudinal wave field and forward transverse wave field, and respectively taking envelope to obtain longitudinal wave forward envelope field
And transverse wave forward envelope field->
g. Calculating the difference between the simulated envelope data and the observed envelope data to obtain an accompanying source, and reversely transmitting the accompanying source to obtain an accompanying envelope field; decomposing the wave field mode of the adjoint envelope field to obtain a longitudinal wave adjoint envelope field
And transverse wave concomitant envelope field->
h. Zero-delay cross-correlation is carried out on the longitudinal wave forward envelope field and the longitudinal wave accompanying envelope field to obtain the longitudinal wave velocity gradient of the elastic wave direct envelope inversion, and the longitudinal wave velocity gradient is shown in a formula (2):
in the formula, v p The shear wave velocity gradient, which represents the velocity of the longitudinal wave and ρ represents the density, of the elastic wave direct envelope inversion can be calculated by equation (3):
in the formula, v s Represents transverse wave velocity, & represents dot product, μ represents shear modulus;
i. selecting proper step length, adopting the steepest descent method to update the longitudinal wave velocity model and the transverse wave velocity model, and assuming that the current iteration number is m, respectively using the longitudinal wave velocity model and the transverse wave velocity model obtained by updating the current iteration
And->
It is shown that the process of applying the anisotropic total variation constraint to both is equivalent to solving the optimization problem shown in equations (4) and (5): />
Wherein J is 1 And J 2 Respectively represent the objective function of applying anisotropic total variation constraint to longitudinal and transverse wave velocity model, alpha 1 And alpha 2 Respectively represent the update step length of the longitudinal wave speed and the transverse wave speed,
and->
Update amounts of longitudinal wave speed and transverse wave speed respectively of mth iteration lambda 1 And lambda (lambda) 2 Weight coefficients for applying anisotropic total variation constraint to longitudinal and transverse wave speeds respectively, the terms represent the two norms ATV The specific expression for calculating the anisotropic total variation norm for the velocity model v is shown in formula (6):
wherein nz and nx respectively represent the number of grid points in the vertical direction and the horizontal direction of the model, solving the optimization problem shown in formulas (4) and (5), and obtaining longitudinal wave velocity models and transverse wave velocity models after the anisotropic total variation constraint as velocity models after the current iteration update;
j. on the updated model, carrying out iteration stop condition judgment; if the stopping condition is not met, taking the updated longitudinal wave speed model as an initial model, and returning to the step e to continue iterative computation; if the stopping condition is met, outputting a result of a strong scattering medium large-scale longitudinal and transverse wave speed structure v after anisotropic total variation constraint pt And v st ；
k. V is set as pt And v st As an initial model, anisotropic total variation constraint bomb is carried outAnd carrying out full waveform inversion on the characteristic waves to obtain a final inversion result, namely a high-precision longitudinal wave speed structure of the strong scattering medium.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210336492.6A CN114779335B (en) | 2022-03-31 | 2022-03-31 | Elastic wave direct envelope inversion method based on anisotropic total variation constraint |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210336492.6A CN114779335B (en) | 2022-03-31 | 2022-03-31 | Elastic wave direct envelope inversion method based on anisotropic total variation constraint |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114779335A CN114779335A (en) | 2022-07-22 |
| CN114779335B true CN114779335B (en) | 2023-05-02 |
Family
ID=82427147
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210336492.6A Active CN114779335B (en) | 2022-03-31 | 2022-03-31 | Elastic wave direct envelope inversion method based on anisotropic total variation constraint |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114779335B (en) |
Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111766628A (en) * | 2020-07-29 | 2020-10-13 | 浪潮云信息技术股份公司 | Preconditioned time domain elastic medium multi-parameter full waveform inversion method |
Family Cites Families (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7158888B2 (en) * | 2001-05-04 | 2007-01-02 | Takeda San Diego, Inc. | Determining structures by performing comparisons between molecular replacement results for multiple different biomolecules |
| CN102368095B (en) * | 2011-09-10 | 2013-03-27 | 吉林大学 | Extraction method for relaxation time spectrum of nuclear magnetic resonance detection signal for underground water by utilizing multi exponent fitting technology |
| US20160349389A1 (en) * | 2015-05-29 | 2016-12-01 | Cgg Services Sa | Method for developing a geomechanical model based on seismic data, well logs and sem analysis of horizontal and vertical drill cuttings |
| CN108345031B (en) * | 2018-01-11 | 2020-01-17 | 吉林大学 | Full waveform inversion method for elastic medium active source and passive source mixed mining seismic data |
| WO2019234469A1 (en) * | 2018-06-08 | 2019-12-12 | Total Sa | Method for generating an image of a subsurface of an area of interest from seismic data |
| CN110927693B (en) * | 2019-12-23 | 2021-07-27 | 航天南湖电子信息技术股份有限公司 | Pulse compression method combining matched filtering with sparse inversion |
| CN111273348B (en) * | 2020-01-21 | 2021-02-05 | 长江大学 | Multipoint geostatistical prestack inversion method based on updated probability ratio constant theory |
| CN111505714B (en) * | 2020-04-16 | 2021-05-25 | 吉林大学 | Elastic wave direct envelope inversion method based on rock physical constraint |
-
2022
- 2022-03-31 CN CN202210336492.6A patent/CN114779335B/en active Active
Patent Citations (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111766628A (en) * | 2020-07-29 | 2020-10-13 | 浪潮云信息技术股份公司 | Preconditioned time domain elastic medium multi-parameter full waveform inversion method |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114779335A (en) | 2022-07-22 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108873066B (en) | Elastic medium wave equation reflected wave travel time inversion method | |
| US11048001B2 (en) | Methods using travel-time full waveform inversion for imaging subsurface formations with salt bodies | |
| CN106526677B (en) | A kind of wideband reverse-time migration imaging method of marine adaptive compacting ghost reflection | |
| US11294088B2 (en) | Methods for simultaneous source separation | |
| US8081540B2 (en) | Rotary subwoofer marine seismic source | |
| CN108549100B (en) | The multiple dimensioned full waveform inversion method of time-domain for opening up frequency based on non-linear high order | |
| RU2007101245A (en) | ATTENUATION OF WAVES-SATELLITES IN THREE-DIMENSIONAL SPACE | |
| CN110687600B (en) | Elastic wave least square reverse time migration method based on acoustic-elastic coupling equation | |
| CN107340540B (en) | Direction wave decomposition method, device and the computer storage medium of elastic wave field | |
| CN109946742B (en) | Pure qP wave seismic data simulation method in TTI medium | |
| CN111505718A (en) | High-resolution underground structure amplitude-preserving imaging method | |
| CN110542928A (en) | Seismic response simulation method based on VTI anisotropic propagation matrix | |
| CN110888159B (en) | Elastic wave full waveform inversion method based on angle decomposition and wave field separation | |
| CN111505714B (en) | Elastic wave direct envelope inversion method based on rock physical constraint | |
| GB2409901A (en) | Determining shear wave velocity from tube wave characteristics | |
| CN111077567A (en) | Method for double-pass wave prestack depth migration based on matrix multiplication | |
| US11199641B2 (en) | Seismic modeling | |
| CN114779335B (en) | Elastic wave direct envelope inversion method based on anisotropic total variation constraint | |
| LI et al. | Multi-scale full waveform inversion based ongradient decomposition in wavenumber domain | |
| Takekawa et al. | An accuracy analysis of a Hamiltonian particle method with the staggered particles for seismic-wave modeling including surface topography | |
| CN109975873B (en) | Method and system for removing low-frequency noise by reverse time migration imaging | |
| CN111175822B (en) | Strong scattering medium inversion method for improving direct envelope inversion and disturbance decomposition | |
| Zhong et al. | Elastic reverse time migration method in vertical transversely isotropic media including surface topography | |
| CN105607123A (en) | Method and apparatus for calculating seismic wave characteristic information of random porous media model | |
| CN114488302A (en) | In-situ anisotropic ground stress field prediction method and system |
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 | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |