CN107193043B - A kind of subsurface structure imaging method of relief surface - Google Patents
A kind of subsurface structure imaging method of relief surface Download PDFInfo
- Publication number
- CN107193043B CN107193043B CN201710337714.5A CN201710337714A CN107193043B CN 107193043 B CN107193043 B CN 107193043B CN 201710337714 A CN201710337714 A CN 201710337714A CN 107193043 B CN107193043 B CN 107193043B
- Authority
- CN
- China
- Prior art keywords
- data
- gradient
- model
- iteration
- offset 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.)
- Expired - Fee Related
Links
- 238000003384 imaging method Methods 0.000 title claims abstract description 46
- 238000002310 reflectometry Methods 0.000 claims abstract description 39
- 238000013508 migration Methods 0.000 claims abstract description 37
- 230000005012 migration Effects 0.000 claims abstract description 37
- 238000005070 sampling Methods 0.000 claims abstract description 37
- 238000000034 method Methods 0.000 claims abstract description 22
- 238000012952 Resampling Methods 0.000 claims abstract description 14
- 238000012986 modification Methods 0.000 claims description 19
- 230000004048 modification Effects 0.000 claims description 19
- 238000004364 calculation method Methods 0.000 claims description 11
- 230000021615 conjugation Effects 0.000 claims description 5
- 230000017105 transposition Effects 0.000 claims description 3
- 230000006378 damage Effects 0.000 abstract description 4
- 238000011478 gradient descent method Methods 0.000 description 5
- 238000002939 conjugate gradient method Methods 0.000 description 4
- 238000010276 construction Methods 0.000 description 4
- 238000004088 simulation Methods 0.000 description 4
- 239000011159 matrix material Substances 0.000 description 3
- 238000010586 diagram Methods 0.000 description 2
- 238000007689 inspection Methods 0.000 description 2
- 238000010606 normalization Methods 0.000 description 2
- 238000011161 development Methods 0.000 description 1
- 235000013399 edible fruits Nutrition 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000036039 immunity Effects 0.000 description 1
- 238000012804 iterative process Methods 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 230000000630 rising effect Effects 0.000 description 1
- 150000003839 salts Chemical class 0.000 description 1
- 239000004576 sand Substances 0.000 description 1
- 238000002945 steepest descent method Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000003325 tomography Methods 0.000 description 1
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/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/34—Displaying seismic recordings or visualisation of seismic data or attributes
- G01V1/345—Visualisation of seismic data or attributes, e.g. in 3D cubes
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V20/00—Geomodelling in general
-
- 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/624—Reservoir parameters
-
- 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/64—Geostructures, e.g. in 3D data cubes
-
- 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/66—Subsurface modeling
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/70—Other details related to processing
- G01V2210/74—Visualisation of seismic data
Landscapes
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- Investigating Or Analysing Materials By Optical Means (AREA)
Abstract
The invention discloses a kind of subsurface structure imaging methods of relief surface, comprising the following steps: S1, input observation data, migration velocity field, while presetting the initial value of reflectivity model；S2 carries out sampling or resampling based on related coefficient to observation data, obtains offset data；S3 obtains data residual error according to current reflectivity model and offset data；S4 determines gradient according to data residual error, carries out reflectivity model iteration based on gradient, obtains final reflectivity model.The method of the present invention is by introducing the sampling policy based on related coefficient, avoid stochastical sampling, it can be effectively prevented from offset data while reducing offset data amount and generate great variety with the number of iterations, and then weaken replacement offset data to the destruction of gradient direction conjugacy and the jitter phenomenon of convergence curve, improve computational efficiency and final image quality when iteration.
Description
Technical field
The present invention relates to a kind of subsurface structure imaging methods of relief surface, belong to oil gas physical prospecting engineering field.
Background technique
With going deep into for China's exploration and development, explores key area and be gradually transferred to western exploratory area.The exploration in western exploratory area
There are following main features for exploitation: (1) surface relief is violent, and near surface construction is complicated；(2) Background Construction is complicated, such as inverse covers
It pushes away and covers and imaging etc. under salt；(3) reservoir buries relatively deep, and signal is weaker.Therefore, how real the Important Problems of western exploratory area exploration are
Efficient, high-precision under the conditions of existing relief surface protect width imaging.
Current conventional inversion imaging method is solved by iterative algorithms such as steepest descent method, conjugate gradient methods, although can
To obtain the imaging results of high quality, but calculation amount is excessively huge.If carrying out stochastical sampling to observation data, updating every time
Calculating section data, can reduce calculation amount when model.However, error function when stochastic sampling strategy leads to inversion imaging
All there is randomness with gradient.It is adjacent when subsurface structure is complicated or surface relief is violent due to observation data stochastical sampling
Correlation is poor between offset data in iteration twice, and acutely shaking often occurs in the convergence curve of stochastic gradient descent method, and is carrying out
Image quality is poor when underground structure is imaged, and convergence rate is slower, and computational efficiency is low.
Summary of the invention
The object of the present invention is to provide a kind of subsurface structure imaging methods of relief surface, it can solve current skill
The problem of art, improves computational efficiency when being imaged under the conditions of relief surface, while not influencing its image quality.
In order to solve the above technical problems, the present invention adopts the following technical scheme that: a kind of subsurface structure of relief surface at
Image space method, comprising the following steps:
S1, input observation data, migration velocity field, while presetting the initial value of reflectivity model；
S2 carries out sampling or resampling based on related coefficient to observation data, obtains offset data；
S3 obtains data residual error according to current reflectivity model and offset data；
S4 determines stochastic gradient and/or random conjugate gradient according to data residual error, based on stochastic gradient and/or total at random
Yoke gradient carries out reflectivity model iteration, obtains final reflectivity model.
Compared with prior art, when being imaged under the conditions of relief surface with method of the invention, reducing offset data amount
While can effectively reduce offset data with the variation of the number of iterations, and then weaken replacement offset data and gradient direction be conjugated
The destruction of property and the jitter phenomenon of convergence curve, improve the convergence stability and its image quality when iteration.In addition, of the invention
Relief surface reverse-time migration will be introduced in Least squares inversion thought, can suppress imaging noise, balanced imaging amplitude, improve at
As resolution ratio, near surface can be constructed and infrastructure high-precision is imaged.
Detailed description of the invention
Fig. 1 is the flow diagram of one embodiment of the present of invention；
Fig. 2 is the migration velocity field of one embodiment of the present of invention；
Fig. 3 is the relief surface elevation information of one embodiment of the present of invention；
Fig. 4 is the observation seismic data of one embodiment of the present of invention；
Fig. 5 is the carrying out observation data based on obtaining after related coefficient big gun data sampling of one embodiment of the present of invention
Offset data；
Fig. 6 is the relief surface reverse-time migration result of one embodiment of the present of invention；
Fig. 7 is the relief surface least square reverse-time migration result based on stochastic gradient descent method；
Fig. 8 is the relief surface least square reverse-time migration result of one embodiment of the present of invention；
Fig. 9 is one embodiment of the invention and the single track amplitude comparison of other two kinds of imaging method results；
Figure 10 is the normalization data residual error convergence curve of one embodiment of the invention and other two kinds of imaging method results；
Figure 11 is one embodiment of the invention and the calculating time comparison of other two kinds of imaging methods；
Figure 12 is one embodiment of the present of invention flow diagram.
The present invention is further illustrated with reference to the accompanying drawings and detailed description.
Specific embodiment
The embodiment of the present invention 1: a kind of subsurface structure imaging method of relief surface, as shown in Figure 1, including following step
It is rapid:
S1, input observation data, migration velocity field, while presetting the initial value of reflectivity model；Here sight
Measured data includes but is not limited to relief surface elevation, observation system, field inspection record (i.e. big gun data) and can be with preset inclined
Shifting parameter etc., wherein offset parameter includes migration velocity field transverse direction sampled point nxAnd longitudinal sampled point nz, spatial sampling interval dx,
Time sampling interval dt, time sampling points nt, dominant frequency f0, big gun number N in the offset data that obtains after samplingsAnd model modification
The threshold value of direction controlling factor o and iteration ends.Common, the initial value of reflectivity model is arranged to 0, migration velocity field
It is to be obtained after carrying out conventional velocity analysis to observation data.
S2 carries out sampling or resampling based on related coefficient to observation data, obtains offset data；Here correlation
Coefficient calculation method are as follows:
Wherein, dj(x, t) andFor the single-shot data and its mean value in last iteration hour offset data, di(x, t) andFor
Single-shot data and its mean value in remaining data.In formula (1), first to di(x,t),dj(x, t) carries out average value processing,
Weaken influence of the Non-zero Mean noise to related coefficient, therefore the formula has certain noise immunity.Then, d is calculated separatelyi
The related coefficient of each single-shot data in (x, t) and last iteration hour offset data.Finally, summing to all related coefficients
To diThe related coefficient of (x, t) and offset data.By formula (1)-(3) it is found that related coefficient can reflect a certain single-shot data with
The correlation of last iteration hour offset data, using related coefficient as with reference to so that develop based on related coefficient such as down-sampling
Strategy:
It is exactly in the 1st iteration, according to fixed interval sampling N when sampling for the first timesIt is a uniform along offset distance direction
The observation data of distribution are as offset data, and wherein the sampling interval is greater than big gun interval and (samples and tie when the sampling interval is equal to big gun interval
Fruit is entire data volume；Sampled result is only part big gun data when only the sampling interval is greater than big gun interval).When the number of iterations is greater than
1 and when meeting resampling condition, the related coefficient of every big gun data and last iteration hour offset data in remaining data is calculated,
Choose the maximum N of related coefficientsOffset data of a single-shot data as current iteration；When the number of single-shot data in remaining data
Amount is less than NsWhen, insufficient part is randomly selected from calculated data.It, directly will be upper when being unsatisfactory for resampling condition
Offset data of an iteration hour offset data as current iteration.It should be noted that remaining data here refers to observation
After the offset data sampled when removing former all iteration in data, remaining data.Resampling condition are as follows: when iteration time
Number k divided by model modification direction controlling factor o remainder be 1 when, to observation data carry out resampling, replace offset data；It is no
Then, without to observation data carry out resampling, directly using last iteration when offset data as the offset data of current iteration.
In with machine gun data sampling, offset data changes at random with the number of iterations, when offset numbers in adjacent iteration twice
When poor according to correlation, convergence curve often generates violent shake, and convergence rate is slower, and computational efficiency is low.And it is based on phase relation
Several big gun data sampling strategies can guarantee that the difference of offset data in adjacent iteration twice is minimum, that is, the sampling policy is not
It is influenced by subsurface structure or surface conditions, offset data can be effectively reduced with the variation of the number of iterations, and then weaken replacement
Offset data improves convergence stability when iteration to the destruction of gradient direction conjugacy and the jitter phenomenon of convergence curve,
And then improve computational efficiency.Meanwhile joined resampling decision condition in the big gun data sampling strategy based on related coefficient,
When being unsatisfactory for resampling condition, offset data is constant, ensures that the conjugacy of gradient direction is set up completely in this way, improves iteration
When convergence efficiency and image quality.
S3 obtains data residual error according to current reflectivity model and offset data；Specifically, for reflection coefficient
Model is simulated by inverse time inverse migration and obtains model data P (x), calculation formula used in inverse time inverse migration simulation model data
Are as follows:
P (x)=ω2∫m(x')f(ω)G0(x'；xs)G0(x；x')dx' (4)
Wherein, P (x) represents the model data that inverse time inverse migration obtains, and m (x') represents reflectivity model, and f (ω) is shake
Source, ω are frequency, G0(x'；xs),G0(x；It x') is the corresponding Green's function of background velocity, xsFor shot position, x, x' are respectively
Underground any point position.
Formula (4) can be solved by two step forward simulations:
Wherein, P0(x') in migration velocity v0(x) the background wave field acquired under, by the product of itself and reflectivity model
It is excited again as focus, the model data of inverse time inverse migration simulation can be acquired.Wave field is carried out according to formula (5) and formula (6) to prolong
It is every to require to apply free boundary condition to wave field by a continuation when opening up:
PΠ(x, z)=0 (7)
Wherein, Π represents the position of relief surface, i.e., the wave field on relief surface is assigned a value of 0.
Based on above description, formula (4) can be simplified shown as with matrix operator form:
D'=Lm (8)
Wherein, m is the matrix form of reflectivity model, and d' is the matrix form of model data, and L is the calculation of inverse time inverse migration
Son, i.e. positive operator in the reverse-time migration of relief surface least square.That is reverse-biased by doing the inverse time to reflectivity model
Model data P (x) can be obtained in movement calculation, in the methods of the invention, it is corresponding only to need to calculate offset data in each iteration
Model data, indicate are as follows:
d’s=Lsm (9)
Wherein, d's,LsThe respectively corresponding model data of offset data and inverse time inverse migration operator, obtain model data
Afterwards, subtract each other with current offset data, obtain data residual error WithRespectively kth time iteration
When offset data and its corresponding inverse time inverse migration operator, m(k-1)Obtained reflectivity model is iterated to calculate for kth -1 time.
S4 determines stochastic gradient and/or random conjugate gradient according to data residual error, based on stochastic gradient and/or total at random
Yoke gradient carries out reflectivity model iteration, obtains final reflectivity model；Specifically, by data residual computations with
The calculation formula of machine gradient are as follows:
Wherein,ForConjugate transposition, i.e. reverse-time migration operator.
The calculation formula of random conjugate gradient modifying factor are as follows:
The then calculation formula of random conjugate gradient are as follows:
Then, reflectivity model iteration is carried out based on stochastic gradient and/or random conjugate gradient, obtains final reflection
Modulus Model, specifically includes the following steps:
S41 judges whether the number of iterations k divided by the remainder of model modification direction controlling factor o is 1, if so, selection with
Machine gradient is model modification gradient, the direction as reflectivity model iteration；If it is not, then selecting random conjugate gradient for model
Update gradient, the direction as reflectivity model iteration.
S42 enters step S44, threshold value here is can be according to warp when the model modification gradient is less than threshold value
The requirement adjustment setting with available accuracy is tested, when the model modification gradient is greater than threshold value, according to model modification gradient
It calculates and updates step-length, update the calculation formula of step-length are as follows:
Wherein,For model modification gradient, including stochastic gradient and/or random conjugate gradient, which is selected in specific iteration
A parameter determines in step S41.
S43, based on step-length is updated, iteration updates reflectivity model, the more new formula of reflectivity model are as follows:
And return step S2,
S44 determines that the reflectivity model obtained at this time is final reflectivity model, by final reflection coefficient mould
Type obtains subsurface structure imaging results.It is avoided using method of the invention by introducing the sampling policy based on related coefficient
Stochastical sampling when being imaged under conditions of relief surface, can be effectively prevented from offset data while reducing offset data amount
Great variety is generated with the number of iterations, and then weakens replacement offset data to the destruction of gradient direction conjugacy and convergence curve
Jitter phenomenon improves computational efficiency when iteration.In addition, the present invention will in Least squares inversion thought introduce relief surface it is inverse
Hour offset can suppress imaging noise, balanced imaging amplitude, improve imaging resolution, high near surface construction and infrastructure
Precision imaging；Using mixing gradient in iterative process, random conjugate gradient has been used to solve in a part of iteration, relative to
The convergence of machine gradient is more stable, and imaging effect is more preferable.
Embodiment 2
In order to further illustrate method of the invention, illustrate side of the invention by taking Canadian overthrust fault model as an example
Method.Detailed process is as shown in figure 12.Input observation data, including migration velocity field (as shown in Figure 2), relief surface elevation are (such as
Shown in Fig. 3), observation system, field inspection data (as shown in Figure 4), the threshold value of iteration ends and default bias parameter, and be arranged
The value of initial reflection Modulus Model is 0, observation system distribution are as follows: be uniformly distributed 278 big guns, every big gun in relief surface with 30 meters of intervals
All be that 556 geophone stations receive, be divided into 15 meters between geophone station, offset parameter is as follows: migration velocity field transverse direction sampled point be 556 and
Longitudinal sampled point is 250, and spatial sampling interval is 15 meters, and time sampling interval is 0.5 millisecond, and time sampling points are 4000
A, dominant frequency is 25 hertz, and the big gun number of offset data is 40 after sampling, and the model modification direction controlling factor is 5；To observation data into
Big gun data sampling of the row based on related coefficient, the observation data after sampling are as shown in Figure 5；For current reflectance model, lead to
Inverse time inverse migration simulation model data are crossed, are subtracted each other with offset data, data residual error is calculated；According to data residual computations stochastic gradient
Descent direction；Judge whether to need to correct stochastic gradient direction, if meeting Rule of judgment, stochastic gradient is modified to random conjugation
Gradient, as model modification direction, otherwise, using stochastic gradient as model modification direction；Whether judgment models more new direction is full
Sufficient threshold condition updates step-length according to model modification direction calculating and updates reflectivity model if being unsatisfactory for threshold condition, then
It is secondary that big gun data sampling, inverse time inverse migration, computation model more new direction based on related coefficient are carried out to observation data, until model
More new direction meets threshold condition, if meeting threshold condition, exports final reflectivity model, i.e., final imaging results
(as shown in Figure 8).
From imaging results, using method of the invention, relative to relief surface reverse-time migration result (such as Fig. 6 institute
Show), imaging noise has been suppressed, the imaging energy in deep is compensated for, near surface structure imaging and medium and deep tomography are more clear, mention
High imaging resolution；(such as relative to the relief surface least square reverse-time migration imaging results based on stochastic gradient descent method
Shown in Fig. 7), imaging noise is weaker, and imaging amplitude is more balanced, is more clear to the imaging of near surface construction.
Fig. 9 be from based on conjugate gradient method relief surface least square reverse-time migration imaging results, be based on stochastic gradient
Distance=1000 meters in the relief surface least square reverse-time migration imaging results of descent method and imaging results of the invention
Locate the single track amplitude curve extracted, imaging results amplitude of the invention and the relief surface least square based on conjugate gradient method are inverse
Hour offset imaging results amplitude is close, relative to the relief surface least square reverse-time migration imaging based on stochastic gradient descent method
As a result amplitude, which has, preferably protects width.
Figure 10 be three kinds of imaging methods normalization data residual error convergence curve, data residual error of the invention with based on conjugation
The data residual error of the relief surface least square reverse-time migration method of gradient method can stable convergence, and based under stochastic gradient
The data residual error of the relief surface least square reverse-time migration imaging method of drop method has obvious fluctuation, and convergent speed is slow
In the method for the present invention.
Figure 11 is that the runing time of several method compares, CPU model Intel (R) Xeon (R) CPU used in imaging test
E5-2650v2@2.60GHz.The calculating time of the image method of the invention is far smaller than the relief surface based on conjugate gradient method most
Small two multiply the calculating time of reverse-time migration method, but the two image quality is suitable, relative to rising based on stochastic gradient descent method
It is close that table least square reverse-time migration method of throwing oneself on the ground calculates the time, but image quality is higher.Therefore, the method for the present invention is opposite and shows
There is technical method, combine higher computational efficiency and image quality, is of great significance to the exploration in western exploratory area.
Claims (5)
1. a kind of subsurface structure imaging method of relief surface, which comprises the following steps:
S1, input observation data, migration velocity field, while presetting the initial value of reflectivity model；
S2 carries out sampling or resampling based on related coefficient to observation data, obtains offset data, specifically include:
When the number of iterations k is 1, according to the fixed sampling interval, N is selectedsIt is a along the equally distributed observation data in offset distance direction
As offset data, and the sampling interval is greater than big gun interval；When the number of iterations k is greater than 1, judge the number of iterations k divided by model more
Whether the remainder of new direction controlling elements o is 1, if so, carrying out resampling to observation data, otherwise, resampling, is not obtained
Offset data, wherein the specific method of resampling include:
The related coefficient for calculating every big gun data and last iteration hour offset data in remaining data, it is maximum to choose related coefficient
NsOffset data of a single-shot data as current iteration；As the lazy weight N of single-shot data in remaining datasWhen, insufficient portion
Divide and is randomly selected from calculated data, its calculation method of the related coefficient are as follows:
Wherein, dj(x, t) andFor the single-shot data and its mean value in last iteration hour offset data, di(x, t) andFor residue
Single-shot data and its mean value in data, NsFor the quantity of single-shot data in offset data, nxFor space sampling number in single-shot data, ntFor time sampling points, i, j represents big gun number；
S3 obtains data residual error according to current reflectivity model and offset data；
S4 determines stochastic gradient and/or random conjugate gradient according to data residual error, based on stochastic gradient and/or random conjugation ladder
Degree carries out reflectivity model iteration, obtains final reflectivity model.
2. a kind of subsurface structure imaging method of relief surface according to claim 1, which is characterized in that the step S4
In the circular of stochastic gradient is determined by data residual error are as follows:WhereinPoint
Not Wei kth time iteration when offset data and its corresponding inverse time inverse migration operator, k be the number of iterations；ForConjugation
Transposition, i.e. reverse-time migration operator；m(k-1)Obtained reflectivity model is iterated to calculate for kth -1 time.
3. a kind of subsurface structure imaging method of relief surface according to claim 2, which is characterized in that the step S4
The circular of the middle random conjugate gradient of determination are as follows:
Wherein random conjugate gradient modifying factorWhereinStochastic gradient is represented,Represent random conjugation ladder
Degree, k represent the number of iterations.
4. a kind of subsurface structure imaging method of relief surface according to claim 3, which is characterized in that the step S4
In reflectivity model iteration carried out based on stochastic gradient and/or random conjugate gradient, obtain final reflectivity model, have
Body the following steps are included:
S41 is based on decision condition, and selecting one in stochastic gradient or random conjugate gradient is model modification gradient, as anti-
It penetrates the direction of Modulus Model iteration, wherein the decision condition specifically includes: judging the number of iterations k divided by model modification direction
Whether the remainder of controlling elements o is 1, if so, selecting stochastic gradient for model modification gradient, if it is not, the then random conjugation of selection
Gradient is model modification gradient；
S42 enters step S44 when the model modification gradient is less than threshold value, when the model modification gradient is greater than threshold
When value, is calculated according to model modification gradient and update step-length；
S43, based on step-length is updated, iteration updates reflectivity model, and return step S2；
S44 determines that the reflectivity model obtained at this time is final reflectivity model.
5. a kind of subsurface structure imaging method of relief surface according to claim 4, which is characterized in that the step
Update step-length in S43, iteration update reflectivity model, circular are as follows:Wherein
m(k)Obtained reflectivity model is iterated to calculate for kth time,α(k)It represents and updates step-length, it is thereinRepresentative model updates gradient, including stochastic gradient or random conjugate gradient, []*Conjugate transposition is represented,Represent kth time repeatedly
For when the corresponding inverse time inverse migration operator of offset data.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710337714.5A CN107193043B (en) | 2017-05-15 | 2017-05-15 | A kind of subsurface structure imaging method of relief surface |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710337714.5A CN107193043B (en) | 2017-05-15 | 2017-05-15 | A kind of subsurface structure imaging method of relief surface |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN107193043A CN107193043A (en) | 2017-09-22 |
| CN107193043B true CN107193043B (en) | 2019-03-29 |
Family
ID=59873489
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710337714.5A Expired - Fee Related CN107193043B (en) | 2017-05-15 | 2017-05-15 | A kind of subsurface structure imaging method of relief surface |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107193043B (en) |
Families Citing this family (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108562937A (en) * | 2018-03-15 | 2018-09-21 | 东北石油大学 | A kind of seismic imaging method |
| CN108845355A (en) * | 2018-09-26 | 2018-11-20 | 中国矿业大学(北京) | Seismic migration imaging method and device |
| CN110531429A (en) * | 2019-08-02 | 2019-12-03 | 中国科学院电子学研究所 | A kind of time-domain electromagnetic data object inversion method based on supervision descent method |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101082676A (en) * | 2007-07-11 | 2007-12-05 | 成都理工大学 | Earthquake post-stack forward method of undulating surface |
| CN101840000A (en) * | 2009-03-20 | 2010-09-22 | 中国石油天然气集团公司 | Plane wave pre-stack depth migration method |
| CN102116869A (en) * | 2011-02-12 | 2011-07-06 | 中国石油大学(华东) | High-precision prestack domain least square migration seismic imaging technology |
| WO2014028030A1 (en) * | 2012-08-17 | 2014-02-20 | Landmark Graphics Corporation | Systems and methods for imaging seismic data |
| CN105301641A (en) * | 2015-10-30 | 2016-02-03 | 中国石油天然气集团公司 | Azimuthal anisotropy speed inversion method and apparatus |
| CN106033124A (en) * | 2016-06-29 | 2016-10-19 | 中国石油化工股份有限公司 | Multi-seismic resource sticky sound least square reverse time migration method based on stochastic optimization |
Family Cites Families (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20120051176A1 (en) * | 2010-08-31 | 2012-03-01 | Chevron U.S.A. Inc. | Reverse time migration back-scattering noise removal using decomposed wavefield directivity |
| US9128205B2 (en) * | 2012-11-13 | 2015-09-08 | Total E&P Usa, Inc. | Process for creating image gathers |
-
2017
- 2017-05-15 CN CN201710337714.5A patent/CN107193043B/en not_active Expired - Fee Related
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101082676A (en) * | 2007-07-11 | 2007-12-05 | 成都理工大学 | Earthquake post-stack forward method of undulating surface |
| CN101840000A (en) * | 2009-03-20 | 2010-09-22 | 中国石油天然气集团公司 | Plane wave pre-stack depth migration method |
| CN102116869A (en) * | 2011-02-12 | 2011-07-06 | 中国石油大学(华东) | High-precision prestack domain least square migration seismic imaging technology |
| WO2014028030A1 (en) * | 2012-08-17 | 2014-02-20 | Landmark Graphics Corporation | Systems and methods for imaging seismic data |
| CN105301641A (en) * | 2015-10-30 | 2016-02-03 | 中国石油天然气集团公司 | Azimuthal anisotropy speed inversion method and apparatus |
| CN106033124A (en) * | 2016-06-29 | 2016-10-19 | 中国石油化工股份有限公司 | Multi-seismic resource sticky sound least square reverse time migration method based on stochastic optimization |
Non-Patent Citations (1)
| Title |
|---|
| "地震偏移成像技术研究现状与发展趋势";李振春 等;《石油地球物理勘探》;20140228;第49卷(第1期);第1-21页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN107193043A (en) | 2017-09-22 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN104597490B (en) | Multi-wave AVO reservoir elastic parameter inversion method based on accurate Zoeppritz equations | |
| FR3031210B1 (en) | SEISMIC ELASTIC WAVE SIMULATION FOR TRANSVERSALLY INCLINED ISOTROPIC ENVIRONMENT USING DECADED LEBEDEV GRID | |
| EP0750203B1 (en) | Subsurface modeling from seismic data and secondary measurements | |
| CN110031896B (en) | Seismic random inversion method and device based on multi-point geostatistics prior information | |
| CN107193043B (en) | A kind of subsurface structure imaging method of relief surface | |
| CN108645994B (en) | Geological random inversion method and device based on multipoint geostatistics | |
| NO340762B1 (en) | Method for building velocity models for pre-stack depth migration via simultaneous joint inversion of seismic, gravitational and magnetotelluric data | |
| KR20090075843A (en) | Iterative inversion of data from simultaneous geophysical sources | |
| CN106033124B (en) | A kind of viscous sound least square reverse-time migration method of more focus based on stochastic optimization | |
| SG184803A1 (en) | Artifact reduction in method of iterative inversion of geophysical data | |
| KR20110000634A (en) | An efficient method for inversion of geophysical data | |
| CN107462924B (en) | A kind of absolute wave impedance inversion method independent of well-log information | |
| EA018815B1 (en) | Forming a geological model | |
| NO326598B1 (en) | Three-dimensional geological modeling | |
| CN106324662B (en) | A kind of full waveform inversion method and system for destination layer | |
| CN109917454B (en) | True earth surface prestack depth migration imaging method and device based on double reference surfaces | |
| AU2015411591A1 (en) | History matching of hydrocarbon production from heterogenous reservoirs | |
| US10324228B2 (en) | Gridless simulation of a fluvio-deltaic environment | |
| CN107229071B (en) | A kind of subsurface structure inversion imaging method | |
| CN106959467B (en) | Seimic wave velocity inversion method and device | |
| CN108318921A (en) | A kind of quick earthquake stochastic inversion methods based on lateral confinement | |
| CN113740901A (en) | Land seismic data full-waveform inversion method and apparatus based on complex undulating surface | |
| RU2570827C2 (en) | Hybrid method for full-waveform inversion using simultaneous and sequential source method | |
| CN104316961B (en) | Method for obtaining geological parameters of weathered layer | |
| CN111597753A (en) | Data depth change characteristic self-adaptive two-dimensional resistivity inversion 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 | ||
| CB03 | Change of inventor or designer information |
Inventor after: Li Chuang Inventor after: Huang Jianping Inventor after: Li Zhenchun Inventor before: Li Chuang Inventor before: Huang Pingjian Inventor before: Li Zhenchun |
|
| CB03 | Change of inventor or designer information | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20190329 Termination date: 20200515 |
|
| CF01 | Termination of patent right due to non-payment of annual fee |