CN102867302B - Core fracture identification method based on three-dimensional image information processing - Google Patents

Abstract

The invention relates to a core fracture identification method based on three-dimensional image information processing, belonging to the technical field of image segmentation. The method comprises the following steps of: scanning a core to be processed by using a CT (computerized tomography) machine, storing a scanned core CT sequence image in a computer, reading the well scanned core CT sequence image by using core three-dimensional image fracture identification software, performing three-dimensional image reconstruction on the core CT sequence image by utilizing a marching cubes algorithm, and performing a series of calculation, getting and screening processed by utilizing the differences among pores, random noise and fracture in morphological features in a three-dimensional data field to achieve the purpose of identifying a core fracture. According to the method provided by the invention, the core fracture can be fast and accurately identified; and the problems of causing error and even mistake due to negligence of communication information of the core fracture between layers in the existing isolate analysis of each core CT sequence image can be solved.

Description

Based on the rock core crack identification method of three-dimensional image information process

Technical field

The present invention relates to a kind of image processing techniques, particularly a kind of rock core crack identification method based on three-dimensional image information process, belongs to technical field of image segmentation.

Background technology

In field of petroleum geology, crack identification is the important content of oil reservoir and reservoir research, significant to the raising of oil and gas development and oil recovery factor.Rock produces mechanical destruction under effect of stress, without the rift structure of obvious displacement in crack.Due to crannied existence, and fracture development, connective good, mud shale and carbonate reservoir just can be made to become valuable reservoir.

The research work of China to rock core crack is carried out very early, the later stage fifties in last century, the geologist of China's oil gas field just starts the feature such as morphological feature and tomography position according to rock core, and identified the crack existed in rock core by experience, this is a kind of recognition methods based on qualitative analysis.Since the eighties in last century, mainly through appearing, rock core, well logging, earthquake, well testing, formation testing pilot production, drilling well, the data such as well logging, comprehensively analyze the crack identified in rock core by means of determination method, well test analysis method, crack Statistics Method and well-log information method etc. in geology way of qualitative analysis, core chamber.By this method, rock core crack can be analyzed qualitatively; Obtain the fracture number in rock core according to statistics, the sxemiquantitative information such as fracture length simultaneously.But said method is all crack existence and the distribution situation of being analyzed whole rock core inside by the feature on some surfaces of observation rock core, does not obtain the real structure of rock core inside.From the end of the nineties in last century so far, increasingly mature along with computer image treatment means, moves towards the pattern recognition stage to the identification in rock core crack from the numerical simulation calculation stage.Present oil geology research department widely uses computer assisted tomography technology, namely also known as computed tomography, CT, by CT machine to scanning containing crannied rock core, obtain rock core CT sequence image, rock core CT sequence image described above just can be utilized when harmless rock core to observe, analyze the real structure of this rock core inside.The features of different material on rock core CT sequence image such as crack, hole and rock there are differences, crack and hole can be extracted according to image segmentation algorithm, the crack in rock core is identified again according to morphological feature, not only can analyze the crack on core surface, and the situation of change of its crack of observation, rock core inside in different aspects can be goed deep into.

Wang Renyi, Ao pioneer (based on rock core crack image zooming-out [J] of image processing techniques. Xinjiang Geology, 2006, (04) .) etc. utilize the Binarization methods of maximum variance between clusters to extract crack from core profile; Superb, Wang Zhengyong (based on the rock core crack extract algorithm [J] of Beamlet. information and electronic engineering, 2010, (02) .) a kind of rock core crack pattern converted based on multiple dimensioned Beamlet is proposed as extraction algorithm; Analogy continue industry (rock core CT sequence specific primers-polymerase chain reaction some problem research [D]. Sichuan University's Dissertation Database, 2011) propose extract crack based on phase place unification algorism.But, said method has plenty of analyzes core surface image, have plenty of and analyze each Zhang Yanxin CT sequence image to identify its crack isolatedly, with only rock core two-dimensional image information, and all CT sequence chart are not combined analysis, have ignored the communication information of rock core crack between the rock core CT sequence image of different two or more than two, so be easy to become long and narrow linear in two dimensional image, in 3 d data field, just isolated and disconnected hole or random noise are mistaken for crack mutually, or will be subcircular in two dimensional image, should be that the target of the part in crack is mistaken in 3 d data field not crack, the recognition result in rock core crack is so often caused to occur deviation.Therefore object of the present invention is just to provide a kind of method, can identify rock core crack accurately and rapidly.

Summary of the invention

Object of the present invention is to overcome defect existing in prior art and deficiency just, provides a kind of rock core crack identification method based on three-dimensional image information process that can identify rock core crack quickly and accurately.The method first uses CT machine to scan pending rock core, obtains rock core CT sequence image; According to marching cube (Marching Cubes, MC) algorithm, three-dimensional image reconstruction is carried out to the rock core CT sequence image obtained again; Then utilize the difference that the feature of the different materials such as crack, hole and rock on rock core CT sequence image exists, by a series of calculating, acquisition and screening process, thus reach the object identifying rock core crack.

For achieving the above object, the present invention adopts following technical scheme to realize.

The present invention is based on the rock core crack identification method of three-dimensional image information process, its invention main points mainly carry out three-dimensional image reconstruction by scanning the rock core CT sequence image obtained, and carry out crack identification and calculating in the 3-D view rebuild, and its ultimate principle is as follows:

Except containing except rock, crack and hole information in CT rock core sequence image, also inevitably comprise a large amount of random noises.Wherein, the feature of rock is that gray level is lower; The feature of hole is that gray level is very high, generally shows as irregular polygon, ellipse, sub-circular; Although the gray level in crack is also very high, it is that rock texture loses a successional plane, so different from hole, it is long and narrow linear or dendritic distribution in two dimensional image.The gray level that random noise is then unfixing, generally has area little, the features such as discrete distribution in two dimensional image.In order to obtain quantitative information from rock core CT sequence image, an important step will be split rock core CT sequence image, exactly to isolating hole target and crack target from the rock core CT sequence image of reality.But inevitably picture noise is extracted simultaneously based on the partitioning algorithm of half-tone information.Current method usually differently to be distinguished according to hole and random noise and crack are modal in two dimensional image.But in the target of long and narrow linear in two dimensional image, probably just isolated and disconnected hole or random noise mutually in the 3 d data field of reality, so it is not crack.The target of described sub-circular is also likely the part in crack in 3 d data field.Therefore, analyze each Zhang Yanxin CT sequence image to carry out crack identification and probably occur deviation even mistake isolatedly.

Described crack and hole morphological feature difference in 3 d data field is comparatively large, learns by observing rock core CT sequence image meanwhile, and for real noise, then it does not have continuity between the rock core CT sequence image of different two or more than two; For crack, although be isolated impact point from two dimensional image, between the rock core CT sequence image of different two or more than two, very large probability is had to be interconnection.Therefore, the present invention carries out rock core crack identification according to the difference of crack, hole and random noise morphological feature in 3 d data field.Namely use CT machine to scan pending rock core, obtain rock core CT sequence image and store in a computer, carrying out three-dimensional image reconstruction by marching cube (Marching Cubes) algorithm.One or more than one target is obtained in the 3 d data field that three-dimensional image reconstruction is formed.Each target in described 3 d data field contains the information of the crack of rock core inside, hole and random noise.

The present invention is based on the rock core crack identification method of three-dimensional image information process, it is characterized in that comprising following concrete operation step successively:

Step 1: use CT machine to scan pending rock core, obtains rock core CT sequence image and stores in a computer;

Step 2: read in the rock core CT sequence image that step 1 scanned with computer software, and utilize marching cubes algorithm to carry out three-dimensional image reconstruction to rock core CT sequence image, in the 3 d data field that three-dimensional image reconstruction is formed, obtain one or more than one target;

Step 3: each target that in scanning step 2, three-dimensional image reconstruction obtains, obtains the grown form parameter of each target, comprises target volume, target surface area; Calculate the value of the 3D shape factor F of each target, the relational expression of F and described target volume and target surface area is following formula:

$F=\frac{36\mathrm{\π}{V^2}_p}{{S^3}_p}---\left(1\right)$

Wherein, V _pfor target volume, S _pfor target surface area, calculate the value of the F of each target according to formula (1); Again F value is screened, if meet F < 0.05, illustrate that this target is planar distribution, retain this area target, otherwise the non-area target of this target is described, then delete this non-area target;

Step 4: each area target retained in scanning step 3, obtains the minimal circumscribed sphere radius R of each area target _minwith equivalent spheres radius R _eq, calculate the ratio of minimal circumscribed sphere radius and equivalent spheres radius according to ratio screen again, if meet illustrate that this area target possesses ductility simultaneously, then retain the area target that this possesses ductility, otherwise illustrate that this area target does not possess ductility, then delete the area target that this does not possess ductility;

Step 5: each reservation in calculation procedure 4 possesses the length of the approximate minimum circumscribed rectangular body of the area target of ductility, compares between any two, obtain the longest edge D among them to length and width and height _maxmost minor face D _min, according to ratio screen, if meet illustrate that the malleable area target of this tool is crack, it can be used as the crack identified to retain, if the malleable area target of this tool does not meet described condition, illustrate that the malleable area target of this tool is not crack, then delete this area target.

Marching cube (Marching Cubes) algorithm described in step 2 of the present invention is the classic algorithm in iso-surface patch algorithm, and it is also referred to as isosurface extraction (Isosurface Extraction).Be applied in the present invention, be actually and each width figure in rock core CT sequence image is merged together, regard a 3 d data field as, therefrom the voxel with same grayscale value is extracted, couple together with topological form, form contour surface.The method extracting contour surface is each voxel processed one by one in above-mentioned 3 d data field, and the structural form of the inner contour surface of this voxel is decided according to the gray-scale value of each each angle point of voxel, the contour surface of all voxel inside is combined the most at last.

Described voxel generally has two kinds of definition: a kind of is similar with the definition of pixel in common two dimensional image, directly using the sampled point in 3 d data field as voxel; Another kind of then be eight adjacent sampled points in 3 d data field, i.e. summit on the square of unit length size, the region comprised is defined as voxel, and what the present invention adopted is the second definition.

Described marching cube (Marching Cubes) algorithm, its reference is as follows:

Qian Feng, horse is beautiful, Yang Shengqi, Wan Wanggen. the research of marching cubes algorithm and improvement [J] computer engineering and application, 2010, (34).

Neilson G M.On marching cubes[J].IEEE Transactions on Visualization andComputer Graphics，2003，9(3)：283-297.

Miao Bin and, Deng Yuanmu, Huang Feizeng. based on corresponding point matching faultage image three-dimensional interpolation method [J] Chinese medicine physics magazine, 2000,17 (1): 14-16.

Ginger Yun. based on the image three-dimensional visualization algorithm research [D] of iso-surface patch. Chengdu: University of Electronic Science and Technology, 2007.

The minimal circumscribed sphere of target described in step 4 of the present invention makes all voxels of this target all be positioned at the minimum ball of the radius of ball, and the radius of this ball is exactly the minimal circumscribed sphere radius R of described target _min; Described equivalent spheres radius R _eqit is the radius of the ball equal with described target volume.

Calculate the method that the area target possessing ductility is similar to the length of minimum circumscribed rectangular body described in step 5 of the present invention to comprise the following steps:

(1) height of the minimum circumscribed rectangular body of area target of ductility is possessed described in calculating

Adopt least square method to find the fit Plane of described target, its ultimate principle makes the quadratic sum of deviation minimum, is exactly specifically find a plane, makes every bit in described target minimum to the distance sum of plane; The plane in three-dimensional space formula adopting cartesian coordinate system to determine is as follows:

Ax+By+Cz+D＝0 (2)

Wherein (x, y, z) is a bit in plane, and A, B, C are the coefficient of x, y, z, and D is constant factor.

Any point (t in three dimensions ₁, t ₂, t ₃) as follows to the range formula of plane Ax+By+Cz+D=0:

$d=\frac{|{\mathrm{At}}_1+{\mathrm{Bt}}_2+{\mathrm{Ct}}_3+D|}{\sqrt{A^2+B^2+C^2}}---\left(3\right)$

The error formula calculating data point and plane in three dimensions is as follows:

e＝Ax+By+Cz+D (4)

According to least square ratio juris, the formula of plane least square fitting is as follows:

$\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{e_i}^2=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{Ax}}_i+{\mathrm{By}}_i+{\mathrm{Cz}}_i+D)}^2---\left(5\right)$

Wherein (x _i, y _i, z _i) for needing a bit in fit Plane, e _irepresent the error of this data point and plane.

For formula (2), if D=1, obtain plane Ax+By+Cz=-1, according to formula (5), if point (x _i, y _i, z _i) in the plane, then Ax _i+ By _i+ Cz _i=-1, e _i=Ax _i+ By _i+ Cz _i+ 1, plane least square fitting, namely $\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{e}_i^2=\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{Ax}}_i+{\mathrm{By}}_i+{\mathrm{Cz}}_i+1)}^2$ Minimum, respectively partial differential is got to A, B, C:

$\begin{array}{c}\frac{\&PartialD;}{\&PartialD;A}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{Ax}}_i+{\mathrm{By}}_i+{\mathrm{Cz}}_i+1)}^2=2\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}({\mathrm{Ax}}_i+{\mathrm{By}}_i+{\mathrm{Cz}}_i+1)x_i=2\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}({\mathrm{Ax}}_i^2+{\mathrm{By}}_i{x}_i+{\mathrm{Cz}}_i{x}_i+x_i)=0\\ \frac{\&PartialD;}{\&PartialD;B}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{Ax}}_i+{\mathrm{By}}_i+{\mathrm{Cz}}_i+1)}^2]=2\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}({\mathrm{Ax}}_i+{\mathrm{By}}_i+{\mathrm{Cz}}_i+1)y_i=2\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}({\mathrm{Ax}}_i{y}_i+{\mathrm{By}}_i^2+{\mathrm{Cz}}_i{y}_i+y_i)=0\\ \frac{\&PartialD;}{\&PartialD;B}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{({\mathrm{Ax}}_i+{\mathrm{By}}_i+{\mathrm{Cz}}_i+1)}^2=2\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}({\mathrm{Ax}}_i+{\mathrm{By}}_i+{\mathrm{Cz}}_i+1)Z_i=2\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}({\mathrm{Ax}}_i{z}_i+{\mathrm{By}}_i{z}_i+{\mathrm{Cz}}_i^1+z_i)=0\end{array}---\left(6\right)$

System of equations is obtained after arrangement:

$\left\{\begin{array}{c}A\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i^2+B\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i{x}_i+C\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{z}_i{x}_i+\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i=0\\ A\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i{y}_i+B\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i^2+C\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{z}_i{y}_i+\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i=0\\ A\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i{z}_i+B\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i{z}_i+C\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{z}_i^2+\underset{i=0}{\overset{n}{\mathrm{\&Sigma;}}}{z}_i=0\end{array}\right.---\left(7\right)$

Solve an equation:

$\left[\begin{array}{c}A\\ B\\ C\end{array}\right]{\left[\begin{array}{ccc}\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i^2& \underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i{y}_i& \underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i{y}_i\\ \underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i{y}_i& \underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i^2& \underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i{z}_i\\ \underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i{z}_i& \underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i{z}_i& \underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{z}_i^2\end{array}\right]}^{-1}\left[\begin{array}{c}-\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{x}_i\\ -\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{y}_i\\ -\underset{i=1}{\overset{n}{\mathrm{\&Sigma;}}}{z}_i\end{array}\right]---\left(8\right)$

According to (7) formula, by the computing to point set coordinate, the parameters of fit Plane can be obtained, and described target is formed by a series of voxel and three-dimensional point set, and the distribution in planar, therefore, the method for Selection utilization least square plane matching asks its fit Plane.

According to formula (2), if the equation of fit Plane is Ax+By+Cz+D=0, some P _i(x _i, y _i, z _i) be any point in described target, then its subpoint P in fit Plane _jcoordinate be P _j(x _i-Ak, y _i-Bk, z _i-Ck), wherein can invocation point P according to formula (3) _idistance to plane is travel through each point in described target, each some projection point set PointArray in the plane can be obtained, and each point is to the distance set M of plane.Each numerical value relatively in set M, judges k < 0 or k>=0, obtains being positioned at the ultimate range d of point to plane of plane both sides respectively in above-mentioned target _max1, d _max2, calculating above-mentioned target, to be similar to the formula of the height H of minimum circumscribed rectangular body as follows:

H＝d _max1+d _max2(9)

(2) length of the minimum circumscribed rectangular body of area target of ductility is possessed described in calculating

Any two points p is asked in known three dimensions ₁(a ₁, b ₁, c ₁), p ₂(a ₂, b ₂, c ₂) the formula of spacing as follows:

$f=\sqrt{{(a_1-a_2)}^2+{(b_1-b_2)}^2+{(c_1-c_2)}^2}---\left(10\right)$

To project point set PointArray as a convex closure, and use rotates the algorithm that gets stuck and obtains 2 maximum P of distance on convex closure ₁(x ₁, y ₁, z ₁) and P ₂(x ₂, y ₂, z ₂), the length that described target is similar to minimum circumscribed rectangular body is P ₁to P ₂distance, can obtaining described target according to (10) formula, to be similar to the formula of the length L of minimum circumscribed rectangular body as follows:

$L=\sqrt{{(x_1-x_2)}^2+{(y_1-y_2)}^2+{(z_1-z_2)}^2}---\left(11\right)$

(3) width of the minimum circumscribed rectangular body of area target of ductility is possessed described in calculating

Ask vector to vector the formula of projection as follows:

$n=\frac{|\stackrel{\&RightArrow;}{C}\&times;\stackrel{\&RightArrow;}{D}|}{\left|\stackrel{\&RightArrow;}{D}\right|}---\left(12\right)$

Set up an office P _k(x _k, y _k, z _k) be any point on convex closure, vector vector can invocation point P by formula (12) _kto P ₁, P ₂the distance d of the straight line determined _kformula as follows:

$d_k=\frac{|\stackrel{\&RightArrow;}{A}\&times;\stackrel{\&RightArrow;}{B}|}{\left|\stackrel{\&RightArrow;}{B}\right|}---\left(13\right)$

Each point on traversal convex closure, can obtain each point on convex closure and arrive P ₁and P ₂the distance set of the straight line determined.Each numerical value relatively in distance set, judges value be less than 0 or be more than or equal to 0, obtain respectively convex closure is positioned at by P ₁and P ₂the point of the straight line both sides determined is to the ultimate range d of this straight line _max3and d _max4, the width that above-mentioned target is similar to minimum circumscribed rectangular body is two ultimate range sums, and calculating above-mentioned target, to be similar to the formula of the width W of minimum circumscribed rectangular body as follows:

W＝d _max3+d _max4(14)

Through said process, finally calculate described target and be similar to the height of minimal circumscribed sphere, length and wide.

The present invention has following characteristics and useful technique effect:

Rock core crack identification method of the present invention uses CT machine to scan pending rock core, obtains rock core CT sequence image.Carry out three-dimensional image reconstruction by the rock core CT sequence chart obtained scanning, and utilize the difference that hole, random noise and crack in 3 d data field, morphological feature exist, rock core crack can be identified accurately and rapidly; Solve in the past isolated each Zhang Yanxin CT sequence image of analysis, but have ignored the error even problem of mistake that rock core crack communication information between the layers causes.

Accompanying drawing explanation

The rock core CT sequence image that Fig. 1 the present invention is pending; Wherein, (a) is 4 pictures in rock core CT sequence image A, and (b) is 4 pictures in rock core CT sequence image B;

Fig. 2 is the result in Fig. 1 under rock core CT sequence image different points of view after three-dimensional image reconstruction; Wherein, (a) is the three-dimensional image reconstruction result of rock core CT sequence image A, and (b) is the three-dimensional image reconstruction result of rock core CT sequence image B;

Fig. 3 is the result in the three-dimensional rock core crack extracted in rock core CT sequence image in Fig. 2; Wherein, (a) is the extraction result of rock core CT sequence image A, and (b) is the extraction result of rock core CT sequence image B; From then on figure intuitively can clearly see the crack existed in rock core CT sequence image, namely illustrates that rock core crack identification effect of the present invention is accurately.

Embodiment

Below in conjunction with embodiment and result figure, the present invention is described in further detail, but protection scope of the present invention is not limited to the following stated; All " rock core crack identification software based on three-dimensional image information process " mentioned is all by seminar of the present invention independent research below, and its use operation steps is the operation steps for rock core crack identification method.

Embodiment 1

Step 1: use CT machine to scan pending rock core, obtains 400 rock core CT sequence images and stores in a computer, as (a) figure in Fig. 1, observes picture and finds, position to the right, figure center has the target of long and narrow linear;

Step 2: the reading picture function using core three-dimension image crack identification software, read in the rock core CT sequence image that step 1 has scanned, this software carries out three-dimensional image reconstruction according to marching cube (Marching Cubes) algorithm to this CT sequence chart, obtains the 3 d data field be made up of above-mentioned rock core CT sequence image.

Step 3: using the targeted scans function of core three-dimension image crack identification software to carry out stroke scanning to rebuilding the 3 d data field obtained in step 2, often obtaining a target and just being preserved.After having scanned, obtain 53 mutual disjunct targets altogether, as (a) figure in Fig. 2, from then in figure, can be observed the target of long and narrow planar, and granular Small object;

Step: 4: the calculating form factor function using core three-dimension image crack identification software, in calculation procedure 3, the grown form parameter of 53 targets that targeted scans obtains, comprises target volume, target surface area; The value of the 3D shape factor F of each target above-mentioned is calculated by formula (1); After calculating completes, find to only have a goal satisfaction F < 0.05, be area target, this area target retains by this software, and all the other targets is all deleted;

Step 5: the calculating radius of a ball function using core three-dimension image crack identification software, the minimal circumscribed sphere radius R of the area target retained in calculation procedure 4 _minwith equivalent spheres radius R _eq, and calculate the ratio of minimal circumscribed sphere radius and equivalent spheres radius find this goal satisfaction for possessing malleable area target, the area target that this is possessed ductility by this software retains;

Step 6: the digital simulation plane function using core three-dimension image crack identification software, by the fit Plane possessing the area target of ductility retained in formula (2) and formula (8) calculation procedure 5.

Step 7: the calculating minimum circumscribed rectangular height function using core three-dimension image crack identification software, by the approximate minimum circumscribed rectangular height possessing the area target of ductility retained in formula (1) and formula (9) calculation procedure 5.

Step 8: the calculating minimum circumscribed rectangular body length function using core three-dimension image crack identification software, by the approximate minimum circumscribed rectangular body length possessing the area target of ductility retained in formula (11) calculation procedure 5.

Step 9: the calculating minimum circumscribed rectangular body width function using core three-dimension image crack identification software, by the approximate minimum external wide side's height possessing the area target of ductility retained in formula (13) and formula (14) calculation procedure 5.

Step 10: the calculating minimum circumscribed rectangular body length of side ratio function using core three-dimension image crack identification software, to step 7, step 8, obtains length, the width of approximate minimum circumscribed rectangular body respectively and compares highly between any two, obtaining the longest edge D among them in step 9 _maxmost minor face D _min, calculate ratio, find that the area target possessing ductility retained in step 5 meets this software it can be used as the crack identified to retain, as shown in (a) figure in Fig. 3, through identifying that this rock core co-exists in a crack.

Have recorded crack 1 target that extracts through above-mentioned steps and wherein three hole targets---hole 1 in following table 1, the three-dimensional configuration parameter comparison of hole 2 and hole 3, by formula (1), with calculating acquired results as shown in table 1 below.

Table 1 crack 1 and hole 1, the three-dimensional configuration parameter comparison of hole 2 and hole 3

Observation table 1 finds, crack 1 wherein meets step 3, and three screening conditions of step 4 and step 5, are successfully identified; And hole 1, hole 2 and 3, hole do not meet step 3, three screening conditions of step 4 and step 5, finally deleted fall.

Embodiment 2

Step 1: use CT machine to scan pending rock core, obtain 205 rock core CT sequence images and store in a computer, as (b) figure in Fig. 1, observation picture finds, except the target of long and narrow linear, also have sub-circular and granular target in figure, their gray level closely;

Step 2: the reading picture function using core three-dimension image crack identification software, read in the rock core CT sequence image that step 1 has scanned, this software carries out three-dimensional image reconstruction according to marching cube (Marching Cubes) algorithm to this CT sequence chart, obtains the 3 d data field be made up of above-mentioned rock core CT sequence image.

Step 3: using the targeted scans function of core three-dimension image crack identification software to carry out stroke scanning to rebuilding the 3 d data field obtained in step 2, often obtaining a target and just being preserved.After having scanned, obtain 340 mutual disjunct targets altogether, as (b) figure in Fig. 2, learn from figure, have the target of long and narrow planar, and graininess, approximate sphericity Small object;

Step: 4: the calculating form factor function using core three-dimension image crack identification software, in calculation procedure 3, the grown form parameter of 340 targets that targeted scans obtains, comprises target volume, target surface area; The value of the 3D shape factor F of each target above-mentioned is calculated by formula (1); After calculating completes, finding that there is 3 goal satisfaction F < 0.05, is area target, and these 3 area targets retain by this software, and all the other targets are all deleted;

Step 5: the calculating radius of a ball function using core three-dimension image crack identification software, the minimal circumscribed sphere radius R of 3 area targets retained in calculation procedure 4 _minwith equivalent spheres radius R _eq, and calculate the ratio of minimal circumscribed sphere radius and equivalent spheres radius find that there is 2 goal satisfactions for possessing malleable area target, this software possesses ductility area target by these 2 retains, and deletes another one target;

Step 6: the digital simulation plane function using core three-dimension image crack identification software, possesses the fit Plane of the area target of ductility by 2 of reservation in formula (2) and formula (8) calculation procedure 5.

Step 7: the calculating minimum circumscribed rectangular height function using core three-dimension image crack identification software, possesses the approximate minimum circumscribed rectangular height of the area target of ductility by 2 of reservation in formula (1) and formula (9) calculation procedure 5.

Step 8: the calculating minimum circumscribed rectangular body length function using core three-dimension image crack identification software, possesses the approximate minimum circumscribed rectangular body length of the area target of ductility by retain in formula (11) calculation procedure 52.

Step 9: the calculating minimum circumscribed rectangular body width function using core three-dimension image crack identification software, possesses approximate minimum external wide side's height of the area target of ductility by 2 of reservation in formula (13) and formula (14) calculation procedure 5.

Step 10: the calculating minimum circumscribed rectangular body length of side ratio function using core three-dimension image crack identification software, to step 7, step 8, obtains length, the width of approximate minimum circumscribed rectangular body respectively and compares highly between any two, obtaining the longest edge D among them in step 9 _maxmost minor face D _min, calculate ratio, find that 2 area targets possessing ductility retained in step 5 only have 1 to meet this software it can be used as the crack identified to retain, and deletes another one target, as shown in (b) figure in Fig. 3, through identifying that this rock core co-exists in a crack.

Have recorded crack 2 target that extracts through above-mentioned steps and wherein three hole targets---hole 4 in following table 2, the three-dimensional configuration parameter comparison of hole 5 and hole 6, by formula (1), with calculating acquired results as shown in table 2 below.

The three-dimensional configuration parameter comparison of table 2 crack 2 and hole 4, hole 5 and hole 6

	Crack 2	Hole 4	Hole 5	Hole 6
					Form factor	0.00156837	0.205297	0.252467	0.894843
The ratio of two radius	5.528	1.3812	2.6641	1.3476
					The ratio on both sides	4.303	2.18996	1.16752	2.47214

Observation table 2 finds, crack 2 wherein meets step 3, and three screening conditions of step 4 and step 5, are successfully identified; And hole 4, hole 5 and 6, hole do not meet step 3, three screening conditions of step 4 and step 5, finally deleted fall.

Claims (1)

1., based on the recognition methods in the rock core crack of three-dimensional image information process, it is characterized in that comprising following concrete operation step successively:

Step 1: use CT machine to scan pending rock core, obtains rock core CT sequence image and stores in a computer;

Step 2: use core three-dimension image crack identification software reads in the rock core CT sequence image that step 1 has scanned, and utilize marching cubes algorithm to carry out three-dimensional image reconstruction to rock core CT sequence image, in the 3 d data field that three-dimensional image reconstruction is formed, obtain one or more than one target;

Step 3: each target that in scanning step 2, three-dimensional image reconstruction obtains, obtains the grown form parameter of each target, comprises target volume, target surface area; Calculate the value of the 3D shape factor F of each target, the relational expression of F and described target volume and target surface area is following formula:

$F=\frac{36\mathrm{\&pi;}{V^2}_p}{{S^3}_p}---\left(1\right)$

Wherein, V _pfor target volume, S _pfor target surface area, calculate the value of the F of each target according to formula (1); Again F value is screened, if meet F < 0.05, illustrate that this target is planar distribution, retain this area target, otherwise the non-area target of this target is described, then delete this non-area target;

Step 4: each area target retained in scanning step 3, obtains the minimal circumscribed sphere radius R of each area target _minwith equivalent spheres radius R _eq, calculate the ratio of minimal circumscribed sphere radius and equivalent spheres radius according to ratio screen again, if meet illustrate that this area target possesses ductility simultaneously, then retain the area target that this possesses ductility, otherwise illustrate that this area target does not possess ductility, then delete the area target that this does not possess ductility;

Step 5: each reservation in calculation procedure 4 possesses the length of the approximate minimum circumscribed rectangular body of the area target of ductility, compares between any two, obtain the longest edge D among them to length and width and height _maxmost minor face D _min, according to ratio screen, if meet illustrate that the malleable area target of this tool is crack, it can be used as the crack identified to retain, if the malleable area target of this tool does not meet described condition, illustrate that the malleable area target of this tool is not crack, then delete this area target.