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 PDF

Info

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
Application number
CN202210336492.6A
Other languages
Chinese (zh)
Other versions
CN114779335A (en
Inventor
张盼
韩立国
周奕秀
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Jilin University
Original Assignee
Jilin University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Jilin University filed Critical Jilin University
Priority to CN202210336492.6A priority Critical patent/CN114779335B/en
Publication of CN114779335A publication Critical patent/CN114779335A/en
Application granted granted Critical
Publication of CN114779335B publication Critical patent/CN114779335B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30Analysis
    • G01V1/303Analysis for determining velocity profiles or travel times
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/60Analysis
    • G01V2210/62Physical property of subsurface
    • G01V2210/622Velocity, density or impedance
    • G01V2210/6222Velocity; 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

Elastic wave direct envelope inversion method based on anisotropic total variation constraint
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
Figure BDA0003574503380000041
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 +.>
Figure BDA0003574503380000042
Calculating an objective function sigma of elastic wave direct envelope inversion through a formula (1) EDEI
Figure BDA0003574503380000043
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
Figure BDA0003574503380000044
And transverse wave forward envelope field->
Figure BDA0003574503380000045
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
Figure BDA0003574503380000046
And transverse wave concomitant envelope field->
Figure BDA0003574503380000047
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):
Figure BDA0003574503380000048
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):
Figure BDA0003574503380000049
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
Figure BDA00035745033800000410
And->
Figure BDA0003574503380000051
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):
Figure BDA0003574503380000052
Figure BDA0003574503380000053
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,
Figure BDA0003574503380000054
and->
Figure BDA0003574503380000055
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):
Figure BDA0003574503380000056
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
Figure BDA0003574503380000081
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 +.>
Figure BDA0003574503380000082
Calculating an objective function sigma of elastic wave direct envelope inversion through a formula (7) EDEI
Figure BDA0003574503380000083
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
Figure BDA0003574503380000084
And transverse wave forward envelope field->
Figure BDA0003574503380000085
Wherein the wavefield pattern decomposition employs the methods of equations (8) and (9):
Figure BDA0003574503380000086
Figure BDA0003574503380000087
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.
Figure BDA0003574503380000088
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
Figure BDA0003574503380000089
And transverse wave concomitant envelope field->
Figure BDA00035745033800000810
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):
Figure BDA0003574503380000091
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):
Figure BDA0003574503380000092
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
Figure BDA0003574503380000093
And->
Figure BDA0003574503380000094
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):
Figure BDA0003574503380000095
Figure BDA0003574503380000096
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,
Figure BDA0003574503380000097
and->
Figure BDA0003574503380000098
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):
Figure BDA0003574503380000099
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
Figure FDA0004107225950000011
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 +.>
Figure FDA0004107225950000012
Calculating an objective function sigma of elastic wave direct envelope inversion through a formula (1) EDEI
Figure FDA0004107225950000013
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
Figure FDA0004107225950000014
And transverse wave forward envelope field->
Figure FDA0004107225950000015
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
Figure FDA0004107225950000021
And transverse wave concomitant envelope field->
Figure FDA0004107225950000022
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):
Figure FDA0004107225950000023
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):
Figure FDA0004107225950000024
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
Figure FDA0004107225950000025
And->
Figure FDA0004107225950000026
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): />
Figure FDA0004107225950000027
Figure FDA0004107225950000028
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,
Figure FDA0004107225950000029
and->
Figure FDA00041072259500000210
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):
Figure FDA00041072259500000211
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.
CN202210336492.6A 2022-03-31 2022-03-31 Elastic wave direct envelope inversion method based on anisotropic total variation constraint Active CN114779335B (en)

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)

* Cited by examiner, † Cited by third party
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)

* Cited by examiner, † Cited by third party
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

Patent Citations (1)

* Cited by examiner, † Cited by third party
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
CN111505714B (en) Elastic wave direct envelope inversion method based on rock physical constraint
CN110058303A (en) Acoustic anisotropy reverse-time migration mixed method
CN110687600B (en) Elastic wave least square reverse time migration method based on acoustic-elastic coupling equation
CN111505718B (en) High-resolution underground structure amplitude-preserving imaging method
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
CN110542928A (en) Seismic response simulation method based on VTI anisotropic propagation matrix
CN111077567B (en) Method for double-pass wave prestack depth migration based on matrix multiplication
CN110888159B (en) Elastic wave full waveform inversion method based on angle decomposition and wave field separation
CN114488302A (en) In-situ anisotropic ground stress field prediction method and system
GB2409901A (en) Determining shear wave velocity from tube wave characteristics
US11199641B2 (en) Seismic modeling
CN114779335B (en) Elastic wave direct envelope inversion method based on anisotropic total variation constraint
Takekawa et al. An accuracy analysis of a Hamiltonian particle method with the staggered particles for seismic-wave modeling including surface topography
CN109885906B (en) Magnetic resonance sounding signal sparse noise elimination method based on particle swarm optimization
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
CN104732093A (en) FCT-FDM forward simulation method based on dispersion viscosity wave equation
CN114488286B (en) Amplitude weighting-based towline and seabed seismic data joint waveform inversion method
CN115542387A (en) Method and system for pre-stack migration to produce flare angle migration gathers

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