CN104076391A - Partial angle domain anisotropy offset method based on TTI medium four-order travel-time equation - Google Patents

Abstract

The invention discloses a partial angle domain anisotropy offset method based on a TTI medium four-order travel-time equation. The method comprises the steps that a TTI medium four-order time difference term is introduced in step-out time analysis, simplification is carried out with the help of horizontal reflector hypothesis of pre-stack time offset, and a concise travel-time inversion equation is further given out for achieving anisotropy parameter estimation; on this basis, common offset distance domain earthquake offset is expanded to a partial angle domain offset system, accurate description on seismic wave rays is achieved, and a fracture prediction level is improved. The method essentially improves the level of the anisotropy offset technology, more accurate and quantitative extraction and description on anisotropy characteristics (such as speed and amplitude) are achieved, fracture detection can be served better, and fracture inversion quality and accuracy can be improved. Azimuth information of seismic data can be fully utilized, sampling on azimuth is sufficient, the credibility of finally obtained information is improved, and the method is less influenced by noise.

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, relate in particular to a kind of based on TTI(titled transversely isotropy, Transverse Isotropy with a Tilted axis of symmetry, be called for short TTI) the Local angle domain anisotropy offset method of medium quadravalence travel-time equation, main services is in high precision fracture detection and analysis.

Background technology

Affected by the factor such as nonuniformity and anisotropy, it (is Normal Moveout versus azimuth that normal moveout in geological data processing changes with observed bearing conventionally, being called for short NMOz) research and the utilization of .NMOz phenomenon 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 branch.Be subject to the impact of the factors such as the non-average characteristics of underground medium and anisotropic character, geological data is processed normal moveout and is conventionally presented the feature changing with the variation of observed bearing.The people such as Grachka and Tsvanskin has provided the Azimuthal NMO velocity expression formula in 3-D elastic anisotropic media in 1998

$t^2\left(\mathrm{\α}\right)=t_0^2+A_2{x}^2=t_0^2+x^2/V_{\mathrm{nmo}}^2\left(\mathrm{\α}\right)---\left(1\right)$

Wherein,

T ₀zero-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 equally x ²-t ²the approximate NMO velocity of initial slope of relation,

(2)

Under CMP recording geometry, when (formula 1) walked by upward traveling wave when totally walking and when descending ripple walks, carry out respectively two-dimentional Taylor expansion, remain to second order abbreviation, phase adduction substitution equivalence coordinate represents, obtains normal moveout coefficient A ₂expression formula is

$A_2=V_{\mathrm{nmo}}^{-2}\left(\mathrm{\α}\right)=w_{11}{\cos }^2\left(\mathrm{\α}\right)+2w_{12}\sin \left(\mathrm{\α}\right)\cos \left(\mathrm{\α}\right)+w_{22}{\sin }^2\left(\mathrm{\α}\right)---\left(3\right)$

In formula, w ₁₁, w ₁₂, w ₂₂represent three elliptic anisotropy parameters.Observed bearing in above-mentioned formula moves towards orientation as normative reference (being corresponding α=90 °, fracture azimuth) taking vertically oriented fracture.For the purpose of convenient research, interchangeable work taking fracture azimuth as initial 0,

$A_2=V_{\mathrm{nmo}}^{-2}\left(\mathrm{\α}\right)=w_{11}{\sin }^2\left(\mathrm{\α}\right)+2w_{12}\sin \left(\mathrm{\α}\right)\cos \left(\mathrm{\α}\right)+w_{22}{\cos }^2\left(\mathrm{\α}\right)---\left(4\right)$

For thering is the TTI medium of flat reflector, w _ijmathematical expression form be:

$\left\{\begin{array}{c}{w}_{11}=\frac{1}{V_{p0}^2}[1-2\mathrm{\δ}+2{\mathrm{\ϵ}\sin }^{2}v-14(\mathrm{\ϵ}-\mathrm{\δ}){\sin }^2{v\cos }^{2}v]\\ w_{12}=0\\ w_{22}=\frac{1}{V_{p0}^2}[1-2\mathrm{\δ}-2(\mathrm{\ϵ}-\mathrm{\δ}){\sin }^{2}v(1+{\cos }^{2}v)]\end{array}\right.---\left(5\right)$

In formula,

ε and δ represent Thomsen parameter,

V represents fracture dip,

V _p0representative is along the axial velocity of longitudinal wave of symmetry.

Affected by the factor of many possibilities such as rotation of coordinate, in application, conventionally adopt the normal moveout coefficient A2 expression formula (formula 4) of three.

At present, industry member aspect treatment media anisotropic character generally taking point azimuth deviation as main, skew before by prestack CMP(Common Middle Point, being common midpoint gather) road collection data are divided into several sectors and carry out respectively poststack inverting according to gathering orientation, form orientation Dao Ji with this and carry out AVOz(Amplitude Versus Offset and Azimuth, amplitude offset distance and orientation change) crack inverting.The shortcoming of this sector offset maximum is to make full use of the azimuth information of geological data, causes orientation undersampling, reduces the credibility of final acquired information and affected by noise serious.For this reason, having started to advocate gradually in recent years both at home and abroad utilizes the pre-stack time migration of protecting orientation to realize anisotropy skew, this method realizes the continuous utilization to full spectrum information, and in a large amount of practices, prove that the method has retained by anisotropy feature preferably at compacting noise, architectonic impact, thereby can serve better Crack Detection.The technical characterictic of this method is that the time difference of automatic Picking isotropy skew output road collection sets up equation and realize the inverting of 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 that anisotropic parameters is incorporated in skew calculates, and then realizes anisotropy skew.

At present, analyze as main taking second order NMOz in anisotropy skew both at home and abroad, using ground offset distance and the determine standard of collection orientation as output road collection.But, current most NMOz algorithms be to hourage difference Second Order Hyperbolic approximate, ignored the caused non-double curve time difference effect of large offseting distance.In fact, even just can cause whilst on tour to depart from hyperbolic time difference rule under offset distance and reflection horizon depth ratio reach 1 situation.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 conventionally cannot be with the seismic ray approximate representation exciting on ground and receive.If observation offset distance and true incident angle and the position angle of the received seismic event of observed bearing angle of transformation in underground medium with earth's surface will inevitably cause gross error, so that mislead follow-up interpretation and evaluation work.Say in essence, the ground that the cautious orientation of offset distance and big gun is seismic ray characterizes, and Local angle domain is the best scale of weighing underground seismic ray real features.

Summary of the invention

The technical problem to be solved in the present invention is in current industry member, generally to use point azimuth deviation and the skew of common offset anisotropy, this type of offset method all, using orientation, ground and offset distance size as the criteria for classifying, has been lost the accurate description at underground propagation to seismic ray.Meanwhile, in anisotropy skew, use travel-time equation taking second order accuracy as main, the large offseting distance that is difficult to be applicable to day by day increase gathers data.The shortcoming and defect of this two aspect, 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 carry out abbreviation by the flat reflector hypothesis of pre-stack time migration, and then provide simple and clear Travel Time Inversion equation and realize anisotropic parameters and estimate; On this basis, conventional migration technique is extended to Local angle domain skew system apart from territory seismic migration, realization is accurately portrayed seismic ray, improves FRACTURE PREDICTION level.

The method concrete steps are as follows:

Step 1: introduce quadravalence moveout term in anisotropic medium travel-time equation, improve 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 extractive technique, set up Travel Time Inversion equation solution 6 each anisotropic parameters bodies;

Step 4: obtaining on the basis of anisotropic parameters body, calculating four angles of two ray slowness vectors and Local angle domain;

Step 5: by pre-stack time migration framework using angle as the common image gather of determining standard output anisotropy skew.

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 is introduced TTI medium quadravalence moveout term and is carried out abbreviation by the flat reflector hypothesis of pre-stack time migration owing to having adopted 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, from promoting in essence the level of anisotropy migration technology, realize more accurately, quantitatively extracting and portraying anisotropic character (as speed, amplitude etc.), serve better Crack Detection, improve crack inverting quality and order of accuarcy.Can make full use of the azimuth information of geological data, sufficient to orientation sampling, improve the credibility of final acquired information and affected by noise less.

2, the technical program is introduced TTI medium quadravalence moveout term owing to having adopted in step 1, so, solve traditional second order travel-time equation and cannot accurately summarize the shortcoming of large offseting distance rule while walking, again due to the quadravalence moveout equation abbreviation to 3 that has adopted pre-stack time migration, to the hypothesis of flat reflector, script is contained to 20 triangulo operations in step 2, greatly simplify counting yield, significantly improve the practical value in pre-stack time migration system.

3, the technical program is carried out space differentiate owing to having adopted in step 4 to quadravalence travel-time equation, set up incident wave ray vector and reflection wave ray vector, and then by four formation angles of space vector analysis acquisition subsurface imaging place Local angle domain, as the imaging section and the technological means that is total to image gather of determining the skew of standard output Local angle domain anisotropy, so, so, make migration imaging and output road 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:

During about the walking of large offseting distance, Changing Pattern can be traced back to the 60-70 age in last century the earliest, and in order to improve the energy accumulating of velocity spectrum, the people such as Tanner have proposed famous whilst on tour formula

t ²=A ₀+A ₂x ²+A ₄x ⁴+... (8)

Wherein,

T is overall whilst on tour size,

X is geophone offset,

A ₀zero-offset double-pass reflection whilst on tour,

A ₂with A ₄represent respectively second order and quadravalence time difference coefficient.

Above-mentioned formula has been considered seismic event bending deviation effect in the air to a certain extent, has improved the whilst on tour simulation precision of VTI medium, is also one of basic principle of most curved rays Kirchhoff migration algorithm.For anisotropic medium, the quadravalence time difference coefficient A in formula (8) ₄pure mathematics be expressed as

$A_4=\underset{x\→0}{\lim }\frac{1}{2}\frac{d}{d\left(x^2\right)}\left[\frac{d\left(t^2\right)}{d\left(x^2\right)}\right]---\left(9\right)$

The introducing of quadravalence moveout term has improved the simulation precision of equation to non-double curve whilst on tour part while walking, and in this respect, Pech has provided TTI medium A first ₄concrete mathematical form

$A_4^{\mathrm{TTI}}=-2\mathrm{\ηF}(\mathrm{\α},\mathrm{\φ},v)/\left(t_{p0}^2{V}_{p0}^2\right)---\left(10\right)$

In formula,

η=(ε-δ)/(1+2δ)，

T _p0represent the two-way time of zero-offset,

It is initial CMP line position angle that α represents to be inclined to inclined-plane,

φ represents reflecting surface inclination angle,

V represents fracture dip, and

$\begin{array}{c}F(\mathrm{\α},\mathrm{\φ},v)=\frac{1}{128}[18-24\cos 2\mathrm{\α}+6\cos 4\mathrm{\α}+8\cos (6\mathrm{\α}-4v)+4\cos (\mathrm{\α}-2v)\\ -4\cos (4\mathrm{\α}-2v)+24\cos 2(\mathrm{\φ}-2v)+12\cos 2(\mathrm{\α}+\mathrm{\φ}-2v)+8\cos 2(\mathrm{\α}+2\mathrm{\φ}-2v)\\ +4\cos (\mathrm{\α}+3\mathrm{\φ}-2v)+\cos 4(\mathrm{\α}-v)+32\cos 2(\mathrm{\φ}-v)+32\cos 4(\mathrm{\φ}-v)\\ -16\cos 2(\mathrm{\α}+\mathrm{\φ}-v)+8\cos \left(2v\right)+6\cos \left(4v\right)+\cos 4(\mathrm{\α}+v)-4\cos 2(2\mathrm{\α}+v)\\ +4\cos 2(\mathrm{\α}-3\mathrm{\φ}+2v)+8\cos 2(\mathrm{\α}-2\mathrm{\φ}+2v)+12\cos 2(\mathrm{\α}-\mathrm{\φ}+2v)]\end{array}---\left(11\right)$

Can find out the quadravalence time difference coefficient A in TTI medium ₄be by reflecting surface inclination angle, the combined influence of fracture dip and observed azimuth three aspects: forms, and need to carry out trigonometric function operation 20 times, assesses the cost too high, and this is the direct factor that directly limits that it is applied in anisotropy earthquake migration technology.From improving the demand of calculations of offset efficiency, the hypothesis in conjunction with pre-stack time migration to flat reflector, the present invention supposes φ=0 ° on the basis of formula (11), again derives and abbreviation

$F(\mathrm{\α},0,v)=-\frac{2\mathrm{\η}}{t_{p0}^2{V}_{p0}^2}[\frac{14\cos \left(2\mathrm{\α}\right)+{\cos }^2\left(2\mathrm{\α}\right)+17}{32}\cos 4v-\frac{{\cos }^4\mathrm{\α}-1}{2}\cos 2v+\frac{3}{8}{\sin }^4\mathrm{\α}]---\left(12\right)$

So far the anisotropy travel-time equation form that, obtains horizontal layer TTI medium is

$t^2\left(\mathrm{\α}\right)=t_0^2+[c_1\left(\mathrm{\α}\right)w_{11}+c_2\left(\mathrm{\α}\right)w_{12}+c_3\left(\mathrm{\α}\right)w_{12}]x^2+[c_4\left(\mathrm{\α}\right)m_1+c_5\left(\mathrm{\α}\right)m_2+c_6\left(\mathrm{\α}\right)m_3]x^4---\left(13\right)$

In equation, c _1～6represent only relevant to observed bearing coefficient, w ₁₁, w ₁₂, w ₂₂, m ₁, m ₂, m ₃for 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 form, and form altogether 4 angles: the open-angle γ being made up of incident vector Scattering of Vector ₁with open position angle γ ₂, the inclination angle v of ray to normal vector ₁with position angle v ₂.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 ₁,v ₂,γ ₁,γ ₂) (14)

Wherein,

S, R represents ground shot point and geophone station coordinate,

M is the underground coordinate of imaging point,

T is time index;

V ₁, v ₂, γ ₁, γ ₂for four angles of Local angle domain.

Complete Local angle domain mapping process will produce 7 dimension data, need huge computational resource as support.For realize Local angle domain skew under existing design conditions, conventionally Local angle domain to be decomposed into reflected field and direction territory, reflected field is by open-angle γ ₁with open position angle γ ₂composition, direction territory is by inclination angle v ₁with position angle v ₂composition.Wherein the mathematics projected forms of Local angle domain anisotropy skew is

U(S,R,t)→I _γ(M,γ ₁,γ ₂) (15)

From equation (13)s, get respectively the derivative of outward journey time double offset distance and the degree of depth, abbreviation obtains

$\left\{\begin{array}{c}{p}_s=\frac{\∂t_e}{{\∂x}_e}=\frac{1}{t_e}(A_2{x}_e+8A_4{x}_e^3)\\ q_s=\frac{{\∂t}_e}{\∂z}=\frac{t_0}{2t_e{v}_{\mathrm{int}}}\end{array}\right.---\left(16\right)$

In formula,

T _ewith x _erepresent respectively outward journey time and half offset distance,

T ₀for the self excitation and self receiving time,

V _intfor impact point interval velocity.

A ₂with A ₄represent respectively second order and quadravalence time difference coefficient.

Orientation α in conjunction with shot point to imaging point ground subpoint line _s, obtain the incident slowness vector p in local angle character in Fig. 1 _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 character _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 ground subpoint line _s.

Open-angle γ ₁with open position angle γ ₂can utilize space vector further to obtain.

Brief description of the drawings

Fig. 1 is wave ray and the angle schematic diagram of Local angle domain.

Fig. 2 is quadravalence time difference in different offset distance/depth ratio situations schematic diagram that affects on total whilst on tour.

Fig. 3 is that physical model data is introduced schematic diagram.Wherein, (a) model three-dimensional design figure; (b) crack bottom boundary structural map; (c) data position angle

Distribution plan.

Fig. 4 is that seismic migration CRP road, traditional offset distance territory set pair compares schematic diagram.Wherein, (a) isotropy deflection graph; (b) travel based on second order

Time inverting anisotropy deflection graph; (c) the anisotropy deflection graph based on quadravalence Travel Time Inversion.

Fig. 5 is that Local angle domain seismic migration CRP road set pair compares schematic diagram.Wherein, (a) isotropy deflection graph; (b) based on second order whilst on tour

The anisotropy deflection graph of inverting; (c) the anisotropy deflection graph based on quadravalence Travel Time Inversion.

Fig. 6 is the results of fracture prediction schematic diagram of the crack physical model based on the skew of traditional offset distance territory anisotropy.Wherein, (a) crack inverting

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 inverting is tied

Fruit figure; (b) fracture density statistic histogram; (c) fracture azimuth statistic histogram.

Embodiment

The theory advantage that first the present invention is offset by Fig. 1 image explanation Local angle domain anisotropy.As shown in Figure 1, seismic wave propagation angle, the cautious orientation of terrestrial gun and offset distance are only the ground characterization attributes of seismic ray, are difficult to accurately portray the propagation characteristic of seismic ray.In fact (be, incident vector p by two slowness vectors of imaging point _swith Scattering of Vector p _r) Local angle domain that forms is only and weighs the yardstick of underground seismic ray real features, this is also the starting point that solves and distinguish multiray routing problem in traditional earthquake migration technology.

Fig. 2 is for illustrating the necessity of research and development quadravalence travel-time equation.Suppose to have horizontal layer TTI medium parameter to be: ε=0.1, δ=-0.1, fracture azimuth is 40 °, fixing fracture dip is 90 °.Fig. 2 explanation is greater than 1.5 at the ratio of offset distance and the zone of interest degree of depth, and the quadravalence time difference just cannot directly be cast out and do approximate processing.

Crack physical model is an important evidence that proves advantage of the present invention.This physical model is that geophysical survey circle has one of model of extensive popularity aspect crack treatment and inverting.Fig. 3 represents this crack physical model three-dimensional design figure, wherein has three layers of medium, and 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) has shown crack bottom boundary structural map, can see a dome structure, 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 to process and skew thinking, conventionally can only utilize the data in 0-40 ° and 140-180 ° of orientation, the data in other orientation is because degree of covering is rejected compared with major general, and this will be unfavorable for obtaining the azimuthal anisotropy information of fracture medium.And comprehensive skew can make full use of the data in all orientation, be more conducive to keep and disclose the azimuthal anisotropy feature of fracture medium.

The prestack CRP road collection that Fig. 4 (a) extracts at dome structure place for isotropy skew, is affected by fracture medium, and it is different that Slit bottom lineups can observe obvious travel-time difference.Utilize Travel Time Inversion to enter to analyze time difference variation, can obtain the necessary equivalent anisotropic parameters of anisotropic skew.When the TTI quadravalence that equation and the present invention propose in the time utilizing respectively traditional second order to walk 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 that respectively Fig. 4 (b) is with shown in 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 energy focusing is more outstanding.In fact, when traditional second order is walked, equation is generally applicable to the little offset distance situation that is less than 1 with offset distance and the ratio of object layer depth, is difficult to describe medium and feature during compared with the walking of large offseting distance.And the crack bottom boundary degree of depth at this CRP place is in about 1500m, skew output maximum offset reaches 2200m, so residue when equation will inevitably cause walking on large offseting distance while utilizing traditional second order to walk affects energy focus features.

Be analogous to Fig. 4, Fig. 5 further illustrates the impact of travel-time equation precision on the skew of Local angle domain anisotropy.As shown in Figure 5, affected by anisotropy, CRP road collection is extracted in reflected field isotropy skew, and (Fig. 5 a) can observe the seismic event feature that similar " sine " changes.After anisotropy skew, the signal to noise ratio (S/N ratio) of CRP road collection is improved significantly, and waveform character is specification more.Anisotropy skew road collection based on second order Travel Time Inversion (Fig. 5 b) flatness of lineups and focusing than the anisotropy migration result based on quadravalence Travel Time Inversion, (Fig. 5 is c) inferior to some extent.So far, Fig. 4 and Fig. 5 have realized and having absolutely proved the advantage of quadravalence travel-time equation.

Crack inverting based on amplitude characteristic is the importance of evidence Local angle domain anisotropy skew advantage, uses Fig. 6 and Fig. 7 to describe for this reason.The road collection that the two cover crack inversion results of Fig. 6 and Fig. 7 are obtained by the offset distance territory anisotropy skew based on quadravalence travel-time equation and the skew of the Local angle domain anisotropy based on quadravalence travel-time equation respectively carries out inverting gained, the former has represented the hierarchy of skill of domestic and international industry member aspect anisotropy skew substantially, and latter is that technology of the present invention embodies.Than this technological invention, the anisotropy skew of tradition offset distance territory makes that still the fracture density by tectonic structure brought out of residual a great deal of is abnormal at main structure place (as tomography, dome structure place), and from statistical nature, the focusing of its fracture density and fracture azimuth is all obviously not good enough.So far, the present invention accurately, quantitatively extracts and portrays anisotropic character in realization, and the advantage that improves the aspects such as crack inverting quality and order of accuarcy has obtained quantitative and qualitative analysis two aspects and has all been proven.

1, present embodiment is introduced TTI medium quadravalence moveout term and is carried out abbreviation by the flat reflector hypothesis of pre-stack time migration owing to having adopted 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, from promoting in essence the level of anisotropy migration technology, realize more accurately, quantitatively extracting and portraying anisotropic character (as speed, amplitude etc.), serve better Crack Detection, improve crack inverting quality and order of accuarcy.Can make full use of the azimuth information of geological data, sufficient to orientation sampling, improve the credibility of final acquired information and affected by noise less.

2, present embodiment is introduced TTI medium quadravalence moveout term owing to having adopted in step 1, so, solve traditional second order travel-time equation and cannot accurately summarize the shortcoming of large offseting distance rule while walking, again due to the quadravalence moveout equation abbreviation to 3 that has adopted pre-stack time migration, to the hypothesis of flat reflector, script is contained to 20 triangulo operations in step 2, greatly simplify counting yield, significantly improve the practical value in pre-stack time migration system.

3, present embodiment is carried out space differentiate owing to having adopted in step 4 to quadravalence travel-time equation, set up incident wave ray vector and reflection wave ray vector, and then by four formation angles of space vector analysis acquisition subsurface imaging place Local angle domain, as the imaging section and the technological means that is total to image gather of determining the skew of standard output Local angle domain anisotropy, so, (asking inventor to supplement beneficial effect).

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:

During about the walking of large offseting distance, Changing Pattern can be traced back to the 60-70 age in last century the earliest, and in order to improve the energy accumulating of velocity spectrum, the people such as Tanner have proposed famous whilst on tour formula

t ²=A ₀+A ₂x ²+A ₄x ⁴+... (8)

Wherein,

T is overall whilst on tour size,

X is geophone offset,

A ₀zero-offset double-pass reflection whilst on tour,

A ₂with A ₄represent respectively second order and quadravalence time difference coefficient.

Above-mentioned formula has been considered seismic event bending deviation effect in the air to a certain extent, has improved the whilst on tour simulation precision of VTI medium, is also one of basic principle of most curved rays Kirchhoff migration algorithm.For anisotropic medium, the quadravalence time difference coefficient A in formula (8) ₄pure mathematics be expressed as

$A_4=\underset{x\→0}{\lim }\frac{1}{2}\frac{d}{d\left(x^2\right)}\left[\frac{d\left(t^2\right)}{d\left(x^2\right)}\right]---\left(9\right)$

The introducing of quadravalence moveout term has improved the simulation precision of equation to non-double curve whilst on tour part while walking, and in this respect, Pech has provided TTI medium A first ₄concrete mathematical form

$A_4^{\mathrm{TTI}}=-2\mathrm{\ηF}(\mathrm{\α},\mathrm{\φ},v)/\left(t_{p0}^2{V}_{p0}^2\right)---\left(10\right)$

In formula,

η=(ε-δ)/(1+2δ)，

T _p0represent the two-way time of zero-offset,

It is initial CMP line position angle that α represents to be inclined to inclined-plane,

φ represents reflecting surface inclination angle,

V represents fracture dip, and

$\begin{array}{c}F(\mathrm{\α},\mathrm{\φ},v)=\frac{1}{128}[18-24\cos 2\mathrm{\α}+6\cos 4\mathrm{\α}+8\cos (6\mathrm{\α}-4v)+4\cos (\mathrm{\α}-2v)\\ -4\cos (4\mathrm{\α}-2v)+24\cos 2(\mathrm{\φ}-2v)+12\cos 2(\mathrm{\α}+\mathrm{\φ}-2v)+8\cos 2(\mathrm{\α}+2\mathrm{\φ}-2v)\\ +4\cos (\mathrm{\α}+3\mathrm{\φ}-2v)+\cos 4(\mathrm{\α}-v)+32\cos 2(\mathrm{\φ}-v)+32\cos 4(\mathrm{\φ}-v)\\ -16\cos 2(\mathrm{\α}+\mathrm{\φ}-v)+8\cos \left(2v\right)+6\cos \left(4v\right)+\cos 4(\mathrm{\α}+v)-4\cos 2(2\mathrm{\α}+v)\\ +4\cos 2(\mathrm{\α}-3\mathrm{\φ}+2v)+8\cos 2(\mathrm{\α}-2\mathrm{\φ}+2v)+12\cos 2(\mathrm{\α}-\mathrm{\φ}+2v)]\end{array}---\left(11\right)$

Can find out the quadravalence time difference coefficient A in TTI medium ₄be by reflecting surface inclination angle, the combined influence of fracture dip and observed azimuth three aspects: forms, and need to carry out trigonometric function operation 20 times, assesses the cost too high, and this is the direct factor that directly limits that it is applied in anisotropy earthquake migration technology.From improving the demand of calculations of offset efficiency, the hypothesis in conjunction with pre-stack time migration to flat reflector, the present invention supposes φ=0 ° on the basis of formula (11), again derives and abbreviation

$F(\mathrm{\α},0,v)=-\frac{2\mathrm{\η}}{t_{p0}^2{V}_{p0}^2}[\frac{14\cos \left(2\mathrm{\α}\right)+{\cos }^2\left(2\mathrm{\α}\right)+17}{32}\cos 4v-\frac{{\cos }^4\mathrm{\α}-1}{2}\cos 2v+\frac{3}{8}{\sin }^4\mathrm{\α}]---\left(12\right)$

So far the anisotropy travel-time equation form that, obtains horizontal layer TTI medium is

$t^2\left(\mathrm{\α}\right)=t_0^2+[c_1\left(\mathrm{\α}\right)w_{11}+c_2\left(\mathrm{\α}\right)w_{12}+c_3\left(\mathrm{\α}\right)w_{12}]x^2+[c_4\left(\mathrm{\α}\right)m_1+c_5\left(\mathrm{\α}\right)m_2+c_6\left(\mathrm{\α}\right)m_3]x^4---\left(13\right)$

In equation, c _1～6represent only relevant to observed bearing coefficient, w ₁₁, w ₁₂, w ₂₂, m ₁, m ₂, m ₃for 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 form, and form altogether 4 angles: the open-angle γ being made up of incident vector Scattering of Vector ₁with open position angle γ ₂, the inclination angle v of ray to normal vector ₁with position angle v ₂.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 ₁,v ₂,γ ₁,γ ₂) (14)

Wherein,

S, R represents ground shot point and geophone station coordinate,

M is the underground coordinate of imaging point,

T is time index;

V ₁, v ₂, γ ₁, γ ₂for four angles of Local angle domain.

Complete Local angle domain mapping process will produce 7 dimension data, need huge computational resource as support.For realize Local angle domain skew under existing design conditions, conventionally Local angle domain to be decomposed into reflected field and direction territory, reflected field is by open-angle γ ₁with open position angle γ ₂composition, direction territory is by inclination angle v ₁with position angle v ₂composition.Wherein the mathematics projected forms of Local angle domain anisotropy skew is

U(S,R,t)→I _γ(M,γ ₁,γ ₂) (15)

From equation (13)s, get respectively the derivative of outward journey time double offset distance and the degree of depth, abbreviation obtains

$\left\{\begin{array}{c}{p}_s=\frac{\∂t_e}{{\∂x}_e}=\frac{1}{t_e}(A_2{x}_e+8A_4{x}_e^3)\\ q_s=\frac{{\∂t}_e}{\∂z}=\frac{t_0}{2t_e{v}_{\mathrm{int}}}\end{array}\right.---\left(16\right)$

In formula,

T _ewith x _erepresent respectively outward journey time and half offset distance,

T ₀for the self excitation and self receiving time,

V _intfor impact point interval velocity.

A ₂with A ₄represent respectively second order and quadravalence time difference coefficient.

Orientation α in conjunction with shot point to imaging point ground subpoint line _s, obtain the incident slowness vector p in local angle character in Fig. 1 _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 character _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 ground subpoint line _s.

Open-angle γ ₁with open position angle γ ₂can utilize space vector further to obtain.

List of references

Alkhalifah T.,I.Tsvankin,1995,Velocity analysis for transversely isotropic media:Geophysics,60,1550-1566.

Grechka,V.,L.Tsvankin,1998,3-D description of normal moveout in anisotropic inhomogeneous media:Geophysics,63,1079-1092.

Koren Z.,I.Ravve,2011,Full-azimuth subsurface angle domain wavefield decomposition and imaging Part I:Directional and reflection image gathers,Geophysics,76,S1-S13.

Pech A.,I.Tsvankin,V.Grechka,2003,Quartic moveout coefficient:3D description and application to tilted TI media,Geophysics,68,1601-1610.

Sun Z.,D.Wang,X.Zhou,2013,3D amplitude-preserved full-azimuth anisotropic imaging and CRP gathering:An appealing method for integrated fracture detection:77th SEG Annual Meeting,Expanded Abstracts.

Taner M.T.,F.Koehler,1969,Velocity spectra-digital computer derivation and applications of velocity fuctions,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 is protected width skew He Chou road set algorithm and the application in fracture detection thereof: the 5th of geologic geophysical comprehensive rearch centre technology nd Annual Meeting collection.

Claims (3)

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 carry out abbreviation by the flat reflector hypothesis of pre-stack time migration, and then provide simple and clear Travel Time Inversion equation and realize anisotropic parameters and estimate; On this basis, conventional migration technique is extended to Local angle domain skew system apart from territory seismic migration.

2. the Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation according to claim 1, is characterized in that: the step of the method is as follows:

Step 1: introduce quadravalence moveout term in anisotropic medium travel-time equation, improve 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 extractive technique, set up 6 anisotropic parameters bodies of Travel Time Inversion equation solution;

Step 4: obtaining on the basis of anisotropic parameters body, calculating four angles of two ray slowness vectors and Local angle domain;

Step 5: by pre-stack time migration framework using angle as the common image gather of determining standard output anisotropy skew.

3. the Local angle domain anisotropy offset method based on TTI medium quadravalence travel-time equation according to claim 2, it is characterized in that: in step 4, quadravalence travel-time equation is carried out to space differentiate, set up incident wave ray vector and reflection wave ray vector, and then by four formation angles of space vector analysis acquisition subsurface imaging place Local angle domain, as imaging section and the common image gather of determining the skew of standard output Local angle domain anisotropy.