CN101738639B - Method for improving computing precision of rock fracture parameters - Google Patents
Method for improving computing precision of rock fracture parameters Download PDFInfo
- Publication number
- CN101738639B CN101738639B CN2008102267774A CN200810226777A CN101738639B CN 101738639 B CN101738639 B CN 101738639B CN 2008102267774 A CN2008102267774 A CN 2008102267774A CN 200810226777 A CN200810226777 A CN 200810226777A CN 101738639 B CN101738639 B CN 101738639B
- Authority
- CN
- China
- Prior art keywords
- fracture
- coefficient
- density
- rank correlation
- orientation
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Images
Landscapes
- Investigating Strength Of Materials By Application Of Mechanical Stress (AREA)
Abstract
The invention provides a method for improving computing the precision of rock fracture parameters, for performance evaluation and drilling of oil and gas reservoirs in oil exploration. The method comprises the following steps of: working out fracture directions and the fracture density by using more than two seismic attributes for uniformization treatment, computing the coefficients of rank correlation of the fracture directions in pairs to obtain the rank correlation and fracture density coefficient, and if the correlation coefficient is over 0.9, determining that the density or direction is the density or direction with highest precision; otherwise, computing absolute error values of the fracture directions corresponding to all attributes at the well point, and sorting the absolute error values according to an order from small to large; computing the coefficients of rank correlation of the fracture density corresponding to the attributes and the coherency attributes, sorting the coefficients of rank correlation according to an order from small to large, and finally computing the fracture directions with the highest precision. The method can ensure the leading position of basic attributes and correct the final result to varying degrees according to the data quality so as to improve the forecast precision by using the existing information.
Description
Technical field
The present invention relates to petroleum and natural gas exploration, development technique, specifically is a kind of calculating of the assessment of oily reservoir performance, reserves and method that probing provides the raising computing precision of rock fracture parameters of parameter of can be.
Background technology
In the exploration of oil and natural gas, exploitation industry, the layer position of assembling underground oil and gas calls reservoir, stops the make progress plugged zone of loss of oil gas to be called cap rock above the reservoir.Crack in the reservoir can increase the connectedness between the reservoir hole, the permeability of raising reservoir, also can enlarge reservoir space simultaneously; Educate all the more in the crack, and storage and collection performance is good more.There is the crack then can cause the oil gas loss in the cap rock; Educate all the more in the crack, poor more to the shut-off capacity of oil gas.
Therefore in oil-gas exploration, performance history, clear and definite understanding need be arranged to the fracture development status of reservoir and cap rock, could carry out rational evaluation cap rock shutoff condition, reservoir storage and collection performance; Could reasonably calculate the size of reservoir space and reserves, output; Could be to the position of straight well, the direction of horizontal well, the deployment of water injection well etc. reasonably designs.In a word, FRACTURE PREDICTION is extremely important; Mainly comprise the fracture development direction and grow two CALCULATION OF PARAMETERS of density.
Present FRACTURE PREDICTION mainly obtains through earthquake information being analyzed extract.Earthquake information is called seismic properties again, is meant the value with certain physical significance that can from seismic data, obtain, like speed, amplitude, frequency, phase place etc.According to the difference of earthquake flow chart of data processing, can be divided into two kinds to seismic data: a kind of is the seismic data that superposes entirely, is meant all seismic datas are handled simultaneously, finally obtains the folded inclined to one side seismic data of a sets of data body; Another kind is a branch orientation stack seismic data; Be meant according to shot point---the azimuthal difference of geophone station; Geological data is divided into several sectors handles respectively, finally obtain the pairing folded inclined to one side seismic data in different sectors with many sets of data body of different azimuth corner characteristics.
There are two types to utilize earthquake information to carry out the fracture parameters Calculation Method.One type is based on full stack seismic data, is the basis with the coherence analysis, predicts fracture and micro-fracture through the otherness of analyzing between the seismic trace, thus the method for indirect predictions fracture development band position and fracture density.The another kind of branch orientation stack seismic data that is based on; With the azimuthal anisotropy theory is guidance; Through difference and the distribution characteristics of anisotropy earthquake attribute on different azimuth such as analysis speed, amplitude, frequencies, come directly to predict the direction in crack and the method for density.
First kind theoretical method is ripe, and is simple and quick, the correct macroscopic law of Prediction of fracture, but direction that can not Prediction of fracture.Second class methods can be predicted the direction and the density of fracture development on each earthquake bin theoretically, have very big superiority; But because implementation procedure relative complex; And it is high to the quality requirements of seismic data; Receive easily in seismic data acquisition, the processing procedure to produce and the error of accumulation influences, thereby the reliability that predicts the outcome instability, the result that different data obtains has bigger difference; Therefore up to the present still be in the stage of trying to explore to study, do not have large-scale application in commercial production.
Utilize second class methods need set up the mathematical model of subterranean strata.The anisotropic medium that will comprise one group of vertical crack that is arranged in parallel, has horizontal symmetry axis is called Method in Transverse Isotropic Medium (HTI medium).The HTI medium has two planes of symmetry, and one was the vertical plane of axis of symmetry, and one is the plane vertical with axis of symmetry.On previous plane, show as anisotropic character, on a back plane, show as the isotropy characteristic.Because seismic line are deployed on the surface level, add that vertically oriented fracture is the most important in the oil-gas exploration and development, analyze so utilize seismic data to carry out generally the stratum being reduced to the HTI medium when calculate in the crack.
The anisotropic character of HTI medium can be stated as simply: because the crack is to the absorption and the attenuation of primary seismic wave, compressional wave is when the different directions of HTI medium is propagated, and its speed and frequency will change along with the variation of direction.Be parallel to fracture orientation, absorption that compressional wave is suffered and attenuation are minimum, thereby speed is maximum, frequency is the highest; Perpendicular to fracture orientation, absorption that compressional wave is suffered and attenuation are maximum, thereby speed is minimum, frequency is minimum.Simultaneously, because the variation of speed causes the variation of formation wave impedance, with the variation that brings the reflected P-wave amplitude.In general, these attributes with azimuthal variation characteristic aggregate performance for perpendicular or parallel be the ellipse of major axis or short-axis direction in fracture orientation.
Owing to confirm that a centre coordinate is (0,0), become the ellipse of any angle theta need confirm 3 parameters with the X axle, therefore major axis radius a, minor axis radius b and θ need 3 known point coordinates at least.Utilize certain anisotropy property value (speed, amplitude, frequency or other) of the seismic data volume of 3 or 3 above different azimuth; Through directly separating the positive definite system of equations or passing through least square solution overdetermined equation group; The anisotropy that just can obtain this attribute is oval, thus the direction of calculating fracture and density.
Present this fracture parameters computing method generally only adopt a kind of attribute, improve the precision that fracture parameters calculates through the precision that in the Data Processing process, improves this attribute.Even utilized more than one attribute, also be after calculating crack and direction respectively, carry out comprehensive distinguishing through the comparative analysis of achievement map.The weak point of this method is: put undue emphasis on a certain attribute to anisotropic reflection ability; And ignored in collection, the processing procedure inevitably, because various hypothesis and actual production conditions limit the various errors that produce; And do not utilize out of Memory that these errors are remedied, thereby might lose the chance of further raising precision of prediction.
" FRACTURE PREDICTION technology in the Song Nan work area effect analyze " that Zhu Chenghong calendar year 2001 delivers mentioned utilizing and predicted the crack with longitudinal wave propagation speed, whilst on tour, reflection amplitude and AVO slope etc. with the Changing Pattern of observed bearing.Follow-up research is also many; But mainly concentrate on speed and the reflection amplitude information utilized, like " P ripple data inverting crack method and instance " (Malaysia and China are high, 2003); " the crack prediction method research of direction anisotropy medium " (Du Qizhen etc.; 2003), " compact reservoir P ripple azimuthal anisotropy Crack Detection " (Yang Zhenwu, 2004) etc.Occurred in recent years will branch orientation earthquake information and conventional earthquake information unite and carry out the Study on Forecast achievement; (Zhang Zhi lets etc. like " complicated reservoirs earthquake prediction comprehensive solution "; 2006), " use " (leaf is calm, 2007) etc. based on the earthquake Crack Detection technological synthesis of compressional wave data.But above-mentioned all articles all are to the research that analyses item by item of various attributes, propose to utilize multiple information to carry out the technical scheme of predicted correction simultaneously.
Summary of the invention
The object of the present invention is to provide and a kind ofly utilize multiple anisotropy earthquake attribute to come the density and the direction in computing rock crack simultaneously, and the method for the raising computing precision of rock fracture parameters that optimum is wherein proofreaied and correct.
The present invention provides following technical scheme, comprises the steps:
(1) utilizes two or more anisotropy earthquake attribute, adopt usual way to obtain each self-corresponding fracture orientation and fracture density;
Described anisotropy earthquake attribute is meant the seismic properties with azimuthal anisotropy characteristic, comprises that speed, amplitude, frequency, amplitude change (AVO) gradient and AVO intercept with geophone offset.
Said usual way is meant the method for carrying out fracture orientation and density calculation through match anisotropy ellipse.
(2) adopting usual way that the fracture density that obtains is carried out normalization handles;
Described normalization is handled and is meant the data area unanimity that realizes the fracture density that different attribute obtains through arithmetic.
(3) to the pairing fracture orientation of various attributes, carry out coefficient of rank correlation in twos and calculate, try to achieve n (n-1)/2 a fracture orientation coefficient of rank correlation; To the pairing fracture density of various attributes, carry out coefficient of rank correlation in twos and calculate, try to achieve n (n-1)/2 a fracture density coefficient of rank correlation;
Said n is for participating in the attribute number of computing.
Said coefficient of rank correlation is the statistical study index of two kinds of phenomenon rank correlation degree of reflection in the statistics.
(4) respectively fracture orientation coefficient of rank correlation and fracture density coefficient of rank correlation are compared as follows:
If any coefficient of rank correlation is greater than 0.9, then directly the density of any rock fracture of trying to achieve in the pairing two kinds of seismic properties of this coefficient of rank correlation or direction as precision the highest density or direction; If all coefficient of rank correlations carry out next step all less than 0.9 and greater than 0.4;
Described relatively fracture orientation and fracture density are independently of one another.
(5) be standard with known position angle well-log information, the pairing fracture orientation of each attribute that calculating obtains in step (1) is in the absolute error value at place, well point, and according to absolute error value rank order from small to large, what absolute error value was minimum ranked first;
If described rank order has the position angle well-log information of many mouthfuls of wells in the research work area, then press the absolute error value sum ordering of the fracture orientation at each place, well point.
The coherence properties of (6) trying to achieve with known full stack seismic data is a standard; The pairing fracture density of each attribute that calculating obtains in step (2) and the coefficient of rank correlation of coherence properties; According to coefficient of rank correlation rank order from big to small, what coefficient of rank correlation was maximum ranked first;
(7) have the highest fracture orientation S of precision through following formula calculating:
S＝S
0+∑A
i·(S
i-S
0)，i＝1，2，3，……n-1
And A is arranged
1∈ [0,0.5], A
I+1≤A
i, i=1,2,3 ... N-2.(1)
S
0Fracture orientation for ordering the 1st in the step (5);
S
1Fracture orientation for ordering the 2nd in the step (5);
……；
S
N-1Fracture orientation (n is the number of the attribute of participation computing) for preface n in the step (5).
A
1Be the 1st fracture orientation correction coefficient;
A
2Be the 2nd fracture orientation correction coefficient;
……
A
iBe i fracture orientation correction coefficient.
Described definite fracture orientation correction coefficient A
i(i=1,2,3 ... N-1) method is:
1) according to A
1, A
2... A
N-1Order confirm A successively
i
2) confirm A with Substitution method
iValue, make i item A
i(S
i-S
0) adding can make S
0+ ∑ A
i(S
i-S
0), the absolute error value that calculates according to the method for step (5) is for minimum;
3) adding when the i item can not further reduce S
0+ ∑ A
i(S
i-S
0) absolute error value the time, with A
i, A
I+1, A
I+2... A
N-1Value be made as 0.
(8) have the highest fracture density D of precision through following formula calculating:
D＝D
0+∑B
i·(D
i-D
0)，i＝1，2，3，……n-1
And B is arranged
1∈ [0,0.5], B
I+1≤B
i, i=1,2,3 ... N-2.(2)
D
0Fracture density for ordering the 1st in the step (6);
D
1Fracture density for ordering the 2nd in the step (6);
……；
D
N-1Fracture density (n is the number of the attribute of participation computing) for preface n in the step (6).
B
1Be the 1st fracture density correction coefficient;
B
2Be the 2nd fracture density correction coefficient;
……
B
iBe i fracture density correction coefficient.
Described definite fracture density correction coefficient B
i(i=1,2,3 ... N-1) method is:
1) according to B
1, B
2... B
N-1Order confirm B successively
i
2) confirm B with Substitution method
iValue, make i item B
i(D
i-D
0) adding can make D
0+ ∑ B
i(D
i-D
0), the coefficient of rank correlation that calculates according to the method for step (6) is for maximum;
3) as i item B
i(D
i-D
0) adding can not further improve D
0+ ∑ is B not
i(D
i-D
0) coefficient of rank correlation the time, with B
i, B
I+1, B
I+2... B
N-1Value be made as 0.
One aspect of the present invention has guaranteed the leading position of base attribute, can carry out correction in various degree to end result according to data quality again simultaneously, reaches the purpose of utilizing existing information to improve precision of prediction.
The present invention confirms that through the fracture density of certain district's ORDOVICIAN CARBONATE reservoir is calculated the present invention has improved the precision of final calculation result really.Accompanying drawing 1 is the coherence properties planimetric map that utilizes full stack seismic data to obtain.Accompanying drawing 2 and accompanying drawing 3 are respectively to utilize in the branchs orientation stack seismic data anisotropic character of amplitude and speed to try to achieve fracture density, and the coefficient of rank correlation of they and Fig. 1 is for being respectively 0.85 and 0.74.According to step (6), the amplitude attribute ranked first, the speed attribute ranked second.Contrast accompanying drawing 2 and accompanying drawing 1, the compact district of visible accompanying drawing 2 is more clear, explains that the reflection in 2 pairs of cracks of accompanying drawing is more accurate.With accompanying drawing 3 and accompanying drawing 1 and accompanying drawing 2 contrasts, certain that visible velocity anisotropy can embody is unusual respectively, on the amplitude anisotropy with coherence map on all not obvious, explain that the velocity anisotropy is best in this place's effect.Utilize formula (2), get A
1Value be 0.35, trying to achieve final fracture density figure is accompanying drawing 4.The coefficient of rank correlation of it and accompanying drawing 1 is 0.89, makes a farfetched comparison Fig. 2 and slightly improves to some extent.Contrast four figure, visible accompanying drawing 4 combines the advantage of accompanying drawing 2 and 3 two kinds of attributes of accompanying drawing, and has many places anomaly ratio accompanying drawing 1 also obviously, thereby makes that the precision of accompanying drawing 4 is the highest.
Description of drawings
Accompanying drawing 1 conventional seismic data volume coherence properties planimetric map;
The fracture density prediction planimetric map that accompanying drawing 2 utilizes the amplitude anisotropy to obtain;
The fracture density prediction planimetric map that accompanying drawing 3 utilizes the velocity anisotropy to obtain;
The fracture density prediction planimetric map that accompanying drawing 4 adopts technical scheme of the present invention to obtain.
Embodiment
Given data: the full stack of certain district's one cover geological data, the branch orientation stack seismic data that a cover divides four orientation to handle, and the layer bit data T of expression zone of interest position
0
Embodiment one: utilize two kinds of anisotropy attributes to improve computing precision of rock fracture parameters
(1) utilizes two kinds of anisotropy earthquake attributes, adopt usual way to obtain each self-corresponding fracture orientation and fracture density;
The amplitude attribute A and the speed attribute V that choose zone of interest utilize the ellipse fitting method to try to achieve each self-corresponding fracture orientation and fracture density as the anisotropy earthquake attribute of participating in computing.And be designated as S respectively
A, D
A, S
V, D
V, wherein:
S
A---the fracture orientation of obtaining by amplitude attribute A
D
A---the fracture density of obtaining by amplitude attribute A
S
V---the fracture orientation of obtaining by speed attribute V
D
V---the fracture density of obtaining by speed attribute V
Amplitude attribute A is directly by the T of different azimuth
0Data were extracted in corresponding the branch in the orientation stack geological data; Speed attribute V is then directly by T
0The data replacement (because speed=distance/time, and for underground same point, it is identical that the distance of seismic event can be thought, so T
0Size can representation speed size; T0 is big more, and speed is more little; T0 is more little, and speed is big more).
All interpretation softwares with anisotropy explanation function all have and utilize the ellipse fitting method to try to achieve the module of fracture orientation and fracture density, also can design program voluntarily and calculate.
(2) adopting usual way that the fracture density that obtains is carried out normalization handles;
The fracture density D that calculates in the step (1)
AData area be 200～645, D
VData area be 0～500.To D
ACarry out following computing: D
A=(D
A-200) * 1.124, make D
AData area become 0～500, with D
VConsistent.
(3) two kinds of pairing fracture orientations of attribute are carried out coefficient of rank correlation and calculate, try to achieve fracture orientation coefficient of rank correlation 0.43; Two kinds of pairing fracture densities of attribute are carried out coefficient of rank correlation calculate, try to achieve fracture density coefficient of rank correlation 0.51;
For n kind attribute, can obtain n (n-1)/2 coefficient of rank correlation.When n=2, obtain 1 coefficient of rank correlation.
Coefficient of rank correlation is the statistical study index of two kinds of phenomenon rank correlation degree of reflection in the statistics.The formula of computed rank related coefficient is:
Wherein: x
i, y
iBe respectively the ranking compositor size of two kinds of phenomenons,
Be respectively the mean value of the ranking compositor size of two kinds of phenomenons.The span of r is-1 to 1, and absolute value is big more, and the similarity of two kinds of phenomenons is good more.In this example, during the calculating fracture direction, utilized two kinds of attributes of seismic amplitude and speed (phenomenon) to calculate respectively, obtained different values, the coefficient of rank correlation of these two values is 0.43; During calculating fracture density, in like manner, the coefficient of rank correlation of two values is 0.51.
The module that the computed rank related coefficient is all arranged in common seismic data interpretation software.
(4) respectively fracture orientation coefficient of rank correlation and fracture density coefficient of rank correlation are compared as follows:
If any coefficient of rank correlation is greater than 0.9, then directly the density of any rock fracture of trying to achieve in the pairing two kinds of seismic properties of this coefficient of rank correlation or direction as precision the highest density or direction; If all coefficient of rank correlations carry out next step all less than 0.9 and greater than 0.4;
Fracture orientation coefficient of rank correlation (0.43) and fracture density coefficient of rank correlation (0.51) are all less than 0.9 and greater than 0.4, so can carry out next step.
(5) be standard with known position angle well-log information, the pairing fracture orientation of each attribute that calculating obtains in step (1) is in the absolute error value at place, well point, and according to absolute error value rank order from small to large, what absolute error value was minimum ranked first;
If described rank order has the position angle well-log information of many mouthfuls of wells in the research work area, then press the absolute error value sum ordering of the fracture orientation at each place, well point.
2 mouthfuls of wells are arranged in the study area, and the position angle well-log information shows that zone of interest is 42 ° at the angle, fracture azimuth of Jing1Chu, is 138 ° at the position angle of Jing2Chu.S
AValue at Jing1Chu is 48, is 142 in the value of Jing2Chu; S
VIn the value at Jing1Chu is 52, is 133 in the value of Jing2Chu.So S
AThe absolute error value sum of the fracture orientation of locating in all well points is (48-42)+(142-138)=10; S
VThe absolute error value sum of the fracture orientation of locating in all well points is (52-42)+(138-133)=15.So S
AThe 1st, S
VArrange the 2nd.
The coherence properties of (6) trying to achieve with known full stack seismic data is a standard; The pairing fracture density of each attribute that calculating obtains in step (2) and the coefficient of rank correlation of coherence properties; According to coefficient of rank correlation rank order from big to small, what coefficient of rank correlation was maximum ranked first;
In any a seismic data interpretation software, can be easy to obtain the coherence data body, and then obtain the coherence properties of zone of interest through full stack seismic data is carried out coherent computing.
According to the coefficient of rank correlation computing formula, obtain D
AWith the coefficient of rank correlation of coherence properties be 0.85, obtain D
VWith the coefficient of rank correlation of coherence properties be 0.74.So D
AArrange the 1st, D
VArrange the 2nd.
(7) have the highest fracture orientation S of precision through following formula calculating:
S＝S
0+∑A
i·(S
i-S
0)，i＝1，2，3，……n-1
And A is arranged
1∈ [0,0.5], A
I+1≤A
i, i=1,2,3 ... N-2.(1)
S
0Fracture orientation for ordering the 1st in the step (5);
S
1Fracture orientation for ordering the 2nd in the step (5);
……；
S
N-1Fracture orientation (n is the number of the attribute of participation computing) for preface n in the step (5).
A
1Be the 1st fracture orientation correction coefficient;
A
2Be the 2nd fracture orientation correction coefficient;
……
A
iBe i fracture orientation correction coefficient.
Described definite fracture orientation correction coefficient A
i(i=1,2,3 ... N-1) method is:
1) according to A
1, A
2... A
N-1Order confirm A successively
i
2) confirm A with Substitution method
iValue, make i item A
i(S
i-S
0) adding can make S
0+ ∑ A
i(S
i-S
0), the absolute error value that calculates according to the method for step (5) is for minimum;
3) adding when the i item can not further reduce S
0+ ∑ A
i(S
i-S
0) absolute error value the time, with A
i, A
I+1, A
I+2... A
N-1Value be made as 0.
In this example, n=2 has S=S
0+ A
1(S
1-S
0)
Ranking results according to step (5) has S
0=S
A, S
1=S
V
That is: S=S
A+ A
1(S
V-S
A) (2)
Work as A
1=0.44 o'clock, S was 7.78 in the absolute error of Jing1Chu, was 0 in the error of Jing2Chu, and absolute error and be 7.78 reaches minimum; And respectively less than S
AAnd S
VAbsolute error and 10 and 15.
(8) have the highest fracture density D of precision through following formula calculating:
D＝D
0+∑B
i·(D
i-D
0)，i＝1，2，3，……n-1
And B is arranged
1∈ [0,0.5], B
I+1≤B
i, i=1,2,3 ... N-2.(3)
D
0Fracture density for ordering the 1st in the step (6);
D
1Fracture density for ordering the 2nd in the step (6);
……；
D
N-1Fracture density (n is the number of the attribute of participation computing) for preface n in the step (6).
B
1Be the 1st fracture density correction coefficient;
B
2Be the 2nd fracture density correction coefficient;
……
B
iBe i fracture density correction coefficient.
Described definite fracture density correction coefficient B
i(i=1,2,3 ... N-1) method is:
1) according to B
1, B
2... B
N-1Order confirm B successively
i
2) confirm B with Substitution method
iValue, make i item B
i(D
i-D
0) adding can make D
0+ ∑ B
i(D
i-D
0), the coefficient of rank correlation that calculates according to the method for step (6) is for maximum;
3) as i item B
i(D
i-D
0) adding can not further improve D
0+ ∑ is B not
i(D
i-D
0) coefficient of rank correlation the time, with B
i, B
I+1, B
I+2... B
N-1Value be made as 0.
In this example, n=2 has D=D
0+ B
1(D
1-D
0)
Ranking results according to step (6) has D
0=D
A, D
1=D
V
That is: D=D
A+ B
1(D
V-D
A) (4)
When B1=0.35, the coefficient of rank correlation of D and coherence properties is 0.89, reaches maximum, and respectively greater than D
AAnd D
VCoefficient of rank correlation 0.85 and 0.74.
Description through the foregoing description is visible:
(1) one aspect of the present invention has guaranteed the leading position of base attribute, can carry out correction in various degree to end result according to data quality again simultaneously, reaches the purpose of utilizing existing information to improve precision of prediction.This is one of innovative point of the present invention.To use two kinds of attributes to ask the situation of the fracture density D with full accuracy to be example, formula (4) is converted into:
D＝(1-B
1)D
A+B
1D
V (5)
It is thus clear that D is actual is by D
AAnd D
VWeighted sum obtains.D
ABigger with the coefficient of rank correlation of coherence properties, D is described
ACompare D
VMore near truly separating; B
1∈ [0,0.5] has guaranteed D
AContribution to D is not less than D
VContribution to D.Step (4) has guaranteed D simultaneously
AWith D
VHave the similarity greater than 0.4, therefore (5) formula has realized truly separating under the condition of unknown, and the constraint through multiparameter further approaches truly to separate, and has guaranteed that D is unlikely to differ too big with truly separating.
(2) the present invention can be according to B
1Variation control the dynamics of correction, and can be to initial calculation D as a result
AAnd D
VPrecision (by data quality decision) carry out quantitative evaluation.This is two of an innovative point of the present invention.Know by (5) formula: B
1More little, D
BContribution to D is more little.At D
A≠ D
VSituation under, when β=0, D=D
A, promptly think attribute D
AEnough good, need not to proofread and correct.Work as B
1=0.5 o'clock, D=0.5D
A+ 0.5D
V, think that promptly attribute A and attribute B have same confidence level, predicting the outcome in fact at this moment is least reliable.So B
1Can be used as initial predicted result's quantitative evaluation standard.B
1More little, explain that initial value is accurate more, corresponding data quality is also good more; Otherwise detail file inferior quality then, the initial value precision is lower.And before this, we can only judge the result of calculation of certain attribute or bad qualitatively through absolute error or coefficient of rank correlation.
(3) last, other advantage of the present invention is: 1) no matter improved through the precision of prediction of which kind of mode with certain attribute, this method can directly be utilized its result, and might make to predict the outcome through other attributes and be further improved.2) the present invention has carried out quality control to each attribute in step (4), has guaranteed the reliability that predicts the outcome.3) computing of direction S and density D is relatively independent among the present invention, can not influence each other, and can bring into play the advantage of each attribute in particular aspects (direction or density), thereby guarantees that net result is a best results among all results.
Embodiment two: utilize three kinds of anisotropy attributes to improve computing precision of rock fracture parameters
(1) utilize three kinds of anisotropy earthquake attributes: amplitude attribute A, speed attribute V and frequency attribute F adopt the ellipse fitting method to obtain each self-corresponding fracture orientation and fracture density S
A, D
A, S
V, D
V, S
F, D
F, wherein:
S
F---the fracture orientation of obtaining by frequency attribute F
D
F---the fracture density of obtaining by frequency attribute F
(other parameter is with embodiment one, down together)
(2) adopting usual way that the fracture density that obtains is carried out normalization handles;
The fracture density D that calculates in the step (1)
FData area be 100～420, it is carried out following computing: D
F=(D
F-100) * 1.563, make D
FData area become 0～500, with D
VConsistent.
(3) three kinds of pairing fracture orientations of attribute are carried out coefficient of rank correlation in twos and calculate, try to achieve 3 related coefficient: S
AWith S
VBe 0.43, S
AWith S
FBe 0.4, S
VWith S
FBe 0.36; Three kinds of pairing fracture densities of attribute are carried out coefficient of rank correlation in twos calculate, try to achieve 3 fracture density coefficient of rank correlation D
AWith D
VBe 0.51, D
AWith D
FBe 0.44, D
VWith D
FBe 0.42.
(4) respectively fracture orientation coefficient of rank correlation and fracture density coefficient of rank correlation are compared as follows:
If any coefficient of rank correlation is greater than 0.9, then directly the density of any rock fracture of trying to achieve in the pairing two kinds of seismic properties of this coefficient of rank correlation or direction as precision the highest density or direction; If all coefficient of rank correlations carry out next step all less than 0.9 and greater than 0.4;
Fracture orientation S
VWith S
FCoefficient of rank correlation 0.36 less than 0.4, think that this study area is not suitable for adopting this three kinds of attribute calculating fracture directions.
Three fracture density coefficient of rank correlations (0.51,0.44,0.42) can carry out next step all less than 0.9 and greater than 0.4.
(5) according to the coefficient of rank correlation computing formula, obtain D
FWith the coefficient of rank correlation of coherence properties be 0.68.So D
A(0.85) arranges the 1st, D
V(0.74) arranges the 2nd, D
F(0.68) arranges the 3rd.
(6) have the highest fracture density D of precision through following formula calculating:
D＝D
0+∑B
i·(D
i-D
0)，i＝1，2，3，……n-1
And B is arranged
1∈ [0,0.5], B
I+1≤A
i, i=1,2,3 ... N-2.(1)
In this example, n=2 has D=D
0+ B
1(D
1-D
0)+B
2(D
2-D
0)
Ranking results according to step (5) has D
0=D
A, D
1=D
V, D
2=D
F
That is: D=D
A+ B
1(D
V-D
A)+B
2(D
F-D
A) (2)
Work as B
1=0.35, B
2=0.25 o'clock, the coefficient of rank correlation of D and coherence properties was 0.9, reached maximum, and respectively greater than D
A, D
VAnd D
FCoefficient of rank correlation 0.85,0.74 and 0.68.
Comparative example one is visible with the fracture density result of calculation of embodiment two: in embodiment two; The adding of frequency factor only makes the coefficient of rank correlation of D and coherence properties bring up to 0.9 by 0.89, so frequency attribute is very little to the raising contribution of final fracture density computational accuracy.And the step among the embodiment two (4) shows, frequency can not have a positive effect to the raising of crack direction calculating precision.Therefore, consider Effects of Noise, it is just passable to think that this study area only adopts amplitude and two kinds of attributes of speed to improve final fracture parameters computational accuracy.
Claims (7)
1. a method that improves the rock fracture parameters precision is gathered the seismic properties with azimuthal anisotropy characteristic, comprises that speed, amplitude, frequency, amplitude change (AVO) gradient and AVO intercept with geophone offset, is characterized in that also comprising the steps:
(1) utilizes two or more anisotropy earthquake attribute, adopt usual way to obtain each self-corresponding fracture orientation and fracture density;
(2) adopting usual way that the fracture density that obtains is carried out normalization handles;
(3) to the pairing fracture orientation of various attributes, carry out coefficient of rank correlation in twos and calculate, try to achieve n (n-1)/2 a fracture orientation coefficient of rank correlation; To the pairing fracture density of various attributes, carry out coefficient of rank correlation in twos and calculate, try to achieve n (n-1)/2 a fracture density coefficient of rank correlation; Wherein n is for participating in the seismic properties number of calculating;
(4) respectively fracture orientation coefficient of rank correlation and fracture density coefficient of rank correlation are compared as follows:
If any coefficient of rank correlation is greater than 0.9, then directly the density of any rock fracture of trying to achieve in the pairing two kinds of seismic properties of this coefficient of rank correlation or direction as precision the highest density or direction; If all coefficient of rank correlations carry out next step all less than 0.9 and greater than 0.4;
(5) be standard with known position angle well-log information, the pairing fracture orientation of each attribute that calculating obtains in step (1) is in the absolute error value at place, well point, and according to absolute error value rank order from small to large, what absolute error value was minimum ranked first;
The coherence properties of (6) trying to achieve with known full stack seismic data is a standard; The pairing fracture density of each attribute that calculating obtains in step (2) and the coefficient of rank correlation of coherence properties; According to coefficient of rank correlation rank order from big to small, what coefficient of rank correlation was maximum ranked first;
(7) have the highest fracture orientation S of precision through following formula calculating:
S＝S
0+∑A
i·(S
i-S
0)，i＝1，2，3，……n-1
And A is arranged
1∈ [0,0.5], A
I+1≤A
i, i=1,2,3 ... N-2.(1)
S
0Fracture orientation for ordering the 1st in the step (5);
S
1Fracture orientation for ordering the 2nd in the step (5);
……；
S
N-1Fracture orientation (n is the number of the attribute of participation computing) for preface n in the step (5).
A
1Be the 1st fracture orientation correction coefficient;
A
2Be the 2nd fracture orientation correction coefficient;
……
A
iBe i fracture orientation correction coefficient.
(8) have the highest fracture density D of precision through following formula calculating:
D＝D
0+∑B
i·(D
i-D
0)，i＝1，2，3，……n-1
And B is arranged
1∈ [0,0.5], B
I+1≤B
i, i=1,2,3 ... N-2.(2)
D
0Fracture density for ordering the 1st in the step (6);
D
1Fracture density for ordering the 2nd in the step (6);
……；
D
N-1Fracture density (n is the number of the attribute of participation computing) for preface n in the step (6).
B
1Be the 1st fracture density correction coefficient;
B
2Be the 2nd fracture density correction coefficient;
……
B
iBe i fracture density correction coefficient.
2. the method for raising rock fracture parameters precision according to claim 1 is characterized in that the said usual way of step (1) is meant the method for carrying out fracture orientation and density calculation through match anisotropy ellipse.
3. the method for raising rock fracture parameters precision according to claim 1 is characterized in that it is to realize that through arithmetic the data area of the fracture density that different attribute obtains is consistent that the described normalization of step (2) is handled.
4. the method for raising rock fracture parameters precision according to claim 1 is characterized in that described relatively fracture orientation of step (4) and fracture density are independently of one another.
5. the method for raising rock fracture parameters precision according to claim 1; It is characterized in that step (5) if described rank order has the position angle well-log information of many mouthfuls of wells in the research work area, then press the absolute error value sum ordering of the fracture orientation at each place, well point.
6. the method for raising computing precision of rock fracture parameters according to claim 1 is characterized in that the described definite fracture orientation correction coefficient A of step (7)
i(i=1,2,3 ... N-1) method is:
1) according to A
1, A
2... A
N-1Order confirm A successively
i
2) confirm A with Substitution method
iValue, make i item A
i(S
i-S
0) adding can make S
0+ ∑ A
i(S
i-S
0), the absolute error value that calculates according to the method for step (5) is for minimum;
3) adding when the i item can not further reduce S
0+ ∑ A
i(S
i-S
0) absolute error value the time, with A
i, A
I+1, A
I+2... A
N-1Value be made as 0.
7. the method for raising computing precision of rock fracture parameters according to claim 1 is characterized in that the described definite fracture density correction coefficient B of step (8)
i(i=1,2,3 ... N-1) method is:
1) according to B
1, B
2... B
N-1Order confirm B successively
i
2) confirm B with Substitution method
iValue, make i item B
i(D
i-D
0) adding can make D
0+ ∑ B
i(D
i-D
0), the coefficient of rank correlation that calculates according to the method for step (6) is for maximum;
3) as i item B
i(D
i-D
0) adding can not further improve D
0+ ∑ is B not
i(D
i-D
0) coefficient of rank correlation the time, with B
i, B
I+1, B
I+2... B
N-1Value be made as 0.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN2008102267774A CN101738639B (en) | 2008-11-24 | 2008-11-24 | Method for improving computing precision of rock fracture parameters |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN2008102267774A CN101738639B (en) | 2008-11-24 | 2008-11-24 | Method for improving computing precision of rock fracture parameters |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN101738639A CN101738639A (en) | 2010-06-16 |
| CN101738639B true CN101738639B (en) | 2012-05-30 |
Family
ID=42462357
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN2008102267774A Active CN101738639B (en) | 2008-11-24 | 2008-11-24 | Method for improving computing precision of rock fracture parameters |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN101738639B (en) |
Families Citing this family (21)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102033242B (en) * | 2010-10-22 | 2012-09-26 | 中国石油化工股份有限公司 | Deep inclined fractured reservoir earthquake amplitude prediction method |
| CN102012521B (en) * | 2010-10-27 | 2013-03-20 | 中国石油化工股份有限公司 | Method for detecting pre-stack cracks in seismic reservoir prediction |
| CN102253415B (en) * | 2011-04-19 | 2013-03-20 | 中国石油大学(华东) | Method for establishing earthquake response mode based on fracture equivalent medium model |
| CN102749651B (en) * | 2012-07-12 | 2015-07-08 | 恒泰艾普石油天然气技术服务股份有限公司 | Calibration method for fracture logging of reservoir |
| CN103675913B (en) * | 2012-09-04 | 2017-01-25 | 中国石油天然气集团公司 | Detection method for fractured reservoirs |
| CN103399345B (en) * | 2013-08-16 | 2016-08-24 | 中国科学院地质与地球物理研究所 | The investigation method of a kind of buried hill fissure distribution and device |
| CN104749618B (en) * | 2013-12-26 | 2018-06-05 | 中国石油化工股份有限公司 | Mud shale low-angle dip crack poststack probability quantitatively characterizing method |
| CN103869358A (en) * | 2014-02-18 | 2014-06-18 | 中国石油天然气股份有限公司 | Histogram equalization-based fault identification method and equipment |
| CN104714247B (en) * | 2014-04-24 | 2017-02-15 | 中国石油化工股份有限公司 | Pre-stack and post-stack linkage attribute interpretation method |
| CN103926618B (en) * | 2014-04-25 | 2016-09-28 | 中国科学院地质与地球物理研究所 | Based on cylinder matching 3-d seismic data set fracture detection method |
| CN104914465B (en) * | 2015-06-02 | 2017-06-13 | 中国石油天然气股份有限公司 | Volcanic Rock quantitative forecasting technique and device |
| CN105425292B (en) * | 2015-11-20 | 2019-04-12 | 中国石油天然气集团公司 | A kind of oil gas forecasting method and device |
| CN105758361B (en) * | 2016-02-01 | 2018-05-22 | 绍兴文理学院 | A kind of structural plane roughness coefficient anisotropy quantitative evaluation method |
| CN107728204B (en) * | 2016-08-11 | 2019-03-29 | 中国石油化工股份有限公司 | Based on the anisotropic crack prediction method of prestack longitudinal wave and system |
| CN107145633B (en) * | 2017-04-07 | 2018-04-17 | 中国地质大学（武汉） | A kind of Forecasting Methodology of rock fracture network occurrence three-dimensional statistical distribution |
| CN107678063B (en) * | 2017-09-25 | 2019-11-29 | 中国石油集团东方地球物理勘探有限责任公司 | A kind of multi-component converted wave crack prediction method based on Rank correlation |
| CN107956472B (en) * | 2017-11-10 | 2021-01-01 | 中国石油天然气股份有限公司 | Method, device and medium for determining completion mode of fractured carbonate reservoir |
| CN108957530B (en) * | 2018-05-23 | 2019-08-23 | 电子科技大学 | A kind of crack automatic testing method based on Acceleration Algorithm in Seismic Coherence Cube slice |
| CN110702014B (en) * | 2019-11-26 | 2021-04-06 | 苏交科集团股份有限公司 | Crack monitoring device and method for road slope |
| CN111506861B (en) * | 2020-04-24 | 2020-11-13 | 成都捷科思石油天然气技术发展有限公司 | Method for calculating crack strength of favorable region of target layer |
| CN113970785B (en) * | 2020-07-22 | 2024-05-17 | 中国石油化工股份有限公司 | Method and system for predicting development of crack of underground river, storage medium and electronic equipment |
-
2008
- 2008-11-24 CN CN2008102267774A patent/CN101738639B/en active Active
Also Published As
| Publication number | Publication date |
|---|---|
| CN101738639A (en) | 2010-06-16 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN101738639B (en) | Method for improving computing precision of rock fracture parameters | |
| CN105093313B (en) | A kind of Karst-type oil reservoir individual well oil and gas productivity prediction method and device | |
| WO2016041189A1 (en) | Method for evaluating shale gas reservoir and seeking desert area | |
| CN109425896A (en) | Dolomite oil and gas reservoir distribution forecasting method and device | |
| CN105528522A (en) | Method and device for calculating quantity of resources of continuous oil and gas reservoir based on space grids | |
| CN108680951A (en) | A method of judging that Enriching Coalbed Methane depositional control acts on based on earthquake information | |
| CA2823710C (en) | Methods and systems regarding models of underground formations | |
| CN104678434A (en) | Method for predicting storage layer crack development parameters | |
| CN102466815A (en) | Triassic clastic rock petroleum reservoir identification method | |
| CN104047598A (en) | Method for predicating productivity of nonhomogeneity ancient karst carbonate reservoir | |
| CN104834003B (en) | Phase-controlled compression coefficient earthquake prediction method of unconventional tight gas reservoir | |
| CN107807410A (en) | Reservoir geophysical sensitive parameter method for optimizing and system | |
| CN112034521B (en) | Method for predicting overpressure of under-compacted and hydrocarbon-production mixed formation | |
| CN103869362A (en) | Method and equipment for obtaining body curvature | |
| CN107728205B (en) | A kind of Formation pressure prediction method | |
| CN102914799A (en) | Forward modeling method and device for nonequivalent wave field | |
| CN109283577B (en) | Seismic horizon calibration method | |
| CN111308558B (en) | Shale gas horizontal well longitudinal wave time difference correction method | |
| CN116856845A (en) | Tight sandstone gas reservoir horizontal well geosteering method | |
| CN104749618B (en) | Mud shale low-angle dip crack poststack probability quantitatively characterizing method | |
| CN103278852B (en) | Utilize the method for seismic data volume waveform configuration characteristic model predicting oil/gas | |
| CN104820237A (en) | Method for calculating variance cube of stratum | |
| CN105259576A (en) | Hydrocarbon reservoir identification method by means of seismic statistical characteristic | |
| Li et al. | Three-dimensional reservoir architecture modeling by geostatistical techniques in BD block, Jinhu depression, northern Jiangsu Basin, China | |
| CN110439534A (en) | A kind of seismic data leakage prediction technique |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C10 | Entry into substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| C14 | Grant of patent or utility model | ||
| GR01 | Patent grant |