CN104076391B - Based on the Local angle domain anisotropy offset method of TTI medium quadravalence travel-time equation - Google Patents

Based on the Local angle domain anisotropy offset method of TTI medium quadravalence travel-time equation Download PDF

Info

Publication number
CN104076391B
CN104076391B CN201410150887.2A CN201410150887A CN104076391B CN 104076391 B CN104076391 B CN 104076391B CN 201410150887 A CN201410150887 A CN 201410150887A CN 104076391 B CN104076391 B CN 104076391B
Authority
CN
China
Prior art keywords
quadravalence
time
travel
anisotropy
angle domain
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
Application number
CN201410150887.2A
Other languages
Chinese (zh)
Other versions
CN104076391A (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.)
China University of Petroleum Beijing
Original Assignee
Individual
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 Individual filed Critical Individual
Priority to CN201410150887.2A priority Critical patent/CN104076391B/en
Publication of CN104076391A publication Critical patent/CN104076391A/en
Application granted granted Critical
Publication of CN104076391B publication Critical patent/CN104076391B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Landscapes

  • Geophysics And Detection Of Objects (AREA)

Abstract

The invention discloses a kind of Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation, the method is, in step-out time analysis, introduce TTI medium quadravalence moveout term and by pre-stack time migration flat reflector hypothesis carry out abbreviation, and then provide simple and clear Travel Time Inversion equation realize anisotropic parameters estimate; On this basis, conventional migration technique is extended to Local angle domain skew system apart from territory seismic migration, realizes accurately portraying seismic ray, improve FRACTURE PREDICTION level.The method inherently promotes the level of anisotropy migration technology, realize more accurately, quantitatively extracting and portraying anisotropic character (as speed, amplitude etc.), serve Crack Detection better, improve crack inverting quality and order of accuarcy.The azimuth information of geological data can be made full use of, sufficient to azimuth sample, improve the credibility of final acquisition information and affected by noise less.

Description

Based on the Local angle domain anisotropy offset method of TTI medium quadravalence travel-time equation
Technical field
The present invention relates to petroleum exploration field, particularly relate to a kind of based on TTI (titled transversely isotropy, TransverseIsotropywithaTiltedaxisofsymmetry, be called for short TTI) the Local angle domain anisotropy offset method of medium quadravalence travel-time equation, serve primarily in high precision fracture detection and analysis.
Background technology
Affect by the factor such as nonuniformity and anisotropy, normal moveout in seismic data process changes (i.e. NormalMoveoutversusazimuth with observed bearing change usually, being called for short NMOz) research of .NMOz phenomenon and utilization originate from the beginning of the nineties in last century, and become the important foundation in anisotropy research system at present and be widely applied to many research branches.By the impact of the factors such as the non-average characteristics of underground medium and anisotropic character, seismic data process normal moveout presents the feature changed with the change of observed bearing usually.The people such as Grachka and Tsvanskin gave the expression formula of the Azimuthal NMO velocity in 3-D elastic anisotropic media in 1998
t 2 ( α ) = t 0 2 + A 2 x 2 = t 0 2 + x 2 / V n m o 2 ( α ) - - - ( 1 )
Wherein,
T 0zero-offset double-pass reflection whilst on tour,
T is overall whilst on tour size,
X is geophone offset,
V nmofor NMO velocity,
α represents observed azimuth.
In anisotropic medium, take x equally 2-t 2the initial slope of relation is similar to NMO velocity, namely
Under CMP recording geometry, carry out two-dimentional Taylor expansion respectively to when (formula 1) is walked by upward traveling wave when totally walking and when down going wave walks, remain to second order abbreviation, phase adduction substitutes into equivalent coordinate and represent, obtain normal moveout coefficient A 2expression formula is
A 2 = V n m o - 2 ( α ) = w 11 cos 2 ( α ) + 2 w 12 s i n ( α ) c o s ( α ) + w 22 sin 2 ( α ) - - - ( 3 )
In formula, w 11, w 12, w 22represent three elliptic anisotropy parameters.Observed bearing in above-mentioned formula moves towards orientation for normative reference (i.e. corresponding α=90 °, fracture azimuth) with vertically oriented fracture.For the purpose of convenient research, interchangeable work is initial 0 with fracture azimuth, then
A 2 = V n m o - 2 ( α ) = w 11 sin 2 ( α ) + 2 w 12 s i n ( α ) c o s ( α ) + w 22 cos 2 ( α ) - - - ( 4 )
For the TTI medium with flat reflector, w ijmathematical expression form be:
w 11 = 1 V p 0 2 [ 1 - 2 δ + 2 ϵsin 2 v - 14 ( ϵ - δ ) sin 2 v cos 2 v ] w 12 = 0 w 22 = 1 V p 0 2 [ 1 - 2 δ - 2 ( ϵ - δ ) sin 2 v ( 1 + cos 2 v ) ] - - - ( 5 )
In formula,
ε and δ represents Thomsen parameter,
V represents fracture dip,
V p0representative is along the axial velocity of longitudinal wave of symmetry.
Factor by many possibilities such as rotation of coordinate affects, and usually adopts the normal moveout coefficient A2 expression formula (formula 4) of three in the application.
At present, industry member in treatment media anisotropic character generally based on a point azimuth deviation, namely skew before by prestack CMP (CommonMiddlePoint, i.e. common midpoint gather) road collection data are divided into several sector carry out post-stack inversion respectively according to collection orientation, form orientation Dao Ji with this and carry out AVOz (AmplitudeVersusOffsetandAzimuth, amplitude offset distance and Orientation differences) crack inverting.The maximum shortcoming of this sector offset is the azimuth information that cannot make full use of geological data, causes orientation undersampling, reduces the credibility of final acquisition information and seriously affected by noise.For this reason, having started to advocate gradually in recent years both at home and abroad utilizes the pre-stack time migration protecting orientation to realize anisotropy skew, this method realizes the continuous utilization to full spectrum information, and in a large amount of practice, prove that the method is in compacting noise, architectonic impact and remain by anisotropy feature preferably, thus Crack Detection can be served better.The technical characteristic of this method is that the time difference of automatic Picking isotropy skew Efferent tube collection sets up the inverting that equation realizes three anisotropic parameterses
Wherein,
M represents offset distance number, and n represents position angle number,
X ifor offset distance size, for observed azimuth,
T represents reflection wave whilst on tour.
The whilst on tour incorporated by anisotropic parameters in skew calculates, and then realizes anisotropy skew.
At present, analyzes based on second order NMOz in anisotropy skew both at home and abroad, determine standard using ground offset distance and collection orientation as Efferent tube collection.But, current most NMOz algorithm be to hourage difference Second Order Hyperbolic be similar to, ignore the non-double curve time difference effect caused by large offseting distance.In fact, even if whilst on tour just can be caused to depart from hyperbolic time difference rule when offset distance and reflection horizon depth ratio reach 1.On the other hand, seismic event will experience the effects such as complicated distortion deviation in the communication process of underground medium, makes the ray under arrival point usually cannot with seismic ray approximate representation ground excited and receive.If replace the true incident angle of received seismic event in underground medium and position angle to cause gross error with the observation offset distance on earth's surface and observed azimuth, so that mislead following explanations appraisal.Essentially, orientation examined by offset distance and big gun is that the ground of seismic ray characterizes, and Local angle domain is then the best scale weighing subsurface seismic ray real features.
Summary of the invention
The technical problem to be solved in the present invention generally uses point azimuth deviation and the skew of common offset anisotropy in current industry member, this type of offset method all using orientation, ground and offset distance size as the criteria for classifying, loses the accurate description of seismic ray at underground propagation.Meanwhile, use travel-time equation based on second order accuracy in anisotropy skew, the large offseting distance being difficult to be applicable to day by day increase gathers data.The shortcoming and defect of these two aspects, makes the signal to noise ratio (S/N ratio) of follow-up crack inverting and credibility have a greatly reduced quality to a great extent.
In order to solve the problems of the technologies described above, a kind of Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation of the present invention is, in step-out time analysis, introduce TTI medium quadravalence moveout term and by pre-stack time migration flat reflector hypothesis carry out abbreviation, and then provide simple and clear Travel Time Inversion equation realize anisotropic parameters estimate; On this basis, conventional migration technique is extended to Local angle domain skew system apart from territory seismic migration, realizes accurately portraying seismic ray, improve FRACTURE PREDICTION level.
The method concrete steps are as follows:
Step 1: introduce quadravalence moveout term in anisotropic medium travel-time equation, improves the fitting precision to large offseting distance or wide-angle whilst on tour;
Step 2: towards actual computation demand, further abbreviation quadravalence moveout term, reduces redundant computation amount;
Step 3: adopt automatic time difference receiver technology, sets up each anisotropic parameters body of Travel Time Inversion equation solution 6;
Step 4: on the basis obtaining anisotropic parameters body, calculate four angles of two ray slowness vectors and Local angle domain;
Step 5: by pre-stack time migration framework using angle as the common imaging gather of determining the skew of standard output anisotropy.
Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation of the present invention compared with prior art has following beneficial effect.
1, the technical program introduces TTI medium quadravalence moveout term and flat reflector hypothesis by pre-stack time migration carries out abbreviation owing to have employed in step-out time analysis, and then provides simple and clear Travel Time Inversion equation and realize anisotropic parameters and estimate; On this basis, conventional migration technique is extended to the technological means of Local angle domain skew system apart from territory seismic migration, so, inherently promote the level of anisotropy migration technology, realize more accurately, quantitatively extracting and portraying anisotropic character (as speed, amplitude etc.), serve Crack Detection better, improve crack inverting quality and order of accuarcy.The azimuth information of geological data can be made full use of, sufficient to azimuth sample, improve the credibility of final acquisition information and affected by noise less.
2, the technical program introduces TTI medium quadravalence moveout term in step 1 owing to have employed, so, solve the shortcoming that traditional second order travel-time equation accurately cannot summarize rule when large offseting distance is walked, again owing to have employed pre-stack time migration in step 2 to the hypothesis of flat reflector by the quadravalence moveout equation abbreviation to 3 originally containing 20 triangulo operations, greatly simplify counting yield, significantly improve the practical value in pre-stack time migration system.
3, the technical program carries out space differentiate to quadravalence travel-time equation in step 4 owing to have employed, set up incident wave ray vector and reflection wave ray vector, and then four formation angles of subsurface imaging place Local angle domain are obtained by space vector analysis, as the technological means of the imaging section and common imaging gather of determining the skew of standard output Local angle domain anisotropy, so, make migration imaging and Efferent tube collection quality be improved significantly.
The ultimate principle of the Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation provided by the present invention is as follows:
Can trace back to the 60-70 age in last century the earliest about Changing Pattern during the walking of large offseting distance, in order to improve the energy accumulating of velocity spectrum, the people such as Tanner propose famous whilst on tour formula
t 2=A 0+A 2x 2+A 4x 4+...(8)
Wherein,
T is overall whilst on tour size,
X is geophone offset,
A 0zero-offset double-pass reflection whilst on tour,
A 2with A 4represent second order and quadravalence time difference coefficient respectively.
Above-mentioned formula considers seismic event bending Diffraction dispersion in the air to a certain extent, improves the whilst on tour simulation precision of VTI medium, is also one of basic principle of current most of curved rays Kirchhoff migration algorithm.For anisotropic medium, the quadravalence time difference coefficient A in formula (8) 4pure mathematics be expressed as
A 4 = lim x → 0 1 2 d d ( x 2 ) [ d ( t 2 ) d ( x 2 ) ] - - - ( 9 )
When the introducing of quadravalence moveout term improves away, equation is to the simulation precision of non-double curve whilst on tour part, and in this respect, Pech gives TTI medium A first 4concrete mathematical form
A 4 T T I = - 2 η F ( α , φ , v ) / ( t p 0 2 V p 0 2 ) - - - ( 10 )
In formula,
η=(ε-δ)/(1+2δ),
T p0represent the two-way time of zero-offset,
α representative is initial CMP line position angle to be inclined to inclined-plane,
φ represents reflecting surface inclination angle,
V represents fracture dip, and
F ( α , φ , v ) = 1 128 [ 18 - 24 cos 2 α + 6 cos 4 α + 8 cos ( 6 α - 4 v ) + 4 cos ( α - 2 v ) - 4 cos ( 4 α - 2 v ) + 24 cos 2 ( φ - 2 v ) + 12 cos 2 ( α + φ - 2 v ) + 8 cos 2 ( α + 2 φ - 2 v ) + 4 cos ( α + 3 φ - 2 v ) + cos 4 ( α - v ) + 32 cos 2 ( φ - v ) + 32 cos 4 ( φ - v ) - 16 cos 2 ( α + φ - v ) + 8 cos ( 2 v ) + 6 cos ( 4 v ) + cos 4 ( α + v ) - 4 cos 2 ( 2 α + v ) + 4 cos 2 ( α - 3 φ + 2 v ) + 8 cos ( α - 2 φ + 2 v ) + 12 cos 2 ( α - φ + 2 v ) ] - - - ( 11 )
Can find out, the quadravalence time difference coefficient A in TTI medium 4be by reflecting surface inclination angle, the combined influence of fracture dip and observed azimuth three aspect is formed, and needs to carry out 20 trigonometric function operations, assesses the cost too high, and this is the most direct factor directly limiting that it is applied in anisotropy earthquake migration technology.From the demand improving calculations of offset efficiency, in conjunction with the hypothesis of pre-stack time migration to flat reflector, the present invention supposes φ=0 ° on the basis of formula (11), again derives and abbreviation
F ( α , 0 , v ) = - 2 η t p 0 2 V p 0 2 [ 14 cos ( 2 α ) + cos 2 ( 2 α ) + 17 32 cos 4 v - cos 4 α - 1 2 cos 2 v + 3 8 sin 4 α ] - - - ( 12 )
So far, the anisotropy travel-time equation form obtaining horizontal layer TTI medium is
t 2 ( α ) = t 0 2 + [ c 1 ( α ) w 11 + c 2 ( α ) w 12 + c 3 ( α ) w 12 ] x 2 + [ c 4 ( α ) m 1 + c 5 ( α ) m 2 + c 6 ( α ) m 3 ] x 4 - - - ( 13 )
In equation, c 1 ~ 6represent only relevant to observed bearing coefficient, w 11, w 12, w 22, m 1, m 2, m 3for calculating the whilst on tour parameter of TTI medium.This equation is one of basic theory of the present invention.
On the other hand, Local angle domain is in essence by incident vector p swith Scattering of Vector p rboth are formed, and form 4 angles altogether: the open-angle γ be made up of incident vector Scattering of Vector 1with open position angle γ 2, ray is to the inclination angle v of normal vector 1with position angle v 2.In theory, by the seismic data mapping of ground acquisition to the mathematical notation of Local angle domain be
U(S,R,t)→I(M,v 1,v 212)(14)
Wherein,
S, R represent ground shot point and geophone station coordinate,
M is imaging point underground coordinate,
T is time index;
V 1, v 2, γ 1, γ 2for four angles of Local angle domain.
Complete Local angle domain mapping process, by generation 7 dimension data, needs huge computational resource as support.For realizing Local angle domain skew under existing design conditions, usually Local angle domain will be decomposed into reflected field and directions, namely reflected field is by open-angle γ 1with open position angle γ 2composition, directions is by inclination angle v 1with position angle v 2composition.Wherein the mathe-matical map form of Local angle domain anisotropy skew is
U(S,R,t)→I γ(M,γ 12)(15)
From equation, (13)s, get the derivative of outward journey time double offset distance and the degree of depth respectively, abbreviation obtains
p s = ∂ t e ∂ x e = 1 t e ( A 2 x e + 8 A 4 x e 3 ) q s = ∂ t e ∂ z = t 0 2 t e v int - - - ( 16 )
In formula,
T ewith x erepresent outward journey time and half offset distance respectively,
T 0for the self excitation and self receiving time,
V intfor impact point interval velocity.
A 2with A 4represent second order and quadravalence time difference coefficient respectively.
In conjunction with shot point to the orientation α of imaging point floor projection point line s, obtain the incident slowness vector p in Fig. 1 in local angle feature sdefinition
p s=(p ssinα s,p scosα s,q s)(17)
Similar thinking is derived, and obtains the reflection slowness vector p in local angle feature rbe defined as
p r=(p rsinα r,p rcosα r,q r)(18)
Wherein,
Be combined into the orientation α of picture point to geophone station floor projection point line s.
Open-angle γ 1with open position angle γ 2space vector can be utilized to obtain further.
Accompanying drawing explanation
Fig. 1 is wave ray and the angle schematic diagram of Local angle domain.
Fig. 2 affects schematic diagram to total whilst on tour at quadravalence time difference in different offset distance/depth ratio situation.
Fig. 3 is that physical model data introduces schematic diagram.Wherein, (a) model three-dimensional design figure; (b) crack bottom boundary structural map; (c) data azimuthal distribution figure.
Fig. 4 is that the CRP road set pair of conventional offset distance territory seismic migration compares schematic diagram.Wherein, (a) isotropy deflection graph; B () is based on the anisotropy deflection graph of second order Travel Time Inversion; C () is based on the anisotropy deflection graph of quadravalence Travel Time Inversion.
Fig. 5 is that the CRP road set pair of Local angle domain seismic migration compares schematic diagram.Wherein, (a) isotropy deflection graph; B () is based on the anisotropy deflection graph of second order Travel Time Inversion; C () is based on the anisotropy deflection graph of quadravalence Travel Time Inversion.
Fig. 6 is based on the results of fracture prediction schematic diagram of conventional offset apart from the crack physical model of territory anisotropy skew.Wherein, (a) crack inversion result figure; (b) fracture density statistic histogram; (c) fracture azimuth statistic histogram.
Fig. 7 is the results of fracture prediction schematic diagram of the crack physical model based on the skew of Local angle domain anisotropy.Wherein, (a) crack inversion result figure; (b) fracture density statistic histogram; (c) fracture azimuth statistic histogram.
Embodiment
First the present invention illustrates the theory advantage of Local angle domain anisotropy skew by Fig. 1 image.As shown in Figure 1, seismic wave propagation angle, the ground characterization attributes that orientation and offset distance are only seismic rays examined by terrestrial gun, is difficult to the propagation characteristic accurately portraying seismic ray.In fact, by two slowness vectors (i.e. incident vector p of imaging point swith Scattering of Vector p r) Local angle domain that forms is only the yardstick weighing subsurface seismic ray real features, this is also the starting point solving and distinguish multiray routing problem in traditional earthquake migration technology.
Fig. 2 is for illustration of the necessity of research and development quadravalence travel-time equation.Suppose to have horizontal layer TTI medium parameter be: ε=0.1, δ=-0.1, fracture azimuth is 40 °, and fixing fracture dip is 90 °.Fig. 2 illustrates that the quadravalence time difference does approximate processing with regard to cannot directly cast out when the ratio of offset distance and the zone of interest degree of depth is greater than 1.5.
Crack physical model is the important evidence proving advantage of the present invention.This physical model is geophysical survey circle one of model with extensive popularity in crack treatment and inverting.Fig. 3 represents this crack physical model three-dimensional design figure, and wherein have three layers of medium, the second layer is fracture medium, and fracture strike is designed to 90 °, and all the other are two-layer is isotropic medium.Fig. 3 (b) illustrates crack bottom boundary structural map, can see dome structure, a trap-down and vertical fault structure.Fig. 3 (c) represents the azimuthal distribution feature of physical model data, can find out that this data belongs to narrow orientation data.If adopt conventional point orientation process and skew thinking, usually can only utilize the data in 0-40 ° and 140-180 ° of orientation, the data in other orientation due to degree of covering comparatively major general be rejected, this will be unfavorable for the azimuthal anisotropy information obtaining fracture medium.And comprehensive skew can make full use of the data in all orientation, be more conducive to the azimuthal anisotropy feature keeping and disclose fracture medium.
The prestack CRP road collection that Fig. 4 (a) extracts at dome structure place for isotropy skew, affect by fracture medium, it is different that Slit bottom lineups can observe obvious travel-time difference.Utilize Travel Time Inversion to enter to analyze time difference change, the necessary equivalent anisotropic parameters of anisotropic skew can be obtained.When the TTI quadravalence that equation and the present invention propose when utilizing traditional second order to walk respectively is walked, equation (equation 13) obtains on the basis of anisotropic parameters, two kinds of different anisotropy migration algorithms are applied to physical model data, and the road of extraction is assembled fruit and is seen respectively shown in Fig. 4 (b) Yu Fig. 4 (c).Contrast the integrated fruit in above-mentioned road, two kinds of offset distance territories, can find that Fig. 4 (c) road collection quality is better than Fig. 4 (b), the lineups of its large offseting distance part are more straight, and Voice segment is more outstanding.In fact, when traditional second order is walked, equation is generally applicable to the little offset distance situation being less than 1 with offset distance and the ratio of object layer depth, is difficult to describe medium and compared with feature during the walking of large offseting distance.And the crack bottom boundary degree of depth at this CRP place is at about 1500m, skew exports maximum offset and reaches 2200m, so residue when equation will inevitably cause walking on large offseting distance when utilizing traditional second order to walk, affects Voice segment feature.
Be analogous to Fig. 4, Fig. 5 further illustrates the impact that travel-time equation precision offsets Local angle domain anisotropy.As shown in Figure 5, affect by anisotropy, CRP road collection is extracted in the skew of reflected field isotropy, and (Fig. 5 a) can observe the seismic event feature that similar " sine " change.After anisotropy skew, the signal to noise ratio (S/N ratio) of CRP road collection is improved significantly, and waveform character is specification more.Based on the flatness of anisotropy skew road collection (Fig. 5 b) lineups of second order Travel Time Inversion and focusing more inferior to some extent than the anisotropy migration result (Fig. 5 c) based on quadravalence Travel Time Inversion.So far, Fig. 4 and Fig. 5 has achieved the advantage of quadravalence travel-time equation and has absolutely proved.
Crack inverting based on amplitude characteristic is the importance of evidence Local angle domain anisotropy skew advantage, uses Fig. 6 and Fig. 7 to be described for this reason.The two cover crack inversion results of Fig. 6 and Fig. 7 carry out inverting gained by the road collection of the offset distance territory anisotropy skew based on quadravalence travel-time equation and the Local angle domain anisotropy skew acquisition based on quadravalence travel-time equation respectively, the former substantially represent the hierarchy of skill of domestic and international industry member in anisotropy skew, and latter is that technology of the present invention embodies.Compared to this technological invention, conventional offset apart from the skew of territory anisotropy make major constituents place (as tomography, dome structure place) still residual a great deal of by tectonic structure the fracture density of bringing out abnormal, and from statistical nature, the focusing of its fracture density and fracture azimuth is all obviously not good enough.So far, the present invention is realizing accurately, quantitatively extracting anisotropic character and portraying, and the advantage improving the aspect such as crack inverting quality and order of accuarcy obtains quantitative and qualitative analysis two aspect and is all proven.
1, present embodiment introduces TTI medium quadravalence moveout term and flat reflector hypothesis by pre-stack time migration carries out abbreviation owing to have employed in step-out time analysis, and then provides simple and clear Travel Time Inversion equation and realize anisotropic parameters and estimate; On this basis, conventional migration technique is extended to the technological means of Local angle domain skew system apart from territory seismic migration, so, inherently promote the level of anisotropy migration technology, realize more accurately, quantitatively extracting and portraying anisotropic character (as speed, amplitude etc.), serve Crack Detection better, improve crack inverting quality and order of accuarcy.The azimuth information of geological data can be made full use of, sufficient to azimuth sample, improve the credibility of final acquisition information and affected by noise less.
2, present embodiment introduces TTI medium quadravalence moveout term in step 1 owing to have employed, so, solve the shortcoming that traditional second order travel-time equation accurately cannot summarize rule when large offseting distance is walked, again owing to have employed pre-stack time migration in step 2 to the hypothesis of flat reflector by the quadravalence moveout equation abbreviation to 3 originally containing 20 triangulo operations, greatly simplify counting yield, significantly improve the practical value in pre-stack time migration system.
3, present embodiment carries out space differentiate to quadravalence travel-time equation in step 4 owing to have employed, set up incident wave ray vector and reflection wave ray vector, and then four formation angles of subsurface imaging place Local angle domain are obtained by space vector analysis, as the technological means of the imaging section and common imaging gather of determining the skew of standard output Local angle domain anisotropy, so, make migration imaging and Efferent tube collection quality be improved significantly.
The ultimate principle of the Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation provided by the present invention is as follows:
Can trace back to the 60-70 age in last century the earliest about Changing Pattern during the walking of large offseting distance, in order to improve the energy accumulating of velocity spectrum, the people such as Tanner propose famous whilst on tour formula
t 2=A 0+A 2x 2+A 4x 4+...(8)
Wherein,
T is overall whilst on tour size,
X is geophone offset,
A 0zero-offset double-pass reflection whilst on tour,
A 2with A 4represent second order and quadravalence time difference coefficient respectively.
Above-mentioned formula considers seismic event bending Diffraction dispersion in the air to a certain extent, improves the whilst on tour simulation precision of VTI medium, is also one of basic principle of current most of curved rays Kirchhoff migration algorithm.For anisotropic medium, the quadravalence time difference coefficient A in formula (8) 4pure mathematics be expressed as
A 4 = lim x → 0 1 2 d d ( x 2 ) [ d ( t 2 ) d ( x 2 ) ] - - - ( 9 )
When the introducing of quadravalence moveout term improves away, equation is to the simulation precision of non-double curve whilst on tour part, and in this respect, Pech gives TTI medium A first 4concrete mathematical form
A 4 T T I = - 2 η F ( α , φ , v ) / ( t p 0 2 V p 0 2 ) - - - ( 2 )
In formula,
η=(ε-δ)/(1+2δ),
T p0represent the two-way time of zero-offset,
α representative is initial CMP line position angle to be inclined to inclined-plane,
φ represents reflecting surface inclination angle,
V represents fracture dip, and
F ( α , φ , v ) = 1 128 [ 18 - 24 cos 2 α + 6 cos 4 α + 8 cos ( 6 α - 4 v ) + 4 cos ( α - 2 v ) - 4 cos ( 4 α - 2 v ) + 24 cos 2 ( φ - 2 v ) + 12 cos 2 ( α + φ - 2 v ) + 8 cos 2 ( α + 2 φ - 2 v ) + 4 cos ( α + 3 φ - 2 v ) + cos 4 ( α - v ) + 32 cos 2 ( φ - v ) + 32 cos 4 ( φ - v ) - 16 cos 2 ( α + φ - v ) + 8 cos ( 2 v ) + 6 cos ( 4 v ) + cos 4 ( α + v ) - 4 cos 2 ( 2 α + v ) + 4 cos 2 ( α - 3 φ + 2 v ) + 8 cos ( α - 2 φ + 2 v ) + 12 cos 2 ( α - φ + 2 v ) ] - - - ( 11 )
Can find out, the quadravalence time difference coefficient A in TTI medium 4be by reflecting surface inclination angle, the combined influence of fracture dip and observed azimuth three aspect is formed, and needs to carry out 20 trigonometric function operations, assesses the cost too high, and this is the most direct factor directly limiting that it is applied in anisotropy earthquake migration technology.From the demand improving calculations of offset efficiency, in conjunction with the hypothesis of pre-stack time migration to flat reflector, the present invention supposes φ=0 ° on the basis of formula (11), again derives and abbreviation
F ( α , 0 , v ) = - 2 η t p 0 2 V p 0 2 [ 14 cos ( 2 α ) + cos 2 ( 2 α ) + 17 32 cos 4 v - cos 4 α - 1 2 cos 2 v + 3 8 sin 4 α ] - - - ( 12 )
So far, the anisotropy travel-time equation form obtaining horizontal layer TTI medium is
t 2 ( α ) = t 0 2 + [ c 1 ( α ) w 11 + c 2 ( α ) w 12 + c 3 ( α ) w 12 ] x 2 + [ c 4 ( α ) m 1 + c 5 ( α ) m 2 + c 6 ( α ) m 3 ] x 4 - - - ( 13 )
In equation, c 1 ~ 6represent only relevant to observed bearing coefficient, w 11, w 12, w 22, m 1, m 2, m 3for calculating the whilst on tour parameter of TTI medium.This equation is one of basic theory of the present invention.
On the other hand, Local angle domain is in essence by incident vector p swith Scattering of Vector p rboth are formed, and form 4 angles altogether: the open-angle γ be made up of incident vector Scattering of Vector 1with open position angle γ 2, ray is to the inclination angle v of normal vector 1with position angle v 2.In theory, by the seismic data mapping of ground acquisition to the mathematical notation of Local angle domain be
U(S,R,t)→I(M,v 1,v 212)(14)
Wherein,
S, R represent ground shot point and geophone station coordinate,
M is imaging point underground coordinate,
T is time index;
V 1, v 2, γ 1, γ 2for four angles of Local angle domain.
Complete Local angle domain mapping process, by generation 7 dimension data, needs huge computational resource as support.For realizing Local angle domain skew under existing design conditions, usually Local angle domain will be decomposed into reflected field and directions, namely reflected field is by open-angle γ 1with open position angle γ 2composition, directions is by inclination angle v 1with position angle v 2composition.Wherein the mathe-matical map form of Local angle domain anisotropy skew is
U(S,R,t)→I γ(M,γ 12)(15)
From equation, (13)s, get the derivative of outward journey time double offset distance and the degree of depth respectively, abbreviation obtains
p s = ∂ t e ∂ x e = 1 t e ( A 2 x e + 8 A 4 x e 3 ) q s = ∂ t e ∂ z = t 0 2 t e v int - - - ( 16 )
In formula,
T ewith x erepresent outward journey time and half offset distance respectively,
T 0for the self excitation and self receiving time,
V intfor impact point interval velocity.
A 2with A 4represent second order and quadravalence time difference coefficient respectively.
In conjunction with shot point to the orientation α of imaging point floor projection point line s, obtain the incident slowness vector p in Fig. 1 in local angle feature sdefinition
p s=(p ssinα s,p scosα s,q s)(17)
Similar thinking is derived, and obtains the reflection slowness vector p in local angle feature rbe defined as
p r=(p rsinα r,p rcosα r,q r)(18)
Wherein,
Be combined into the orientation α of picture point to geophone station floor projection point line s.
Open-angle γ 1with open position angle γ 2space vector can be utilized to obtain further.
List of references
AlkhalifahT.,I.Tsvankin,1995,Velocityanalysisfortransverselyisotropicmedia:Geophysics,60,1550-1566.
Grechka,V.,L.Tsvankin,1998,3-Ddescriptionofnormalmoveoutinanisotropicinhomogeneousmedia:Geophysics,63,1079-1092.
KorenZ.,I.Ravve,2011,Full-azimuthsubsurfaceangledomainwavefielddecompositionandimagingPartI:Directionalandreflectionimagegathers,Geophysics,76,S1-S13.
PechA.,I.Tsvankin,V.Grechka,2003,Quarticmoveoutcoefficient:3DdescriptionandapplicationtotiltedTImedia,Geophysics,68,1601-1610.
SunZ.,D.Wang,X.Zhou,2013,3Damplitude-preservedfull-azimuthanisotropicimagingandCRPgathering:Anappealingmethodforintegratedfracturedetection:77thSEGAnnualMeeting,ExpandedAbstracts.
TanerM.T.,F.Koehler,1969,Velocityspectra-digitalcomputerderivationandapplicationsofvelocityfuctions,Geophysics,34,859-881.
Hu Zhiquan, 2009, angle domain pre-stack depth migration and the technical research of Angle Domain Common Image Gather extracting method, master's thesis, Chengdu University of Technology.
Wang Di, Sun Zandong, Yue Hangyu etc., 2013, comprehensive anisotropy Amplitude preserving and take out set algorithm and the application in fracture detection thereof: the 5th, geologic geophysical comprehensive rearch centre technology nd Annual Meeting collection.

Claims (2)

1. the Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation, it is characterized in that: in step-out time analysis, introduce TTI medium quadravalence moveout term and by pre-stack time migration flat reflector hypothesis carry out abbreviation, and then provide simple and clear Travel Time Inversion equation realize anisotropic parameters estimate; On this basis, conventional migration technique is extended to Local angle domain skew system apart from territory seismic migration, and the step of the method is as follows:
Step 1: introduce quadravalence moveout term in anisotropic medium travel-time equation, improves the fitting precision to large offseting distance or wide-angle whilst on tour;
Step 2: towards actual computation demand, further abbreviation quadravalence moveout term, reduces redundant computation amount;
Step 3: adopt automatic time difference receiver technology, sets up Travel Time Inversion equation solution 6 anisotropic parameters bodies;
Step 4: on the basis obtaining anisotropic parameters body, calculate four angles of two ray slowness vectors and Local angle domain;
Step 5: by pre-stack time migration framework using angle as the common imaging gather of determining the skew of standard output anisotropy.
2. the Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation according to claim 1, it is characterized in that: in step 4 space differentiate is carried out to quadravalence travel-time equation, set up incident wave ray vector and reflection wave ray vector, and then obtained four formation angles of subsurface imaging place Local angle domain by space vector analysis, as imaging section and the common imaging gather of determining the skew of standard output Local angle domain anisotropy.
CN201410150887.2A 2014-04-16 2014-04-16 Based on the Local angle domain anisotropy offset method of TTI medium quadravalence travel-time equation Expired - Fee Related CN104076391B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201410150887.2A CN104076391B (en) 2014-04-16 2014-04-16 Based on the Local angle domain anisotropy offset method of TTI medium quadravalence travel-time equation

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201410150887.2A CN104076391B (en) 2014-04-16 2014-04-16 Based on the Local angle domain anisotropy offset method of TTI medium quadravalence travel-time equation

Publications (2)

Publication Number Publication Date
CN104076391A CN104076391A (en) 2014-10-01
CN104076391B true CN104076391B (en) 2015-12-02

Family

ID=51597773

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201410150887.2A Expired - Fee Related CN104076391B (en) 2014-04-16 2014-04-16 Based on the Local angle domain anisotropy offset method of TTI medium quadravalence travel-time equation

Country Status (1)

Country Link
CN (1) CN104076391B (en)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108627871B (en) * 2017-03-24 2020-02-21 中国石油化工股份有限公司 Method for inverting TTI medium crack property parameters
CN111650638B (en) * 2020-05-21 2022-07-05 长江大学 Seismic wave travel time calculation method
CN113917533B (en) * 2020-07-10 2023-04-28 中国石油化工股份有限公司 TI medium double-linkage omnibearing imaging systematic realization method

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102012521A (en) * 2010-10-27 2011-04-13 中国石油化工股份有限公司 Method for detecting pre-stack cracks in seismic reservoir prediction
CN102066978A (en) * 2008-03-24 2011-05-18 雪佛龙美国公司 System and method for migrating seismic data
CN102213769A (en) * 2010-04-07 2011-10-12 中国石油天然气集团公司 Method for determining anisotropic parameters by utilizing data of three-dimensional VSP (Vertical Seismic Profile)
CN102239429A (en) * 2008-12-03 2011-11-09 雪佛龙美国公司 Multiple anisotropic parameter inversion for a tti earth model
CN102759746A (en) * 2011-04-28 2012-10-31 中国石油天然气集团公司 Method for inverting anisotropy parameters using variable offset vertical seismic profile data

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102066978A (en) * 2008-03-24 2011-05-18 雪佛龙美国公司 System and method for migrating seismic data
CN102239429A (en) * 2008-12-03 2011-11-09 雪佛龙美国公司 Multiple anisotropic parameter inversion for a tti earth model
CN102213769A (en) * 2010-04-07 2011-10-12 中国石油天然气集团公司 Method for determining anisotropic parameters by utilizing data of three-dimensional VSP (Vertical Seismic Profile)
CN102012521A (en) * 2010-10-27 2011-04-13 中国石油化工股份有限公司 Method for detecting pre-stack cracks in seismic reservoir prediction
CN102759746A (en) * 2011-04-28 2012-10-31 中国石油天然气集团公司 Method for inverting anisotropy parameters using variable offset vertical seismic profile data

Also Published As

Publication number Publication date
CN104076391A (en) 2014-10-01

Similar Documents

Publication Publication Date Title
CN102540250B (en) Azimuth fidelity angle domain imaging-based fractured oil and gas reservoir seismic exploration method
CN106094029B (en) Utilize the method for offset distance vector piece geological data Predicating Reservoir Fractures
CN102841379B (en) Method for analyzing pre-stack time migration and speed based on common scatter point channel set
WO2016041189A1 (en) Method for evaluating shale gas reservoir and seeking desert area
Husen et al. Induced seismicity during the construction of the Gotthard Base Tunnel, Switzerland: hypocenter locations and source dimensions
CN102841375A (en) Method for tomography velocity inversion based on angle domain common imaging gathers under complicated condition
CN102213769A (en) Method for determining anisotropic parameters by utilizing data of three-dimensional VSP (Vertical Seismic Profile)
CN105093281A (en) Earthquake multi-wave modeling method under inverse framework
CN104216014A (en) Seismic signal frequency division processing method
CN103901465A (en) Design method of holographic three-dimensional seismic prospecting and observing system
CN104459777A (en) Fluid identification method and system based on fluid bulk modulus AVO inversion
CN105277982A (en) Shale total organic carbon content earthquake prediction method
CN103576200A (en) Low signal-to-noise ratio zone shallow wave impedance interface static correction method
CN103257363A (en) Method for detecting inclination angle of fissure in underground fissure type reservoir stratum
CN103149588B (en) Method and system for calculating VTI (Velocity Time Integral) anisotropic parameter by utilizing well seismic calibration
CN104076391B (en) Based on the Local angle domain anisotropy offset method of TTI medium quadravalence travel-time equation
CN102053260B (en) Method for acquiring azimuth velocity of primary wave and method for processing earthquake data
CN104536041B (en) Optimization method of seismological observation system parameters
CN102565852B (en) Angle domain pre-stack offset data processing method aiming to detect oil-gas-bearing property of reservoir
CN103558637B (en) Based on the detection method far away of three component sensor
Zhu et al. Recent applications of turning-ray tomography
CN109085644A (en) True earth's surface imaging method when being walked based on dual-beam
Opara et al. Near-surface seismic velocity model building from first arrival travel-times-a case study from an onshore, Niger Delta Field
Sun et al. Three-Dimensional Urban Subsurface Space Tomography with Dense Ambient Noise Seismic Array
Buia et al. Depth imaging Coil data: Multi azimuthal tomography earth model building and depth imaging the full azimuth Tulip coil project

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
ASS Succession or assignment of patent right

Owner name: SUN XUEKAI

Free format text: FORMER OWNER: SUN ZANDONG

Effective date: 20150605

Owner name: SUN ZANDONG

Effective date: 20150605

C41 Transfer of patent application or patent right or utility model
C53 Correction of patent of invention or patent application
CB03 Change of inventor or designer information

Inventor after: Sun Zandong

Inventor after: Sun Xuekai

Inventor after: Wang Zhaoming

Inventor after: Pan Wenqing

Inventor before: Sun Xuekai

Inventor before: Sun Zandong

Inventor before: Wang Zhaoming

Inventor before: Pan Wenqing

COR Change of bibliographic data

Free format text: CORRECT: INVENTOR; FROM: SUN XUEKAI SUN ZANDONG WANG ZHAOMING PAN WENQING TO: SUN ZANDONG SUN XUEKAI WANG ZHAOMING PAN WENQING

TA01 Transfer of patent application right

Effective date of registration: 20150605

Address after: 102249, China University of Petroleum, 18, Xuefu Road, Beijing, Changping District (Beijing)

Applicant after: Sun Xuekai

Applicant after: Sun Zandong

Address before: 102200 Beijing city Changping District Road No. 18, China University of Petroleum (Beijing)

Applicant before: Sun Zandong

C14 Grant of patent or utility model
GR01 Patent grant
TR01 Transfer of patent right
TR01 Transfer of patent right

Effective date of registration: 20190312

Address after: 102249 18 Fu Xue Road, Changping District, Beijing

Patentee after: China University of Petroleum (Beijing)

Address before: 102249 Xuefu Road, Changping District, Beijing, China University of Petroleum (Beijing)

Co-patentee before: Sun Zandong

Patentee before: Sun Xuekai

CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20151202

Termination date: 20210416