CN101533518A - Method for re-establishing surface of three dimensional target object by unparallel dislocation image sequence - Google Patents

Abstract

The invention provides a method for re-establishing a surface of a three dimensional target object by an unparallel dislocation image sequence, and belongs to the technical field of re-establishment and visualization of a three dimensional image. The method comprises that: the outline of the target object is extracted by cutting the dislocation image; the outline of the target object is subjected to distance conversion, extraction of disc aggregate and removal of redundancy so as to find out a maximum disc aggregate in the outline of the target object; and the disc aggregate is selected out from adjacent dislocation images for interpolation, a vector value of each point, which points to the target object in the image, in the dislocation image is calculated; a distance value of space points between the dislocation images needed for interpolation is determined; and the outline based surface is drawn and displayed. The method realizes the direct re-establishment of the surface of the three dimensional target object by the unparallel dislocation image sequence with high re-establishment accuracy and speed, and can also achieve quite good effect on the surface re-establishment of a three dimensional target object with a complex geometric shape.

Description

A kind of method by non-parallel tomographic sequence reconstruction of three-dimensional destination object surface

(1) technical field

The present invention relates to a kind of image reconstruction and visualization method, relate in particular to a kind of method by non-parallel tomographic sequence reconstruction of three-dimensional destination object surface.

(2) background technology

Some application scenario in three-dimensional image reconstruction and visualization technique field, what people obtained via data acquisition is a series of nonparallel tomographic sequence of objective object, and observe one of objective object most common form is that the surface shows, this just requires us must be by the surface of non-parallel tomographic sequence reconstruction of three-dimensional destination object.Here, the resurfacing based on profile is the method that is most widely used.After normally at first recombinating nonparallel tomographic sequence and reset, the method for surface reconstruction based on profile in the past obtains the volume data space of a rule, in this volume data space, some parallel faultage images are carried out the manual or semi-automatic profile that obtains destination object on it of cutting apart, and then interpolation obtains the profile of the destination object on other faultage image between these parallel fault images, carries out at last showing based on the iso-surface patch of profile.Reorganization rearrangement interpolation processing calculated amount in the previous methods is big, processing procedure engineering practicability comparatively consuming time, method is relatively poor.And previous methods needs interpolation calculation at least twice, has reduced the reconstruction precision of method.

It is investigated, a spot of researcher is also arranged so far to studying by non-parallel tomographic sequence reconstruction of three-dimensional destination object cosmetic issue, in article " Fast surface and volume estimation fromnon-parallel cross-sections; for freehand 3-D ultrasound; 1998 " and " Surfaceinterpolation from sparse cross-sections using region correspondence; 1999 ", realized surface as people such as G.M.Treece, but its process of reconstruction calculated amount is big by non-parallel tomographic sequence reconstruction of three-dimensional destination object based on shape interpolation method, length consuming time.People such as A.L.Bogush adopt cubic spline interpolation method head it off in article " 3D Object Reconstruction from Non-parallelCross-sections; 2004 ", their method can obtain effect preferably to the comparatively simple objective object of geometric shape, but, because their method space geometry relation can not correct calculation profile between very undesirable for the objective object effect that geometric shape is comparatively complicated.

(3) summary of the invention

Length consuming time, precision at the described method for reconstructing of background technology is low, poor practicability, be difficult to adapt to the shortcoming of the objective object of geometric shape complexity, the present invention proposes a kind of method by the direct reconstruction of three-dimensional destination object of original non-parallel tomographic sequence surface, it is directly cut apart by the original faultage image of part and obtains the objective contours of objects, and interpolation obtains showing needed whole profile based on the iso-surface patch of profile again.

A kind of method by non-parallel tomographic sequence reconstruction of three-dimensional destination object surface, step is as follows:

S1) the part faultage image is cut apart, extracted the destination object profile

Earlier the part faultage image is cut apart, extracted the profile of destination object in these faultage images; In the tomographic sequence that extracts the destination object profile, every two adjacent frame faultage images are all carried out following steps S2 successively)-S6);

S2) faultage image is carried out range conversion, obtain range image

The two adjacent frame faultage images that obtain the destination object profile are converted into bianry image; Carry out the Chamfer range conversion then, choose 3 * 3 template, show as Fig. 2; After range conversion finishes, make the distance value of profile external point get original negative, the distance value of profile internal point is constant, to distinguish exterior point in the profile, the image that obtains is the range image of faultage image, each point of destination object profile outside wherein, its pixel value are represented this negative value of some distances recently to the destination object profile; Each inner point, its pixel value are represented this distance value of any recently to the destination object profile; The pixel value of putting on the destination object profile is a null value;

S3) extract the very big disk collection in the destination object profile in the range image, and remove redundant

Extract the disk collection: a series of tangent and fully at the disk collection of destination object profile inside, each disk wherein can not covered by all the other all disks institutes fully with the destination object profile by calculating; These disk collection have been represented the skeleton of destination object and have been moved towards trend;

At first find out the maximum point of interior all pixel values of objective contour in the range image, be the central point of disk: the pixel value to each some p in the destination object profile all will compare judgement with the pixel value of its 8 adjoint points, if satisfy pixel (p)〉pixel (q)-t, then this point is maximum point; Its mid point q is that a p reaches 8 consecutive point on the diagonal up and down, and t is the relative distance value of some q to a p; For a p four adjoint points on the four direction up and down, t=3; For four adjoint points on the p diagonal, t=4; Each point in the range image all will be judged, up to finding out all maximum points; The all corresponding disk of each maximum point, and disc centre is positioned at the maximum point place;

Calculate the radius of the corresponding disk of each maximum point then, each maximum point is carried out following same treatment respectively; The image of at first adjusting the distance carries out initialization, and except that when the maximum point of pre-treatment, all the other some pixel values all are initialized as 0; Then each point in the image is carried out twice scanning similar to range conversion, get 3 * 3 template, as Fig. 2, the pixel value that obtains after each scanning is determined by following formula:

$v_{i,j}^{\mathrm{num}}=\underset{(i,j)\∈\mathrm{mask}}{\max }(v_{i+k,j+l}^{\mathrm{num}-1}-t(k,l))$

Wherein

For on the image (i, j) position is in the pixel value after the scanning the num time, (k l) be the position coordinates at relative its center (0,0) of template, and t (k is a template in that (k l) locates relative distance value with the center l); If point (k, when l) being positioned at four adjoint point places on the four direction up and down of center, t (k, l)=3; If (k, when l) being positioned at four adjoint point places on the diagonal of center, t (k, l)=4; Through the image after twice scanning, have only wherein that the pixel value of point is a nonzero value near the maximum point, and the pixel value of these points is less than the pixel value of maximum point, the pixel value of other point is a null value; Find out minimum non-zero value w wherein, then the radius of this disk is r=v-w, and wherein v is the pixel value at maximum point place, and r is the radius of the corresponding disk of this maximum point; Repeating aforesaid operations is all calculated up to the disc radius of each maximum point correspondence;

The disk collection is gone redundancy: make that M is original very big disk collection, N is redundant very big disk collection; Take out the disk fm of radius maximum among the M, this disk is put into N; Then each remaining among M disk fn is compared with disk fm, get fl and be the distance between two circle disk centers, get fd=radius (fm)-0.5 * radius (fn), if fl ²＜fd ², then fn is removed from M; Repeating aforesaid operations then, is empty up to M; The N that finally obtains is the disk collection after the redundancy;

S4) in faultage image, select to be used for the disk collection of interpolation

During pixel value between the interpolation two adjacent faultage images, wherein interior each the very big disk sm of destination object profile carries out projection (for non-parallel faultage image to another width of cloth image B successively among the piece image A, average normal direction along two adjacent faultage images is carried out projection), as Fig. 4, judge the projection of disk sm and be projected any disk sn that comprises in the destination object profile in the image whether lap is arranged, if lap is arranged, then sm is put into set C ₁, sn puts into set C ₂, if sn is at C ₂In exist, just need not put inward; Repeat said process then, all disks in handling A; Last C ₁With C ₂Be respectively the used disk collection of interpolation in the two adjacent faultage images;

S5) every bit points to the vector value of destination object in the image in the computed tomography image

According to step S4) be used for the disk collection of interpolation in the two adjacent faultage images that obtain, to each the some a in the two adjacent faultage images _iCarry out following calculating,

$\stackrel{\→}{c}=\frac{1}{\underset{d\∈c}{\mathrm{\Σ}}\frac{r_d}{{\left|a_i{\stackrel{\→}{k}}_d\right|}^2}}\underset{d\∈c}{\mathrm{\Σ}}{r}_d\frac{a_i{\stackrel{\→}{k}}_d}{{\left|a_i{\stackrel{\→}{k}}_d\right|}^2}$

Wherein, c represents to be used in the faultage image disk collection of interpolation, k _dBe the center of disk d, r _dBe the radius of disk d, a _iBe the arbitrfary point in the faultage image,

Be a vector value, expression faultage image mid point a _iDirection and distance to destination object; The result promptly obtains and corresponding two width of cloth new images of two adjacent faultage images, wherein the pixel value of every bit is a vector value, has represented the direction and the distance of this some destination object in the image;

S6) determine the distance value of spatial point between the interpolation faultage image

In the range image of two adjacent faultage images wherein a width of cloth A as starting point, traversal is each some a wherein, spend this o'clock and do straight line along the average normal direction of two faultage images and hand over second width of cloth faultage image B, connect a point two vectors corresponding in two faultage images then with b point place in 1 b Obtain new vector i, as shown in Figure 5; Cross the direction of a point along i then, do straight line and hand over faultage image B in a c, line segment ac is the corresponding interpolation line segment of a point;

According to a and 2 distance values of locating of c the point g on the line segment ac is carried out the cosine interpolation, formula is as follows:

$d_g=(1-\frac{\cos \frac{l_1}{{\mathrm{<figure-callout\; id="1"\; label="l"\; filenames="A200910020718E00141.png,A200910020718E00151.png"\; state="\{\{state\}\}">l</figure-callout>}}_1+{\mathrm{<figure-callout\; id="2"\; label="l"\; filenames="A200910020718E00142.png,A200910020718E00151.png"\; state="\{\{state\}\}">l</figure-callout>}}_2}\mathrm{\&pi;}+1}{2}){\mathrm{<figure-callout\; id="1"\; label="d"\; filenames="A200910020718E00141.png,A200910020718E00151.png"\; state="\{\{state\}\}">d</figure-callout>}}_1+\frac{\cos \frac{l_1}{{\mathrm{<figure-callout\; id="1"\; label="l"\; filenames="A200910020718E00141.png,A200910020718E00151.png"\; state="\{\{state\}\}">l</figure-callout>}}_1+{\mathrm{<figure-callout\; id="2"\; label="l"\; filenames="A200910020718E00142.png,A200910020718E00151.png"\; state="\{\{state\}\}">l</figure-callout>}}_2}\mathrm{\&pi;}+1}{2}{\mathrm{<figure-callout\; id="2"\; label="d"\; filenames="A200910020718E00142.png,A200910020718E00151.png"\; state="\{\{state\}\}">d</figure-callout>}}_2$

L wherein ₁And l ₂Be respectively g section along the line ac direction to the distance of some a with some c, d ₁, d ₂Be respectively the distance value of a point and c point place in two range images, d _gBe the distance value at the some g place that obtains of interpolation, as Fig. 6; If when being straight line and faultage image B and not had intersection point by a point, then give up a point, to a bit repeating calculating like this down, after the every bit among the faultage image A all traveled through and finishes, the distance value of spatial point had promptly been come out by correct interpolation between the two adjacent faultage images again; Determine that according to the positive and negative of distance value this point is in profile inside or outside, and zero point is the point on the profile;

All dispose as if every adjacent two frame faultage images in the tomographic sequence that has extracted the destination object profile, then execution in step S7), otherwise continue execution in step S2)-S6);

S7) iso-surface patch based on profile shows

All spatial point distance values are extracted zero point, and carry out showing that based on the iso-surface patch of profile the entire process flow process finishes.

Step S2 wherein) the Chamfer distance in is a transform method classical in the range conversion, and this title is taken from term in the woodworking, i.e. chamfering, and these class methods are successively by these characteristics of two-pass scan of opposite both direction during in order to the expression computed range.

A kind of freedom-arm, three-D ultrasonic image-forming system that uses method for reconstructing of the present invention comprises microcomputer, ultrasound scanner, locating device, image capture device, it is characterized in that ultrasound scanner is directly connected to microcomputer by digital output port; Perhaps be connected to image capture device by analog interface, image capture device is connected to microcomputer by pci interface or USB interface; The moveable part of locating device is fixed on above the probe of ultrasound scanner, and locating device is connected to microcomputer by wired or wireless mode; Microcomputer obtains the image information of ultrasound scanner by digital input port or image capture device, obtains the positional information of ultrasonic scanning instrument probe by locating device.

The course of work of imaging system of the present invention is as follows:

101. collect the tomographic sequence of human organ interested

The scanning probe that the operator moves freely ultrasound scanner obtains the tomographic sequence of human organ interested;

102. the part faultage image is carried out segmented extraction destination object profile

At first cut apart, extract the profile of destination object in the part faultage image; Here, extract profile based on self-adaptation movable contour model (being the T-Snake model);

In the tomographic sequence that extracts the destination object profile, every two adjacent frame faultage images are all carried out following steps 103-107 according to steps in sequence of the present invention;

103. faultage image is carried out range conversion, obtains range image

104. the very big disk collection in the extraction range image in the destination object profile, and remove redundant

Extract the disk collection: can calculate a series of tangent and fully at the disk collection of destination object profile inside, the former dish of wherein each can not covered by all the other all former dishes with the destination object profile; These disk collection have been represented the skeleton of destination object and have been moved towards trend;

105. in faultage image, select to be used for the disk collection of interpolation

106. every bit points to the vector value of destination object in the image in the computed tomography image

107. determine the distance value of spatial point between the interpolation faultage image

All dispose as if every adjacent two frame faultage images in the tomographic sequence that has extracted the destination object profile, then execution in step 108, otherwise continue to carry out 103-107;

108. the iso-surface patch based on profile shows

All spatial point distance values are extracted zero point, and carry out showing that based on the iso-surface patch of profile the entire process flow process finishes.

The present invention realized by the direct reconstruction of three-dimensional destination object of non-parallel tomographic sequence surface, rebuilds the precision height, and speed is fast, the objective object of complex geometry carried out resurfacing also have good effect.

(4) description of drawings

Fig. 1 is the process flow diagram of the inventive method.

Fig. 2 is 3 * 3 Chamfer range conversion mask template synoptic diagram, and it is divided into scan forward template and back to scan templates.

Fig. 3 is the synoptic diagram of a kind of very big disk collection finally obtained among the step S3.

How Fig. 4 selects to be used for the synoptic diagram of very big disk of interpolation among the step S4.Have only a disk among the figure in the adjacent faultage image in the destination object profile, r (sm) wherein, r (sn) is respectively disk sm, the radius of sn.Sl is the distance between latter two circle disk center of projection.D (sm), d (sn) are respectively the distance values at circle center line connecting and profile intersection point place on image A, B.N is the projecting direction of disk, i.e. the average normal direction of two faultage images.

Fig. 5 is the synoptic diagram that calculates the interpolation direction among the step S5.A, b are respectively the intersection point of straight line ab and two fault planes, and c (1) and c (2) are respectively the vectors at a point and b point place, and i is the line of c (1), c (2) two vector terminals.

Fig. 6 is the synoptic diagram that is carried out interpolation among the step S6 by the interpolation direction.Point g is the interpolation point, and d (1) and d (2) were respectively the g point along the straight line of interpolation direction and the intersection point of two faultage images, and l (1) and l (2) are respectively a g through the distance of interpolation direction place straight line to two planes.

(5) embodiment

The present invention will be further described below in conjunction with drawings and Examples, but be not limited thereto.

Embodiment:

S1) the part faultage image is cut apart, extracted the destination object profile

Earlier the part faultage image is cut apart, extracted the profile of destination object in these faultage images; In the tomographic sequence that extracts the destination object profile, every two adjacent frame faultage images are all carried out following steps S2 successively)-S6);

S2) faultage image is carried out range conversion, obtain range image

The two adjacent frame faultage images that obtain the destination object profile are converted into bianry image; Carry out the Chamfer range conversion then, choose 3 * 3 template, show as Fig. 2; After range conversion finishes, make the distance value of profile external point get original negative, the distance value of profile internal point is constant, to distinguish exterior point in the profile, the image that obtains is the range image of faultage image, each point of destination object profile outside wherein, its pixel value are represented this negative value of some distances recently to the destination object profile; Each inner point, its pixel value are represented this distance value of any recently to the destination object profile; The pixel value of putting on the destination object profile is a null value;

S3) extract the very big disk collection in the destination object profile in the range image, and remove redundant

Extract the disk collection: a series of tangent and fully at the disk collection of destination object profile inside, each disk wherein can not covered by all the other all disks institutes fully with the destination object profile by calculating; These disk collection have been represented the skeleton of destination object and have been moved towards trend;

At first find out the maximum point of interior all pixel values of objective contour in the range image, be the central point of disk: the pixel value to each some p in the destination object profile all will compare judgement with the pixel value of its 8 adjoint points, if satisfy pixel (p)〉pixel (q)-t, then this point is maximum point; Its mid point q is that a p reaches 8 consecutive point on the diagonal up and down, and t is the relative distance value of some q to a p; For a p four adjoint points on the four direction up and down, t=3; For four adjoint points on the p diagonal, t=4; Each point in the range image all will be judged, up to finding out all maximum points; The all corresponding disk of each maximum point, and disc centre is positioned at the maximum point place;

Calculate the radius of the corresponding disk of each maximum point then, each maximum point is carried out following same treatment respectively; The image of at first adjusting the distance carries out initialization, and except that when the maximum point of pre-treatment, all the other some pixel values all are initialized as 0; Then each point in the image is carried out twice scanning similar to range conversion, get 3 * 3 template, as Fig. 2, the pixel value that obtains after each scanning is determined by following formula:

$v_{i,j}^{\mathrm{num}}=\underset{(i,j)\&Element;\mathrm{mask}}{\max }(v_{i+k,j+l}^{\mathrm{num}-1}-t(k,l))$

Wherein

For on the image (i, j) position is in the pixel value after the scanning the num time, (k l) be the position coordinates at relative its center (0,0) of template, and t (k is a template in that (k l) locates relative distance value with the center l); If point (k, when l) being positioned at four adjoint point places on the four direction up and down of center, t (k, l)=3; If (k, when l) being positioned at four adjoint point places on the diagonal of center, t (k, l)=4; Through the image after twice scanning, have only wherein that the pixel value of point is a nonzero value near the maximum point, and the pixel value of these points is less than the pixel value of maximum point, the pixel value of other point is a null value; Find out minimum non-zero value w wherein, then the radius of this disk is r=v-w, and wherein v is the pixel value at maximum point place, and r is the radius of the corresponding disk of this maximum point; Repeating aforesaid operations is all calculated up to the disc radius of each maximum point correspondence;

The disk collection is gone redundancy: make that M is original very big disk collection, N is redundant very big disk collection; Take out the disk fm of radius maximum among the M, this disk is put into N; Then each remaining among M disk fn is compared with disk fm, get fl and be the distance between two circle disk centers, get fd=radius (fm)-0.5 * radius (fn), if fl ²＜fd ², then fn is removed from M; Repeating aforesaid operations then, is empty up to M; The N that finally obtains is the disk collection after the redundancy;

S4) in faultage image, select to be used for the disk collection of interpolation

During pixel value between the interpolation two adjacent faultage images, wherein interior each the very big disk sm of destination object profile carries out projection (for non-parallel faultage image to another width of cloth image B successively among the piece image A, average normal direction along two adjacent faultage images is carried out projection), as Fig. 4, judge the projection of disk sm and be projected any disk sn that comprises in the destination object profile in the image whether lap is arranged, if lap is arranged, then sm is put into set C ₁, sn puts into set C ₂, if sn is at C ₂In exist, just need not put inward; Repeat said process then, all disks in handling A; Last C ₁With C ₂Be respectively the used disk collection of interpolation in the two adjacent faultage images;

S5) every bit points to the vector value of destination object in the image in the computed tomography image

According to step S4) be used for the disk collection of interpolation in the two adjacent faultage images that obtain, to each the some a in the two adjacent faultage images _iCarry out following calculating,

$\stackrel{\&RightArrow;}{c}=\frac{1}{\underset{d\&Element;c}{\mathrm{\&Sigma;}}\frac{r_d}{{\left|a_i{\stackrel{\&RightArrow;}{k}}_d\right|}^2}}\underset{d\&Element;c}{\mathrm{\&Sigma;}}{r}_d\frac{a_i{\stackrel{\&RightArrow;}{k}}_d}{{\left|a_i{\stackrel{\&RightArrow;}{k}}_d\right|}^2}$

Wherein, c represents to be used in the faultage image disk collection of interpolation, k _dBe the center of disk d, r _dBe the radius of disk d, a _iBe the arbitrfary point in the faultage image, Be a vector value, expression faultage image mid point a _iDirection and distance to destination object; The result promptly obtains and corresponding two width of cloth new images of two adjacent faultage images, wherein the pixel value of every bit is a vector value, has represented the direction and the distance of this some destination object in the image;

S6) determine the distance value of spatial point between the interpolation faultage image

In the range image of two adjacent faultage images wherein a width of cloth A as starting point, traversal is each some a wherein, spend this o'clock and do straight line along the average normal direction of two faultage images and hand over second width of cloth faultage image B, connect a point two vectors corresponding in two faultage images then with b point place in 1 b Obtain a new vector i as shown in Figure 5; Cross the direction of a point along i then, do straight line and hand over faultage image B in a c, line segment ac is the corresponding interpolation line segment of a point;

According to a and 2 distance values of locating of c the point g on the line segment ac is carried out the cosine interpolation, formula is as follows:

$d_g=(1-\frac{\cos \frac{l_1}{{\mathrm{<figure-callout\; id="1"\; label="l"\; filenames="A200910020718E00141.png,A200910020718E00151.png"\; state="\{\{state\}\}">l</figure-callout>}}_1+{\mathrm{<figure-callout\; id="2"\; label="l"\; filenames="A200910020718E00142.png,A200910020718E00151.png"\; state="\{\{state\}\}">l</figure-callout>}}_2}\mathrm{\&pi;}+1}{2}){\mathrm{<figure-callout\; id="1"\; label="d"\; filenames="A200910020718E00141.png,A200910020718E00151.png"\; state="\{\{state\}\}">d</figure-callout>}}_1+\frac{\cos \frac{l_1}{{\mathrm{<figure-callout\; id="1"\; label="l"\; filenames="A200910020718E00141.png,A200910020718E00151.png"\; state="\{\{state\}\}">l</figure-callout>}}_1+{\mathrm{<figure-callout\; id="2"\; label="l"\; filenames="A200910020718E00142.png,A200910020718E00151.png"\; state="\{\{state\}\}">l</figure-callout>}}_2}\mathrm{\&pi;}+1}{2}{\mathrm{<figure-callout\; id="2"\; label="d"\; filenames="A200910020718E00142.png,A200910020718E00151.png"\; state="\{\{state\}\}">d</figure-callout>}}_2$

L wherein ₁And l ₂Be respectively g section along the line ac direction to the distance of some a with some c, d ₁, d ₂Be respectively the distance value of a point and c point place in two range images, d _gBe the distance value at the some g place that obtains of interpolation, as Fig. 6; If when being straight line and faultage image B and not had intersection point by a point, then give up a point, to a bit repeating calculating like this down, after the every bit among the faultage image A all traveled through and finishes, the distance value of spatial point had promptly been come out by correct interpolation between the two adjacent faultage images again; Determine that according to the positive and negative of distance value this point is in profile inside or outside, and zero point is the point on the profile;

All dispose as if every adjacent two frame faultage images in the tomographic sequence that has extracted the destination object profile, then execution in step S7), otherwise continue execution in step S2)-S6);

S7) iso-surface patch based on profile shows

All spatial point distance values are extracted zero point, and carry out showing that based on the iso-surface patch of profile the entire process flow process finishes.

Claims (1)

1, a kind of method by non-parallel tomographic sequence reconstruction of three-dimensional destination object surface, step is as follows:

S1) the part faultage image is cut apart, extracted the destination object profile

Earlier the part faultage image is cut apart, extracted the profile of destination object in these faultage images; In the tomographic sequence that extracts the destination object profile, every two adjacent frame faultage images are all carried out following steps S2 successively)-S6);

S2) faultage image is carried out range conversion, obtain range image

The two adjacent frame faultage images that obtain the destination object profile are converted into bianry image; Carry out the Chamfer range conversion then, choose 3 * 3 template, after range conversion finishes, make the distance value of profile external point get original negative, the distance value of profile internal point is constant, and to distinguish exterior point in the profile, the image that obtains is the range image of faultage image, each point of destination object profile outside wherein, its pixel value are represented this negative value of some distances recently to the destination object profile; Each inner point, its pixel value are represented this distance value of any recently to the destination object profile; The pixel value of putting on the destination object profile is a null value;

S3) extract the very big disk collection in the destination object profile in the range image, and remove redundant

Extract the disk collection: a series of tangent and fully at the disk collection of destination object profile inside, each disk wherein can not covered by all the other all disks institutes fully with the destination object profile by calculating; These disk collection have been represented the skeleton of destination object and have been moved towards trend;

At first find out the maximum point of interior all pixel values of objective contour in the range image, be the central point of disk: the pixel value to each some p in the destination object profile all will compare judgement with the pixel value of its 8 adjoint points, if satisfy pixel (p)〉pixel (q)-t, then this point is maximum point; Its mid point q is that a p reaches 8 consecutive point on the diagonal up and down, and t is the relative distance value of some q to a p; For a p four adjoint points on the four direction up and down, t=3; For four adjoint points on the p diagonal, t=4; Each point in the range image all will be judged, up to finding out all maximum points; The all corresponding disk of each maximum point, and disc centre is positioned at the maximum point place;

Calculate the radius of the corresponding disk of each maximum point then, each maximum point is carried out following same treatment respectively; The image of at first adjusting the distance carries out initialization, and except that when the maximum point of pre-treatment, all the other some pixel values all are initialized as 0; Then each point in the image is carried out twice scanning similar to range conversion, get 3 * 3 template, the pixel value that obtains after each scanning is determined by following formula:

$v_{i,j}^{\mathrm{num}}=\underset{(i,j)\&Element;\mathrm{mask}}{\max }(v_{i+k,j+l}^{\mathrm{num}-1}-t(k,l))$

Wherein

For on the image (i, j) position is in the pixel value after the scanning the num time, (k l) be the position coordinates at relative its center (0,0) of template, and t (k is a template in that (k l) locates relative distance value with the center l); If point (k, when l) being positioned at four adjoint point places on the four direction up and down of center, t (k, l)=3; If (k, when l) being positioned at four adjoint point places on the diagonal of center, t (k, l)=4; Through the image after twice scanning, have only wherein that the pixel value of point is a nonzero value near the maximum point, and the pixel value of these points is less than the pixel value of maximum point, the pixel value of other point is a null value; Find out minimum non-zero value w wherein, then the radius of this disk is r=v-w, and wherein v is the pixel value at maximum point place, and r is the radius of the corresponding disk of this maximum point; Repeating aforesaid operations is all calculated up to the disc radius of each maximum point correspondence;

The disk collection is gone redundancy: make that M is original very big disk collection, N is redundant very big disk collection; Take out the disk fm of radius maximum among the M, this disk is put into N; Then each remaining among M disk fn is compared with disk fm, get fl and be the distance between two circle disk centers, get fd=radius (fm)-0.5 * radius (fn), if fl ²＜fd ², then fn is removed from M; Repeating aforesaid operations then, is empty up to M; The N that finally obtains is the disk collection after the redundancy;

S4) in faultage image, select to be used for the disk collection of interpolation

During pixel value between the interpolation two adjacent faultage images, wherein interior each the very big disk sm of destination object profile carries out projection to another width of cloth image B successively among the piece image A, judge the projection of disk sm and be projected any disk sn that comprises in the destination object profile in the image whether lap is arranged, if lap is arranged, then sm is put into set C ₁, sn puts into set C ₂, if sn is at C ₂In exist, just need not put inward; Repeat said process then, all disks in handling A; Last C ₁With C ₂Be respectively the used disk collection of interpolation in the two adjacent faultage images;

S5) every bit points to the vector value of destination object in the image in the computed tomography image

According to step S4) be used for the disk collection of interpolation in the two adjacent faultage images that obtain, to each the some a in the two adjacent faultage images _iCarry out following calculating,

$\stackrel{\&RightArrow;}{c}=\frac{1}{\underset{d\&Element;c}{\mathrm{\&Sigma;}}\frac{r_d}{{\left|a_i{\stackrel{\&RightArrow;}{k}}_d\right|}^2}}\underset{d\&Element;c}{\mathrm{\&Sigma;}}{r}_d\frac{a_i{\stackrel{\&RightArrow;}{k}}_d}{{\left|a_i{\stackrel{\&RightArrow;}{k}}_d\right|}^2}$

Wherein, c represents to be used in the faultage image disk collection of interpolation, k _dBe the center of disk d, r _dBe the radius of disk d, a _iBe the arbitrfary point in the faultage image,

Be a vector value, expression faultage image mid point a _iDirection and distance to destination object; The result promptly obtains and corresponding two width of cloth new images of two adjacent faultage images, wherein the pixel value of every bit is a vector value, has represented the direction and the distance of this some destination object in the image;

S6) determine the distance value of spatial point between the interpolation faultage image

In the range image of two adjacent faultage images wherein a width of cloth A as starting point, traversal is each some a wherein, spend this o'clock and do straight line along the average normal direction of two faultage images and hand over second width of cloth faultage image B, connect a point two vectors corresponding in two faultage images then with b point place in 1 b

Obtain new vector i, cross the direction of a point along i then, do straight line and hand over faultage image B in a c, line segment ac is the corresponding interpolation line segment of a point;

According to a and 2 distance values of locating of c the point g on the line segment ac is carried out the cosine interpolation, formula is as follows:

$d_g=(1-\frac{\cos \frac{l_1}{l_1+l_2}\mathrm{\&pi;}+1}{2})d_1+\frac{\cos \frac{l_1}{l_1+l_2}\mathrm{\&pi;}+1}{2}{d}_2$

L wherein ₁And l ₂Be respectively g section along the line ac direction to the distance of some a with some c, d ₁, d ₂Be respectively the distance value of a point and c point place in two range images, d _gIt is the distance value at the some g place that obtains of interpolation; If when being straight line and faultage image B and not had intersection point by a point, then give up a point, to a bit repeating calculating like this down, after the every bit among the faultage image A all traveled through and finishes, the distance value of spatial point had promptly been come out by correct interpolation between the two adjacent faultage images again; Determine that according to the positive and negative of distance value this point is in profile inside or outside, and zero point is the point on the profile;

All dispose as if every adjacent two frame faultage images in the tomographic sequence that has extracted the destination object profile, then execution in step S7), otherwise continue execution in step S2)-S6);

S7) iso-surface patch based on profile shows

All spatial point distance values are extracted zero point, and carry out showing that based on the iso-surface patch of profile the entire process flow process finishes.

Images