CN110703318B - Direct inversion method of pre-stack seismic data - Google Patents
Direct inversion method of pre-stack seismic data Download PDFInfo
- Publication number
- CN110703318B CN110703318B CN201910904075.5A CN201910904075A CN110703318B CN 110703318 B CN110703318 B CN 110703318B CN 201910904075 A CN201910904075 A CN 201910904075A CN 110703318 B CN110703318 B CN 110703318B
- Authority
- CN
- China
- Prior art keywords
- data
- mapping
- inversion
- observation
- observation data
- 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 41
- 238000013507 mapping Methods 0.000 claims abstract description 34
- 238000004364 calculation method Methods 0.000 claims abstract description 17
- 238000004140 cleaning Methods 0.000 claims description 3
- 238000012217 deletion Methods 0.000 claims description 3
- 230000037430 deletion Effects 0.000 claims description 3
- 230000001629 suppression Effects 0.000 claims description 3
- 230000008030 elimination Effects 0.000 claims description 2
- 238000003379 elimination reaction Methods 0.000 claims description 2
- 238000004088 simulation Methods 0.000 description 5
- 238000010586 diagram Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000005457 optimization Methods 0.000 description 2
- 230000002159 abnormal effect Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000007812 deficiency Effects 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 230000000704 physical effect Effects 0.000 description 1
- 230000003252 repetitive effect Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 239000013049 sediment Substances 0.000 description 1
- 238000003786 synthesis reaction Methods 0.000 description 1
- 230000017105 transposition 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/282—Application of seismic models, synthetic seismograms
-
- 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
- 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/306—Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles
-
- 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/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
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
According to the method for directly inverting the pre-stack seismic data, the mapping from the observation data to the model parameters is constructed to be used as the inversion operator, and the constructed inversion operator is used for directly inverting the pre-stack seismic data, so that the calculation efficiency is improved, the calculation amount is reduced, and meanwhile, the problems of dependence on the initial model and reduction of the convergence of the inversion result are effectively avoided. The method comprises the following implementation steps: (1) acquiring a limited group of observation data and model parameter data which is less than m and corresponds to the limited group of observation data, wherein D is used as sample data; (2) sorting the obtained sample data; (3) constructing a mapping from observed data to model parameters(4) And directly inverting the data to be inverted by utilizing the constructed mapping relation.
Description
Technical Field
The invention relates to a method for inverting prestack seismic data, in particular to a method for directly inverting after constructing a mapping relation from the prestack seismic data to model parameters by using sample data, and belongs to the field of seismic exploration.
Background
The seismic data prestack inversion is an effective means for solving the physical property parameters of the underground medium based on the change relation of the reflection coefficient along with the incident angle. The prestack inversion technology takes the Zoeppritz (1919) equation and the approximate equation thereof (Bortfeld, 1961; Aki and Richards, 1980; Shuey,1985) as the theoretical basis, has been proposed since the 80 th century (Ostrander,1982,1984), and is continuously developed and plays an important role in the field of oil and gas resource exploration.
In the prestack seismic data, the variation of the seismic reflection coefficient with the incidence angle includes various physical parameters of the stratum, such as: longitudinal wave velocity, shear wave velocity, density, and the like, and combinations thereof. In the pre-stack seismic inversion, observation data are pre-stack gather data which come from the same reflection point and contain reflection coefficient change information along with incident angles, parameters to be inverted are model parameters containing longitudinal wave velocity, transverse wave velocity, density and the like, and the relation is established between the observation data D and a model parameter vector m through a forward operator f:
D=f(m) (1)
in the above equation (1), the positive operator f is a mapping from the model parameter vector m to the observation data D, and is specifically Zoeppritz (1919) equation in the prestack inversion.
However, the forward operator f is only a mapping from m to D, and the geophysical inversion is a process of solving the model parameter m from the known observation data D, and the process cannot be directly realized by the equation (1). Moreover, the forward modeling process shown in equation (1) is often a multidimensional, highly nonlinear mapping relationship, and a mapping from the observation data D to the model parameter m cannot be directly obtained by inverting f:
m=f-1(D) (2)
therefore, in the inversion problem, an objective function is constructed through the formula (1), and the optimal model parameters are solved as the inversion result by using a forward iteration optimization method. The procedure is shown in the attached FIG. 1, first, initial model parameters m are given, and simulation data D is obtained by forward calculation (formula (1))synRecording actual observation data as DobsThen the least squares objective function can be written as:
determining an objective functionWhether or not the minimum value is reached, ifIf the minimum value is not reached, further adjusting the model parameter m, and obtaining new simulation data D through forward calculation of the adjusted msyn(ii) a Judgment ofWhether the minimum value is reached or not, and if the minimum value is reached, outputting m at the moment as an inversion result m'; if the minimum value is not reached, continuing to adjust m, and continuing to iteratively calculate until the minimum value is reachedReaching a minimum value.
The complexity of the forward operator determines the difficulty of inversion, and the optimization solving process based on forward iteration is often large in calculated amount and low in calculation efficiency, and has the problems of large dependence on an initial model and convergence. The exact Zoeppritz equation is also not physically intuitive due to its mathematically complex form.
To overcome the mathematical complexity and physical non-intuitiveness of the Zoeppritz equation, many scholars approximate the equation from different perspectives. The approximation equation, while overcoming the deficiencies of the exact equation, reduces the accuracy of the calculation.
As summarized above, although the precision of the iterative inversion method based on the precise equation is high, the iterative inversion method is large in calculation amount, low in calculation efficiency, large in dependence on an initial model, and has the problem of convergence. While the linear inversion technology based on the approximate equation has high calculation efficiency, the inversion result has low precision.
In view of this, the present patent application is specifically proposed.
Disclosure of Invention
The direct inversion method of the prestack seismic data is designed to solve the problems in the prior art, and the direct inversion method is provided.
The direct inversion method of the pre-stack seismic data mainly comprises the following implementation steps:
(1) acquiring a limited group of observation data and model parameter data which is less than m and corresponds to the limited group of observation data, wherein D is used as sample data;
(2) sorting the obtained sample data;
(4) And directly inverting the data to be inverted by utilizing the constructed mapping relation.
Further, in the step (1), each set of the primitive data < m, D > needs to satisfy the mapping relationship from the model parameter m to the observation data D within a certain error range, and the sample data can be obtained by theoretical forward calculation or can be obtained from actual observation data. The observation data D is prestack gather data containing the variation of the seismic reflection coefficient with the incidence angle.
In the step (2), the sample data is sorted by adopting data cleaning operations including outlier rejection, noise suppression and repeated sample deletion.
The step (3) specifically includes:
(3-1) regarding the obtained multiple groups of sample data as discrete samples on the mapping F from the observation data to the model parameters to be constructed, and constructing the mapping F from the observation data to the model parameters in the sample data coverage range by utilizing an interpolation or regression method.
wherein,for any pre-stack seismic data to be inverted in the coverage range of a plurality of groups of obtained sample data < m and D,for the inversion result corresponding to the pre-stack data, wiRepresenting their corresponding sub-functions for the weight coefficientsThe proportion in the mapping.
(3-3) weight coefficient w in step (3-2)iSum subfunction giSamples obtained as described aboveThis data is solved:
in the formula, | m-F (D) | represents that the sample data is less than m, D > is the distance between the model parameter m and the observation data D after inversion by mapping F, sigma-represents that the obtained finite group of samples are summed,represents w when the objective function takes the minimum valueiAnd gi。
In summary, compared with the iterative inversion method in the prior art, the direct inversion method has the following advantages:
1. the constructed mapping from the observation data to the model parameters is directly utilized to carry out direct inversion, and an iterative process adopted by solving a target function in the traditional inversion process is not needed, so that the problems of initial values, convergence and the like can be effectively avoided or reduced.
2. The constructed inversion operator is stored for future use, not only for one specific inversion task, thus greatly improving the computational efficiency.
3. On the basis of improving the calculation efficiency, the accuracy of the inversion result cannot be reduced.
Drawings
FIG. 1 is a schematic flow chart of an iterative inversion method adopted in the prior art;
FIG. 2 is a schematic flow diagram of a direct inversion method described herein;
FIG. 3 is a schematic view of the transverse distribution of elastic parameters;
FIG. 4 is a graph of reflection coefficient as a function of angle of incidence;
FIG. 5 is a schematic diagram of the results of a direct inversion;
FIG. 6 is a graph of error comparison of direct inversion and conventional inversion method results;
Detailed Description
Example 1 the present invention will now be described in further detail with reference to the accompanying drawings, which are intended to represent only some, but not all, embodiments of the present application. All other embodiments, which can be obtained by persons skilled in the art without any inventive work, based on the embodiments of the present invention, belong to the scope of the present invention.
As shown in fig. 2, the direct inversion method of pre-stack seismic data mainly includes the following implementation steps:
(1) acquiring a limited group of observation data and model parameter data which is less than m and corresponds to the limited group of observation data, wherein D is used as sample data;
the data of each group is less than m, D is more than the requirement of meeting the mapping relation between the model parameter m and the observation data D within a certain error range, and the sample data can be directly obtained from actual observation data or can be calculated by theoretical forward modeling. In this example, the sample data is obtained by theoretical forward modeling calculation, that is, observation data is obtained by model parameter calculation using the Zoeppritz equation shown in formula (1). Model parameter ml=[αiβjρk]In this example, the longitudinal wave velocity α of the sea bottomiSea bottom transverse wave velocity betajAnd seafloor density ρkRespectively has a value range of 1530-alphai≤1570、180≤βj≤220、1630≤ρkThe value interval is less than or equal to 1670 and is 2.5. From mlForward calculation to obtain observation data DlWherein the observation data DlThe sea bottom reflection coefficient in the range of 1 deg. to 55 deg. is calculated at intervals of 1 deg..
(2) Sorting the obtained sample data;
the sample data is sorted by adopting data cleaning operations of abnormal value elimination, noise suppression, repeated sample deletion and the like. As can be seen from the above forward calculations, i is 1,2, …,17, j is 1,2, …,17, k is 1,2, …, 17; then the model parameter mlThere are a total of 17 × 17 × 17 — 4913 combinations. We have obtained 4913 sets of non-repetitive sample data altogether<ml,Dl>,l=1,2,…,4913。
(3-1) the obtained multiple groups of sample data < ml,DlRegarding as discrete sampling on the mapping F from the observation data to the model parameters to be constructed, and constructing the mapping F from the observation data to the model parameters in the sample data coverage range by utilizing an interpolation or regression method.
(3-3) the weight coefficient w in the step (3-2) is obtained by using the formula (5)iSum subfunction giSubstituting all sample data into < ml,DlMore than get:
M=G·w (7)
To this end, we have obtained the mapping to be constructed(i.e., equation (4)) that the construction of the mapping F from the observation data D to the model parameters m is completed.
(4) And directly inverting the data to be inverted by utilizing the constructed mapping relation.
In order to facilitate the analysis of the accuracy of the inversion result, the data to be inverted used in the embodiment is theoretical forward data, and the longitudinal wave velocity alpha of the sediment at the bottom of the sea is assumed2Transverse wave velocity beta2And density ρ2Distribution along a two-dimensional seismic line, as shown in FIG. 3:
specifically, the elastic parameter distribution is given by equation (8):
in the above formula, n is 1,2, …,200 is CRP gather number.
The sea bottom reflection coefficient in the range of 1 ° to 55 ° was calculated at intervals of 1 ° for each reflection point (each CRP), and random errors were added to the obtained reflection coefficients to simulate data containing noise in practice, as shown in fig. 4, which shows the reflection coefficients at CRP numbers 25, 75, and 150, respectively, as a function of the incident angle.
Thus, a two-dimensional line theoretical simulation data comprising 200 CRP gathers was obtained.
The theoretical simulation data used here only replaces the actual observation data, which is convenient for comparing inversion errors; when the actual data is used for inversion, the simulation synthesis of the seismic data is not needed to be carried out firstly, but the actual data is directly obtained and inverted.
Taking the observation data as the mapping from the observation data to the model parameters constructed in step (3)The obtained model parameter values are the direct inversion results, and the inversion results are shown in fig. 5.
In fig. 5, a thick solid line represents model parameters, and a thin solid line represents inversion results when observed data contains noise; the left side shows the results obtained with the direct inversion method proposed by the present invention, and the right side is the results obtained with the conventional inversion. Obviously, under the influence of noise, the inversion result jumps up and down with a given model value as a center, and the error of the transverse wave velocity is the largest among three parameters of the inversion, and the longitudinal wave velocity and the density are the second order, which is consistent with the traditional inversion result.
As shown in fig. 6, a comparison between the inversion error of the direct inversion method proposed in the present application and the inversion error of the conventional inversion method is shown. In the figure, "+" represents the inversion error of the direct inversion method provided by the invention, and a solid line represents the inversion error obtained by the traditional inversion method. The two errors are obviously basically consistent, and the direct inversion method provided by the invention can achieve the inversion accuracy consistent with that of the traditional inversion method.
In summary, the method provided by the present invention has the advantages described above, and the accuracy of the inversion result is not reduced, so that the method provided by the present invention is feasible.
Claims (4)
1. A direct inversion method of pre-stack seismic data is characterized by comprising the following steps: comprises the following implementation steps of the method,
(1) acquiring a limited group of observation data and model parameter data which is less than m and corresponds to the limited group of observation data, wherein D is used as sample data;
(2) sorting the obtained sample data;
(4) Directly inverting the data to be inverted by utilizing the constructed mapping relation;
wherein,for any pre-stack seismic data to be inverted in the coverage range of a plurality of groups of obtained sample data < m and D,for the inversion result corresponding to the pre-stack data, wiRepresenting their corresponding sub-functions for the weight coefficientsThe proportion in the mapping;
weight coefficient wiSum subfunction giThe calculation of (1):
in the formula, | m-F (D) | represents that the sample data is less than m, D > is the distance between the model parameter m and the observation data D after inversion by mapping F, sigma-represents that the obtained finite group of samples are summed,represents w when the objective function takes the minimum valueiAnd gi。
2. The method of direct inversion of prestack seismic data according to claim 1, characterized by: in the step (1), each group of sample data is less than m, D > needs to satisfy the mapping relationship from the model parameter m to the observation data D within a certain error range, the sample data can be obtained by theoretical forward calculation or obtained from actual observation data, and the observation data D is pre-stack gather data containing the variation relationship of the seismic reflection coefficient with the incident angle.
3. The method of direct inversion of prestack seismic data according to claim 1, characterized by: in the step (2), the sample data is sorted by adopting data cleaning operations including outlier elimination, noise suppression and repeated sample deletion.
4. The method of direct inversion of prestack seismic data according to claim 1, characterized by: in the step (3), the obtained multiple groups of sample data are regarded as discrete samples on the mapping F from the observation data to the model parameters to be constructed, and the mapping F from the observation data to the model parameters in the sample data coverage range is constructed by utilizing an interpolation or regression method.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910904075.5A CN110703318B (en) | 2019-09-24 | 2019-09-24 | Direct inversion method of pre-stack seismic data |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910904075.5A CN110703318B (en) | 2019-09-24 | 2019-09-24 | Direct inversion method of pre-stack seismic data |
Publications (2)
Publication Number | Publication Date |
---|---|
CN110703318A CN110703318A (en) | 2020-01-17 |
CN110703318B true CN110703318B (en) | 2021-06-11 |
Family
ID=69195894
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910904075.5A Active CN110703318B (en) | 2019-09-24 | 2019-09-24 | Direct inversion method of pre-stack seismic data |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110703318B (en) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN113569493B (en) * | 2021-09-26 | 2021-12-17 | 自然资源部第一海洋研究所 | Domain decomposition-based physical driving deep learning inversion method |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103097914A (en) * | 2010-04-06 | 2013-05-08 | 道达尔公司 | A process for characterising the evolution of a reservoir |
CN107092029A (en) * | 2017-04-26 | 2017-08-25 | 中国石油大学(北京) | A kind of seismic inversion method and device |
CN108072892A (en) * | 2016-11-09 | 2018-05-25 | 中国石油化工股份有限公司 | A kind of geological structure constraint chromatography conversion method of automation |
CN108139499A (en) * | 2015-10-02 | 2018-06-08 | 埃克森美孚上游研究公司 | The full wave field inversion of Q- compensation |
Family Cites Families (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US8923093B2 (en) * | 2009-08-25 | 2014-12-30 | Westerngeco L.L.C. | Determining the quality of a seismic inversion |
US10670751B2 (en) * | 2015-03-27 | 2020-06-02 | Cgg Services Sas | Full waveform inversion method for seismic data processing using preserved amplitude reverse time migration |
-
2019
- 2019-09-24 CN CN201910904075.5A patent/CN110703318B/en active Active
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103097914A (en) * | 2010-04-06 | 2013-05-08 | 道达尔公司 | A process for characterising the evolution of a reservoir |
CN108139499A (en) * | 2015-10-02 | 2018-06-08 | 埃克森美孚上游研究公司 | The full wave field inversion of Q- compensation |
CN108072892A (en) * | 2016-11-09 | 2018-05-25 | 中国石油化工股份有限公司 | A kind of geological structure constraint chromatography conversion method of automation |
CN107092029A (en) * | 2017-04-26 | 2017-08-25 | 中国石油大学(北京) | A kind of seismic inversion method and device |
Non-Patent Citations (1)
Title |
---|
非高斯AVO三参数反演算法及其应用研究;刘洋;《中国博士学位论文全文数据库 基础科学辑》;20151115(第11期);第24页 * |
Also Published As
Publication number | Publication date |
---|---|
CN110703318A (en) | 2020-01-17 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN106526674B (en) | Three-dimensional full waveform inversion energy weighting gradient preprocessing method | |
US8406081B2 (en) | Seismic imaging systems and methods employing tomographic migration-velocity analysis using common angle image gathers | |
CN110780351B (en) | Longitudinal wave and converted wave prestack joint inversion method and system | |
AU2012220584A1 (en) | Sensitivity kernel-based migration velocity analysis in 3D anisotropic media | |
CN111025387B (en) | Pre-stack earthquake multi-parameter inversion method for shale reservoir | |
CN109521469B (en) | Regularization inversion method for elastic parameters of submarine sediments | |
US9952341B2 (en) | Systems and methods for aligning a monitor seismic survey with a baseline seismic survey | |
CN112698390B (en) | Pre-stack seismic inversion method and device | |
CN105629299A (en) | Travel-time table and angle table acquisition method for angle domain prestack depth migration and imaging method | |
IL92132A (en) | Homeomorphical imaging method of analyzing the structure of a medium | |
WO2021127382A1 (en) | Full waveform inversion in the midpoint-offset domain | |
CN111025388B (en) | Multi-wave combined prestack waveform inversion method | |
CN110703318B (en) | Direct inversion method of pre-stack seismic data | |
CN111175821B (en) | Step-by-step inversion method for anisotropic parameters of VTI medium | |
CN107976714B (en) | A kind of channel set calculation method of complicated earth surface classification spatial distance weighting | |
CN111538084B (en) | Method and system for converting OVT domain data into azimuth angle domain imaging gather | |
CN117388944A (en) | Multi-physical parameter inversion method of geologic model constraint | |
TANG et al. | Multiple Ray Tracing Within 3‐D Layered Media with the Shortest Path Method | |
CN110208858B (en) | Method and system for directly estimating 'sweet spot' probability based on pre-stack inversion | |
CN112083492A (en) | Full-path compensation primary wave and multiple wave combined imaging method under deep sea environment | |
CN112526605B (en) | Method for simulating and exploring natural gas hydrate by adopting seismic numerical value | |
CN115184986B (en) | Global envelope cross-correlation full waveform inversion method independent of seismic source | |
NL2031350B1 (en) | Method and system for determining geosteering irregular observation system | |
CN118244355B (en) | Reflection wave travel time inversion method based on reconstructed observation seismic data | |
Purba et al. | Improving the accuracy of the expanded anisotropic eikonal equation at larger offsets using Levin T-transformation |
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 |