CN101504767B - Image splitting method based on level set relay - Google Patents

Abstract

The invention discloses a sequential level set-based method for automatically segmenting images, which comprises the steps: (1) inputting an image I to be segmented and automatically initializing a level set function in an image domain by using a morphological method; (2) evolving the level set function by using an iterative approach in a corresponding subdomain to obtain the level set function and a detected image border, and using the inside domain of the edge as a next subdomain; (3) converting an image edge detected in the subdomain into a binary image and storing the binary image; (4) judging whether the area of the next subdomain is zero, and executing a step (5) if the area of the next subdomain is zero, or returning to the step (1) if the area of the next subdomain is not zero; and (5) performing a 'combination' operation of the stored binary image, combining the image edges of all subdomains, and completing the image segmentation. The method has the advantages of detecting more image edges and automatically completing image segmentation, and can be used for segmenting natural images and medical images.

Description

Image partition method based on level set relay

Technical field

The invention belongs to technical field of image processing, relate to level set (Level Set Method) image partition method, can be used for cutting apart medical image and natural image.

Background technology

Image obtains with multi-form and means observation objective world by various recording geometrys, can directly or indirectly act on human eye and and then produce the entity of vision.The human information that obtains from the external world has 75% to obtain from image approximately, and this had both illustrated that amount of image information was huge, shows that also the mankind have higher utilization factor to image information.Along with the development of signal processing theory and computer science and technology, Image Engineering also become one abundant in content and develop subject rapidly.And as the basis of Image Engineering, image segmentation is the committed step from the Flame Image Process to the graphical analysis, also is a kind of basic computer vision technique.

To the history in existing 40 years of Study of Image Segmentation, during emerged various dividing methods based on different theories, however, this field never forms design and the realization that a unified theoretical system instructs new dividing method.In recent years, initiatively consistency profiles showed in the image segmentation field and enlivened, and had produced a lot of new algorithms, as the snake method, and Level Set Method (Level set Method).Different according to the expression of evolution curve and implementation method, initiatively consistency profiles can be divided into:

1) parameter active profile (PAC, Parametric Active Contour), this class methods operation parameter curve method explicitly is represented the curve that develops, and instructs curve to develop in image by constructing corresponding energy function, and finally reaches the purpose of cutting apart.The shortcoming of this method maximum is the change in topology that must solve target in the image, and this problem is to be difficult to very much solve.

2) how much active profile (GAC, Geometric Active Contour), such method is converted into the evolution that is defined in the level set function in the higher dimensional space with the evolution in the plane of evolution curve, and the zero level collection of usage level set function (ZeroLevel Set) is implicitly represented the evolution curve.People such as Suri are at its article " J.S.Suri; K.Liu; S.Singh; S.N.Laxminarayan; X.Zeng and L.Reden.Shape recovery algorithms using level sets in 2-d/3-dmedical imagery:a state-of-the-art review.IEEE Trans.Information Technology in Biomedicine; 2002,6 (1): 8-28. " at length introduce and analyzed the advantage of Level Set Method in; wherein most important is exactly that the develop implicit representation of curve makes that Level Set Method can be simply and the topological deformation of processing target neatly; as division, merge etc.

Level Set Method the earliest is to be proposed in " S.Osher; J.A.Sethian.Fronts propagatingwith curvature-dependent speed:algorithms based on Hamilton-Jacobi formulation.J.Computational Physics; 1988; 79 (1): 12-49. " literary composition by Osher and Sethian, after this people such as Maladi has been used in it cut apart " R.Malladi; J.A.Sethian; and B.C.Vemuri.Shape modeling with frontpropagation:a level set approach.IEEE Trans.Pattern Analysis and Machine Intelligence; 1995,17 (2): 158-175. " of complex target again.There are following three subject matters in these early stage methods:

(1) need artificial initialization evolution curve, promptly manual closed curve of placement on image to be split, thereby ability initialization level set function, this process can not be finished automatically;

(2) instruct the evolution curve to stop on the edge in the image by definition edge indicator function, to reach the purpose of split image, but this indicator function can only make curve rest on the nearest edge of the initial evolution curve of distance, thus omission the initial evolution curve image border far away of distance;

(3) must be in the evolution of level set function add artificially and reinitialize link, make level set function in evolutionary process, remain the symbolic distance function.This link that reinitializes all is by finding the solution a partial differential equation aperiodically in the process of level set function iterative, and calculated amount is bigger.And Gomes also points out the very coordination of theory of this mode and level set in its article " J.Gomes and O.Faugeras.Reconciling distance functions and level sets.J.Visual Communication and ImageRepresentation; 2000,11 (2): 209-223. ".Given this, people such as Li proposed a kind of variation level diversity method " C.Li; C.Xu; C.Gui; and M.D.Fox.Level set evolution without re-initialization:a new variational formulation.IEEEConf.Computer Vision Pattern Recognition; 2005; 1:430-436. " of reinitializing of not needing in 2005 again, this method guarantees that by introduce penalty term in energy function level set function is the symbolic distance function, has avoided reinitializing process in evolution.But still there is following problem:

(1) needing the people is initialization evolution curve, thus ability initialization level set function;

(2) employing is based on the edge indicator function of gradient of image and gray scale, and zero level collection curvilinear motion direction shows as shrinks or expansion, makes this method can only detect the initial nearest edge of evolution curve of distance.

The Level Set Method that not needing of the present invention is based on that people such as Li proposes reinitializes, automatically image is divided into a series of nested subregions according to edge feature, level set function keeps sequential evolution in each subregion, and detect edge in this zone, synthetic at last detected all edges, realized the cutting apart fully of image overcome two shortcomings mentioned above.

Summary of the invention

The object of the invention is to overcome the deficiency of above-mentioned existing level set cutting techniques, a kind of image partition method and device based on level set relay proposed, not need the people be the situation of initialization evolution curve to be implemented in, automatically more edge in the detected image improves the automatization level and the segmentation effect of existing level set dividing method based on the edge.

The technical thought that realizes the object of the invention is: the usage level collection that replaces develops and nested subregion auto-initiation, makes level set function develop in the subregion of correspondence and detects respective edges.When the group region area is zero, synthesize detected edge in each sub regions, get segmentation result to the end.

Proposed a kind of image partition method, comprised the steps: based on level set relay

A. import image I to be split, and in image, utilize morphological method auto-initiation level set function;

B. to the level set function φ after the initialization _{X, y, t=0} ^r, use iterative computation that it is developed in corresponding subregion, obtain the level set function φ of stable state _{X, y} ^rAnd detected image border

And with this edge

Interior zone is as next subregion Ω ^R+1

C. with the edge

Be converted into bianry image, and this bianry image is kept in the internal memory;

D. judge next subregion Ω ^R+1Area whether be zero, if area is zero then execution in step E, otherwise return steps A, stop condition is:

Area(Ω ^r+1)＝0；

E. the bianry image of all preservations is carried out " also " operation,, finish image segmentation to merge the edge of all subregions.

Above-mentioned image partition method based on level set relay, the described morphological method auto-initiation level set function that utilizes of steps A wherein, carry out as follows:

A1. according to stable state evolution curve in the top area, obtain nested subregion Ω ^r, wherein r represents the r sub regions, top layer regions fixedly installs and is Ω ⁰=I;

A2. use morphological method initialization evolution curve

Wherein t=0 represents initially to develop curve, and its initialization formula is:

$\∂{\mathrm{\Ω}}_{t=0}^r=|\▿\mathrm{Erode}({\mathrm{\Ω}}^r,{\mathrm{\ϵ}}_d)|$

In the formula, Erode () expression corrosion operation, ε _dFor positive constant span is 3-5,

Represent gradient, || represent the amplitude of gradient;

A3. according to subregion Ω ^r, initialization evolution curve

With its interior zone Ω _T=0 ^rCalculate initial level set function φ _{X, y, t=0} ^r, computing formula is::

${\mathrm{\φ}}_{x,y,t=0}^r=\left\{\begin{array}{cc}d& (x,y)\∈{\mathrm{\Ω}}_{t=0}^r\\ 0& (x,y)\∈\∂{\mathrm{\Ω}}_{t=0}^r\\ -d& (x,y)\∈{\mathrm{\Ω}}^r\backslash {\mathrm{\Ω}}_{t=0}^r-\∂{\mathrm{\Ω}}_{t=0}^r\end{array},\right.$

In the formula,

Expression evolution curve, Ω _T=0 ^rThe zone of expression evolution curve inside,

Expression subregion Ω ^rIn and the zone of evolution extra curvature, d is that the common value of positive constant is 4.

Above-mentioned image partition method based on level set relay, wherein step C is described with the edge

Being converted into bianry image, is to generate a width of cloth and the onesize image B of image I ^r(x, y), with image B ^r(x, y) middle corresponding edge

Locational pixel assignment be 1, image B ^r(x, y) in other locational pixel assignment be 0, calculate this bianry image according to following formula:

$B^r(x,y)=\left\{\begin{array}{cc}1& (x,y)\∈\∂{\mathrm{\Ω}}^r\\ 0& (x,y)\∉\∂{\mathrm{\Ω}}^r\end{array}\right.$

In the formula,

Represent zero level collection curve in the r sub regions, promptly detected image border.

Above-mentioned image partition method based on level set relay, wherein the described bianry image to all preservations of step e is carried out " also " operation, is by detected whole edges in the following formula computed image:

$\∂\mathrm{\Ω}=\underset{r=1}{\overset{R}{\mathrm{\∪}}}\∂{\mathrm{\Ω}}^r$

In the formula,

Be illustrated in detected all edges in the image I, Represent detected edge in the r sub regions.

Image segmentation device based on level set relay provided by the invention comprises:

The auto-initiation device is used to calculate initial level set function φ _{X, y, t=0} ^r, promptly according to known image region Ω ^r, utilize formula to calculate initial evolution curve earlier: $\∂{\mathrm{\Ω}}_{t=0}^r=|\▿\mathrm{Erode}({\mathrm{\Ω}}^r,{\mathrm{\ϵ}}_d)|,$ Computing formula according to level set function calculates initial level set function φ again _{X, y, t=0} ^r, wherein t=0 represents the level set function in 0 moment, Erode () expression corrosion operation, ε _dBe positive constant, span is 3-5, Represent gradient, || represent the amplitude of gradient; The evolution device is used to carry out the iterative computation that level set function develops, and obtains the level set function φ of stable state _{X, y} ^rAnd detected image border

And with this edge Interior zone is as next subregion Ω ^R+1

Save set is used for the image border Be converted into bianry image B ^r(x y), and is kept at this bianry image in the internal memory, and the conversion formula of bianry image is:

$B^r(x,y)=\left\{\begin{array}{cc}1& (x,y)\∈\∂{\mathrm{\Ω}}^r\\ 0& (x,y)\∉\∂{\mathrm{\Ω}}^r\end{array}\right.$

In the formula, Represent zero level collection curve in the r sub regions, promptly detected image border;

Decision maker is used to judge next subregion Ω ^R+1Area whether be zero, if be zero then carry out the function of edge synthesizer, otherwise return the auto-initiation device;

The edge synthesizer is used for the bianry image of all preservations is passed through formula: $\∂\mathrm{\Ω}=\underset{r=1}{\overset{R}{\mathrm{\∪}}}\∂{\mathrm{\Ω}}^r$ Carry out " also " operation,, finish image segmentation to merge the edge of all subregions, in the formula,

Be illustrated in detected whole edges in the image,

Represent detected edge in the r sub regions.

The present invention compared with prior art has the following advantages:

(1) the present invention has avoided artificial mutual initial evolution curve owing to use morphological method auto-initiation level set function, has improved the automatization level of method;

(2) the present invention has realized cutting apart fully of image because image is used alternatingly auto-initiation and two processes of level set function;

(3) simulation result shows, the more existing level set dividing method based on the edge of the present invention has detected the more images edge.

Technical process of the present invention and effect can describe in detail in conjunction with the following drawings.

Description of drawings

Fig. 1 is a schematic flow sheet of the present invention;

Fig. 2 is a device synoptic diagram of the present invention;

Fig. 3 is the present invention and existing Li and the contrast and experiment figure of Chan-Vese method on the non-homogeneous target image;

Fig. 4 be the present invention with existing Li and Chan-Vese method the contrast and experiment figure on the nested circular image of non-homogeneous;

Fig. 5 is the present invention and existing Li and the contrast and experiment figure of Chan-Vese method on complicated artificial image;

Fig. 6 is the present invention and existing Li and the contrast and experiment figure of Chan-Vese method on the nuclear-magnetism image;

Fig. 7 is the present invention and existing Li and the contrast and experiment figure of Chan-Vese method on natural image.

Embodiment

With reference to Fig. 1, the image segmentation that the present invention is based on level set relay comprises:

Step 1, the auto-initiation level set function.

Import image I to be split, use morphological method auto-initiation evolution curve, and further initialization level set function.The auto-initiation process is divided into determines nested subregion, initialization evolution curve and initialization level set function three sub-steps, that is:

1.1, then image is made as top layer subregion, i.e. Ω if carry out the auto-initiation process for the first time ⁰=I; Otherwise subregion Ω ^rBe the zone of inside, detected edge in the top area, promptly ${\mathrm{\Ω}}^r=\mathrm{Inside}\left({\∂\mathrm{\Ω}}^{r-1}\right);$

1.2 behind the known subregion, use morphological method initialization evolution curve

Computing formula is as follows:

$\∂{\mathrm{\Ω}}_{t=0}^r=|\▿\mathrm{Erode}({\mathrm{\Ω}}^r,{\mathrm{\ϵ}}_d)|$

Wherein Erode () expression corrosion operation, ε _dRepresent a positive constant, value 3-5;

The expression compute gradient, || the amplitude of expression compute gradient;

1.3 according to subregion Ω ^r, initialization evolution curve

With its interior zone Ω _T=0 ^rCalculate initial level set function φ _{X, y, t=0} ^r, computing formula is:

${\mathrm{\φ}}_{x,y,t=0}^r=\left\{\begin{array}{cc}d& (x,y)\∈{\mathrm{\Ω}}_{t=0}^r\\ 0& (x,y)\∈\∂{\mathrm{\Ω}}_{t=0}^r\\ -d& (x,y)\∈{\mathrm{\Ω}}^r\backslash {\mathrm{\Ω}}_{t=0}^r-\∂{\mathrm{\Ω}}_{t=0}^r\end{array}\right.$

In the formula;

Expression evolution curve, Ω _T=0 ^rThe zone of expression evolution curve inside, Expression subregion Ω ^rIn and the zone of evolution extra curvature, d is that the common value of positive constant is 4.

Step 2, evolution level set function φ in subregion _{X, y, t=0} ^r

2.1 input initial level set function, by iterative formula, the steady state solution of calculated level set function:

$\frac{\∂{\mathrm{\φ}}^r}{\∂t}=\mathrm{\μ}[\mathrm{\Δ}{\mathrm{\φ}}^r-\mathrm{div}\left(\frac{\▿{\mathrm{\φ}}^r}{|\▿{\mathrm{\φ}}^r|}\right)]+\mathrm{\λ}{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^r\right)\mathrm{div}\left(g^r\frac{\▿{\mathrm{\φ}}^r}{|\▿{\mathrm{\φ}}^r|}\right)+v{g}^r{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^r\right)$

In the formula, Represent gradient, Δ is a Laplace operator, represents the divergence of gradient, || be the amplitude of gradient, φ ^rBe the level set function in the r sub regions, μ, λ and v are positive constants, δ _ε() is the Dirac function of serialization, and computing formula is as follows:

${\mathrm{\δ}}_{\mathrm{\ϵ}}\left(x\right)=\left\{\begin{array}{cc}0& \left|x\right|>\mathrm{\ϵ}\\ \frac{1}{2}[1+\cos \left(\frac{\mathrm{\πx}}{\mathrm{\ϵ}}\right)]& \left|x\right|\≤\mathrm{\ϵ}\end{array}\right.$

ε is positive constant, g ^rBe edge indicator function, computing formula is as follows:

$g^r=\frac{1}{1+{|\▿G_{\mathrm{\σ}}*I|}^2}$

Wherein

Being to use variance is that the Gaussian filter of σ is to the filtered amplitude response of image I.

As can be seen, the iterative formula of level set function is a partial differential equation, and the present invention uses explicit numerical method to find the solution, and specifically is calculated as follows:

${\mathrm{\φ}}_{x,y,t+1}^r={\mathrm{\φ}}_{x,y,t}^r+\mathrm{\Δt\μ}[\left(\frac{{\mathrm{\φ}}_{x+1,y,t}^r-2{\mathrm{\φ}}_{x,y,t}^r+{\mathrm{\φ}}_{x-1,y,t}^r}{\mathrm{\Δx}}\right)+\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-2{\mathrm{\φ}}_{x,y,t}^r+{\mathrm{\φ}}_{x,y-1,t}^r}{\mathrm{\Δy}}\right)-K^r]$

$+\mathrm{\Δt\λ}{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^r\right)[\frac{\frac{(g_{x+1,y}^r-g_{x-1,y}^r)}{2\mathrm{\Δx}}\left(\frac{{\mathrm{\φ}}_{x+1,y}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}\right)}{\sqrt{{\left(\frac{{\mathrm{\φ}}_{x+1,y}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}\right)}^2+{\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)}^2}}+\frac{\left(\frac{g_{x,y+1}^r-g_{x,y-1}^r}{2\mathrm{\Δy}}\right)\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)}{\sqrt{{\left(\frac{{\mathrm{\φ}}_{x+1,y}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}\right)}^2+{\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)}^2}}+g^r{K}^r]$

$+\mathrm{\Δtv}{g}^r{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^r\right)$

Wherein, φ _{X, y, t} ^rRepresent t level set function constantly in the position (x, the y) value on, φ _{X, y, t+1} ^rRepresent next level set function constantly in the position (x, the y) value on, φ _{X ± 1, y ± 1, t} ^rRepresent t level set function constantly in the position value on (x ± 1, y ± 1), g _{X ± 1, y ± 1} ^rFor edge indicator function in the position value on (x ± 1, y ± 1), Δ t is the time step value, Δ x and Δ y represent the step value on the horizontal direction and vertical direction on the plane respectively, carry out on image owing to develop, and generally select Δ x=Δ y=1, K ^rBe the curvature of evolution curve, calculate by following formula:

$K^r=\mathrm{div}\left(\frac{\▿{\mathrm{\φ}}^r}{|\▿{\mathrm{\φ}}^r|}\right)=\mathrm{div}\left(\frac{\frac{{\mathrm{\φ}}_{x+1,y,t}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}i+\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}j}{\sqrt{{\left(\frac{{\mathrm{\φ}}_{x+1,y,t}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}\right)}^2+{\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)}^2}}\right)$

Divergence is calculated in div () expression in the formula.

2.2 because level set function adopts the method for iteration to calculate,,, obtain steady state solution, otherwise continue iteration up to stable if satisfy condition iteration stopping so need to judge the iteration stopping condition.The iteration stopping condition is:

any(Δφ _1，...，M)＜δ _φ

Just when in M iteration, level set function with respect to the increment of the level set function of last iteration all less than given threshold value δ _φThen iteration stopping has obtained steady state solution φ _{X, y} ^r, otherwise continue iteration;

2.3 obtained the level set function φ of stable state _{X, y} ^rAfter, calculate the zero level collection curve of this level set function correspondence Also just obtained the edge in should the zone, computing formula is as follows:

$\∂{\mathrm{\Ω}}^r=\left\{(x,y)\right|{\mathrm{\φ}}_{x,y}^r(x,y)=0\};$

2.4 because the zero level collection curve of level set function correspondence is a closed curve, thus the present invention circulated as the next one in the zone of this curve inside in the input of auto-initiation step, promptly ${\mathrm{\Ω}}^{r+1}=\mathrm{inside}(\∂{\mathrm{\Ω}}^r).$

Step 3 is kept at detected edge in the subregion.Through after the evolution of level set function, the edge in the subregion is zero level collection curve just, and these edges must be preserved in order to carrying out the synthetic of edge at last.The form of preserving is a bianry image, and computing formula is as follows:

$B^r(x,y)=\left\{\begin{array}{cc}1& (x,y)\∈\∂{\mathrm{\Ω}}^r\\ 0& (x,y)\∉\∂{\mathrm{\Ω}}^r\end{array}\right.$

In the formula,

Represent zero level collection curve in the r sub regions, promptly detected image border.

Step 4 is judged next subregion Ω ^R+1Area whether be zero, whether decision also needs to generate the evolution of nested subregion and executive level set function.Concrete condition is:

Area(Ω ^r+1)＝0.

Wherein Area () represents zoning Ω ^R+1Area.If area is the null representation level set function to be developed on entire image and has detected all edges, also just there is no need to create again nested subregion, then execution in step five, otherwise jump to step 1.

Step 5 after the traversal entire image, thereby need be synthesized the edge that generates entire image in detected edge with level set in each nested subregion.Because the edge in each zone is all preserved in the mode of bianry image, so only these bianry images need be carried out " also " operation here, computing formula is as follows:

$\∂\mathrm{\Ω}=\underset{r=1}{\overset{R}{\mathrm{\∪}}}\∂{\mathrm{\Ω}}^r.$

Wherein,

Be exactly all detected edges in the image,

It is the image border in the r sub regions.

With reference to Fig. 2, the virtual bench of the image partition method based on level set relay of the present invention comprises:

The auto-initiation device is according to known subregion Ω ^rUtilize morphological method to generate initial evolution curve

Definition according to level set function generates initial level set function φ again _{X, y, t=0} ^r, wherein t=0 represents t=0 level set function constantly, i.e. initial level set function.The known subregion Ω of this device input ^r, output initialization level set function.Comprise evolution curve initial module and level set function initialization module.

Evolution curve initial module is input as subregion Ω ^r, be output as initial evolution curve Be formulated as follows:

$\∂{\mathrm{\Ω}}_{t=0}^r=|\▿\mathrm{Erode}({\mathrm{\Ω}}^r,{\mathrm{\ϵ}}_d)|$

In the formula, Erode () expression corrosion operation, ε _dRepresent positive constant, value 3-5; If carry out initialization module for the first time, the curve module that then develops is input as Ω ⁰=I is image area.

The level set function initialization module is input as subregion Ω ^r, initialization evolution curve

With its interior zone Ω _T=0 ^r, be output as initial level set function φ _{X, y, t=0} ^r, be formulated as follows:

${\mathrm{\φ}}_{x,y,t=0}^r=\left\{\begin{array}{cc}d& (x,y)\∈{\mathrm{\Ω}}_{t=0}^r\\ 0& (x,y)\∈\∂{\mathrm{\Ω}}_{t=0}^r\\ -d& (x,y)\∈{\mathrm{\Ω}}^r\backslash {\mathrm{\Ω}}_{t=0}^r-\∂{\mathrm{\Ω}}_{t=0}^r\end{array}\right.$

In the formula,

Expression evolution curve, Ω _T=0 ^rThe zone of expression evolution curve inside, Expression subregion Ω ^rIn and the zone of evolution extra curvature, d is that the common value of positive constant is 4.

The evolution device is used for the iterative computation that level set function develops, and preserving level set function, to develop stable be corresponding zero level collection curve Wherein r represents the r sub regions.This device input initial level set function φ _{X, y, t=0} ^r, export stable zero level collection curve Nested subregion Ω with the next one ^R+1Comprise iterative computation module, zero level collection curve calculation module and subregion generation module.

The iterative computation module is input as initial level set function φ _{X, y, t=0} ^r, output steady-state level set function φ _{X, y} ^r, computing formula is:

${\mathrm{\φ}}_{x,y,t+1}^r=$

${\mathrm{\φ}}_{x,y,t}^r+\mathrm{\Δt\μ}[\left(\frac{{\mathrm{\φ}}_{x+1,y,t}^r-2{\mathrm{\φ}}_{x,y,t}^r+{\mathrm{\φ}}_{x-1,y,t}^r}{\mathrm{\Δx}}\right)+\frac{{\mathrm{\φ}}_{x,y+1,t}^r-2{\mathrm{\φ}}_{x,y,t}^r+{\mathrm{\φ}}_{x,y-1,t}^r}{\mathrm{\Δy}}-K^r]$

$+\mathrm{\Δt\λ}{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^r\right)[\frac{\frac{(g_{x+1,y}^r-g_{x-1,y}^r)}{2\mathrm{\Δx}}\left(\frac{{\mathrm{\φ}}_{x+1,y}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}\right)}{\sqrt{{\left(\frac{{\mathrm{\φ}}_{x+1,y}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}\right)}^2+{\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)}^2}}+\frac{\left(\frac{g_{x,y+1,t}^r-g_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)}{\sqrt{{\left(\frac{{\mathrm{\φ}}_{x+1,y}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}\right)}^2+{\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)}^2}}+g^r{K}^r]$

$+\mathrm{\Δtv}{g}^r{\mathrm{\δ}}_{\mathrm{\ϵ}}\left({\mathrm{\φ}}^r\right)$

In the formula, φ _{X, y, t} ^rRepresent t level set function constantly in the position (x, y) value on; φ _{X, y, t+1} ^rRepresent next level set function constantly in the position (x, y) value on; φ _{X ± 1, y ± 1, t} ^rRepresent t level set function constantly in the position value on (x ± 1, y ± 1); Δ t is the time step value, and Δ x and Δ y represent the step value on the horizontal direction and vertical direction on the plane respectively, carry out on image owing to develop, and generally select Δ x=Δ y=1, K ^rBe the curvature of evolution curve, calculate by following formula:

$K^r=\mathrm{div}\left(\frac{\▿{\mathrm{\φ}}^r}{|\▿{\mathrm{\φ}}^r|}\right)=\mathrm{div}\left(\frac{\frac{{\mathrm{\φ}}_{x+1,y,t}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}i+\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}j}{\sqrt{{\left(\frac{{\mathrm{\φ}}_{x+1,y,t}^r-{\mathrm{\φ}}_{x-1,y,t}^r}{2\mathrm{\Δx}}\right)}^2+{\left(\frac{{\mathrm{\φ}}_{x,y+1,t}^r-{\mathrm{\φ}}_{x,y-1,t}^r}{2\mathrm{\Δy}}\right)}^2}}\right)$

Divergence is calculated in div () expression in the formula.

Zero level collection curve calculation module is input as steady-state level set function φ _{X, y} ^r, be output as zero level collection curve

Computing formula is:

$\∂{\mathrm{\Ω}}^r=\left\{(x,y)\right|{\mathrm{\φ}}_{x,y}^r(x,y)=0\}$

The subregion generation module is input as zero level collection curve Be output as next subregion ${\mathrm{\Ω}}^{r+1}=\mathrm{inside}(\∂{\mathrm{\Ω}}^r).$

Preserve module, be input as and detect the edge, be i.e. zero level collection curve

Be output as bianry image, computing formula is as follows:

$B^r(x,y)=\left\{\begin{array}{cc}1& (x,y)\∈\∂{\mathrm{\Ω}}^r\\ 0& (x,y)\∉\∂{\mathrm{\Ω}}^r\end{array}\right.$

In the formula,

Represent zero level collection curve in the r sub regions, promptly detected image border;

Decision maker is used to judge whether cycling condition satisfies, and curve promptly develops

Whether the inner region area is zero, continues to carry out if condition satisfies, otherwise continues to generate subregion and carry out the level set iteration.The mathematical notation of decision condition is as follows:

Area(Ω ^r)＝0.

Whether the area of promptly judging subregion is zero, wherein the area of Area () expression zoning.

The edge synthesizer is used for the edge in the synthetic entire image in detected edge, each zone.Because the edge is to preserve with the form of bianry image, so adopt " also " operation to synthesize the edge.The edge of this device input picture saves as the bianry image form, is output as the edge of entire image.Be formulated as:

$\∂\mathrm{\Ω}=\underset{r=1}{\overset{R}{\mathrm{\∪}}}\∂{\mathrm{\Ω}}^r.$

Wherein,

Be exactly all detected edges in the image,

It is the image border in the r sub regions.

Below validity and practicality by emulation experiment checking the inventive method.

The emulation content:

(1) employing contrast experiment's form selects two representative level set image segmentation methods to test on same image, to verify validity of the present invention.It is concrete that what select is the method that people such as Li proposes, concrete list of references " C.Li; C.Xu; C.Gui; and M.D.Fox.Level set evolution without re-initialization:a newvariational formulation.IEEE Conf.Computer Vision Pattern Recognition; 2005; 1:430-436. ", method (abbreviating the Chan-Vese method as) with Chan and Vese, concrete list of references " T.F.Chan; L.A.Vese.Active contours without edges.IEEE Trans.Image Processing, 2001,10 (2): 266-277. ");

(2) use artificial image, medical image and natural image carry out emulation experiment, with the segmentation effect of checking the present invention to different images.Concrete simulated conditions sees the description of each experiment for details.

Experiment one, image to be split is 256 * 256 artificial image, and image is made up of three target and backgrounds, and wherein the gray scale of triangular day mark and background gray scale are approaching.Triplex row is represented method, Chan-Vese method and the level set function evolutionary process of the present invention of Li respectively from top to bottom in Fig. 4.Simulation result shows: the method omission of Li the inward flange of annular target, the evolution that in a single day velocity function that reason is based on shade of gray reaches the edge level set function will stop, so the method for Li can only detect the edge nearest apart from initial curve; The present invention is cut apart on a plurality of subregions, so after the method for Li stops, still can continuing to cut apart in the zone of inside, has detected whole edges; The omission of Chan-Vese method triangular day mark, its essential reason is that Chan-Vese algorithm hypothesis image is that target and background by homogeneity constitutes, and with image segmentation linear expression is:

I＝α ₁H(φ)+α ₂(1-H(φ))

In the formula, H () step function, the zone in H (φ) the expression curve; The zone of 1-H (φ) expression extra curvature; φ is a level set function; a ₁And a ₂It is respectively the gray scale barycenter of the pixel of zero level collection inside or outside of curve.Because triangular day mark and background gray scale are approaching,, be used as background so be considered to belong to the part of extra curvature.From top analysis as can be seen, under the situation of target nonuniformity, the present invention has better segmentation effect than side and the Chan-Vese of Li in image, detects whole edges in the image preferably.

Experiment two, image to be split is 256 * 256 artificial image, image is made of the nested circle of three different gray scales, among Fig. 5 from top to bottom triplex row represent method, Chan-Vese method and the level set function evolutionary process of the present invention of Li respectively.Simulation result shows that the method convergent speed of Li is very fast, but only detected the edge of close initial curve, omission inner edge; The Chan-Vese method detects two-layer edge, for experiment in 1 same reason omission the edge of center white object; The present invention is because carry out repeated segmentation at different subregions to image, so detected whole edges in the image.

Experiment three, image to be split is 256 * 256 artificial image, image is made up of the complex target of some difformities and gray scale and the literal of gray scale gradual change, and Fig. 6 triplex row from top to bottom represents method, Chan-Vese method and the level set function evolutionary process of the present invention of Li respectively.Simulation result shows: the method for Li still only detects the edge nearest apart from initial curve, omission than general objective, as the edge of oval and square target and alphabetical p and e inside; The Chan-Vese method is cut apart according to the half-tone information of pixel in the entire image, has detected all than general objective and the less target of part, omission the alphabetical e close and square and ellipse target internal edge with the background gray scale; The present invention is repeated segmentation in different subregions, has considered the local feature of subregion, so detected whole edges in the image.

Experiment four, image to be split is 181 * 217 human brain nuclear magnetic resonance image, and image is the axial bitmap of human brain, and tissue is respectively fat, brainpan and intracranial tissue from outside to inside.Fig. 7 triplex row from top to bottom represents method, Chan-Vese method and the level set function evolutionary process of the present invention of Li respectively.Simulation result shows: the method for Li only detects outermost edge between brain and the background; The Chan-Vese method contains the background of a large amount of black owing in the image, so also only detect the outermost edge of brain; The present invention has detected more internal edge.

Experiment five, image to be split is 192 * 256 natural image, in the image be to put the hand on the desktop and the object of mouse, and uneven illumination is even.Fig. 8 triplex row from top to bottom represents method, Chan-Vese method and the level set function evolutionary process of the present invention of Li respectively.Simulation result shows: though the method for Li detects left hand and mouse pad edge, omission the edge on the right hand and left hand sleeve limit; The Chan-Vese method detects the right hand portion edge, but omission whole left hand, this method is not fine for the even image segmentation effect of uneven illumination; The present invention is complete has detected right-hand man and other edge.

Experimental result shows, invents more existing Level Set Method and can detect more edge, especially can detect more internal edge; The image even to uneven illumination also has tangible segmentation effect; And medical image and natural image all had segmentation effect and practicality preferably.

Claims (5)

1. the image partition method based on level set relay comprises the steps:

A. import image I to be split, and in image I to be split, utilize morphological method auto-initiation level set function,, then image to be split is made as top layer subregion, i.e. Ω if carry out the auto-initiation process for the first time ⁰=I; Otherwise subregion Ω ^rBe the zone of inside, detected edge in the top area, promptly

B. to the level set function after the initialization

, use iterative computation that it is developed in the subregion that the expression current iteration is used, obtain the level set function of stable state

And detected image border

, and with this edge

Interior zone is as next subregion Ω ^R+1

C. with the edge Be converted into bianry image, and this bianry image is kept in the internal memory;

D. judge next subregion Ω ^R+1Area whether be zero, if area is zero then execution in step E, otherwise at regional Ω ^R+1In utilize morphological method auto-initiation level set function, stop condition is:

Area(Ω ^r+1)＝0；

E. the bianry image of all preservations is carried out " also " operation,, finish image segmentation to merge the edge of all subregions.

2. according to right 1 described image partition method based on level set relay, the described morphological method auto-initiation level set function that utilizes of steps A wherein, carry out as follows:

2a., obtain subregion Ω according to stable state evolution curve in the top area ^r, wherein r represents the r sub regions, top layer regions fixedly installs and is Ω ⁰=I;

2b. use morphological method initialization evolution curve

, wherein t=0 represents initially to develop curve, its initialization formula is:

${\∂\mathrm{\Ω}}_{t=0}^r=|\▿\mathrm{Erode}({\mathrm{\Ω}}^r,{\mathrm{\ϵ}}_d)|$

In the formula, Erode () expression corrosion operation, ε _dFor positive constant span is 3-5, Represent gradient, || represent the amplitude of gradient;

2c. according to subregion Ω ^r, initialization evolution curve With its interior zone

Calculate the initial level set function

, computing formula is:

${\mathrm{\φ}}_{x,y,t=0}^r=\left\{\begin{array}{cc}d& (x,y)\∈{\mathrm{\Ω}}_{t=0}^r\\ d& (x,y)\∈\∂{\mathrm{\Ω}}_{t=0}^r\\ -d& (x,y)\∈{\mathrm{\Ω}}^r\backslash {\mathrm{\Ω}}_{t=0}^r-{\∂\mathrm{\Ω}}_{t=0}^r\end{array}\right.,$

In the formula, Expression evolution curve,

The zone of expression evolution curve inside,

Expression subregion Ω ^rIn and the zone of evolution extra curvature, d is that positive constant value is 4.

3. according to right 1 described image partition method based on level set relay, wherein step C is described with the edge

Being converted into bianry image, is to generate a width of cloth and the onesize image B of image I to be split ^r(x, y), with image B ^r(x, y) middle corresponding edge

Locational pixel assignment be 1, image B ^r(x, y) in other locational pixel be 0, calculate its bianry image according to following formula:

$B^r(x,y)=\left\{\begin{array}{cc}1& (x,y)\∈{\∂\mathrm{\Ω}}^r\\ 0& (x,y)\∉{\∂\mathrm{\Ω}}^r\end{array}\right.$

In the formula, Represent zero level collection curve in the r sub regions, promptly detected image border.

4. according to right 1 described image partition method based on level set relay, wherein the described bianry image to all preservations of step e is carried out " also " operation, is by detected all edges in the following formula computed image:

$\∂\mathrm{\Ω}=\underset{r=1}{\overset{R}{\∪}}{\∂\mathrm{\Ω}}^r.$

In the formula, Be illustrated in detected all edges in the image I to be split, Represent detected edge in the r sub regions.

5. image segmentation device based on level set relay comprises:

The auto-initiation device is used to calculate the initial level set function

, promptly earlier according to subregion Ω ^r, calculate initial evolution curve: Wherein t=0 represents the level set function in 0 moment, Erode () expression corrosion operation, ε _dBe positive constant, span is 3-5,

Represent gradient, || represent the amplitude of gradient; Again according to subregion Ω ^r, initialization evolution curve

With its interior zone

Calculate the initial level set function

, computing formula is:

${\mathrm{\φ}}_{x,y,t=0}^r=\left\{\begin{array}{cc}d& (x,y)\∈{\mathrm{\Ω}}_{t=0}^r\\ d& (x,y)\∈\∂{\mathrm{\Ω}}_{t=0}^r\\ -d& (x,y)\∈{\mathrm{\Ω}}^r\backslash {\mathrm{\Ω}}_{t=0}^r-{\∂\mathrm{\Ω}}_{t=0}^r\end{array}\right.,$

In the formula,

Expression evolution curve,

The zone of expression evolution curve inside,

Expression subregion Ω ^rIn and the zone of evolution extra curvature, d is that positive constant value is 4;

The evolution device is used to carry out the iterative computation that level set function develops, and obtains the level set function of stable state

And detected image border

, and with this edge

Interior zone is as next subregion Ω ^R+1

Save set is used for the image border

Be converted into bianry image, and this bianry image is kept in the internal memory, the conversion formula of bianry image is:

$B^r(x,y)=\left\{\begin{array}{cc}1& (x,y)\∈{\∂\mathrm{\Ω}}^r\\ 0& (x,y)\∉{\∂\mathrm{\Ω}}^r\end{array}\right.$

In the formula,

Represent zero level collection curve in the r sub regions, promptly detected image border, image B ^r(x, y) middle corresponding edge Locational pixel assignment be 1, image B ^r(x, y) in other locational pixel assignment be 0;

Decision maker is used to judge next subregion Ω ^R+1Area whether be zero, if be zero then carry out the function of edge synthesizer, otherwise return the auto-initiation device;

The edge synthesizer is used for the bianry image of all preservations is passed through formula:

Carry out " also " operation,, finish image segmentation to merge the edge of all subregions, in the formula,

Be illustrated in detected all edges in the image,

Represent detected edge in the r sub regions.

Images