CN111028201A - Fundus blood vessel image segmentation system and method based on multi-scale level set - Google Patents

Abstract

The invention provides a fundus blood vessel image segmentation system and method based on a multi-scale level set, and relates to the technical field of image processing. The invention comprises an input module, a preprocessing module, a scale construction module, a level set module and an output module; the present invention improves upon existing level set methods based on local features. The fundus blood vessels are different in thickness, and the blood vessel segmentation precision is not high due to the fact that the blood vessel segmentation granularity is controlled by a single scale. Based on the characteristic that the blood vessels are different in thickness, the optimal response scale is found for each pixel point, namely the scales of the blood vessels different in thickness are different, and the problems of missing, wrong and the like of the blood vessels with different thicknesses are solved, so that the robustness of the segmentation method is enhanced, and an accurate image segmentation result can be obtained by applying the method.

Description

Fundus blood vessel image segmentation system and method based on multi-scale level set

Technical Field

The invention relates to the technical field of image processing, in particular to a fundus blood vessel image segmentation system and method based on a multi-scale level set.

Background

The eyeground blood vessel image is the basis for screening eye diseases such as cataract, glaucoma, age-related macular degeneration and the like. The problems of complicated structure of the fundus blood vessels and uneven image gray scale cause the problems of insufficient segmentation and the like of a plurality of algorithms at the adjacent blood vessels. On the other hand, the fundus blood vessels are different in thickness, and the blood vessel segmentation precision is not high due to the fact that the blood vessel segmentation granularity is controlled by a single scale in the existing segmentation method. Therefore, the invention provides a fundus blood vessel image segmentation method based on a multi-scale level set. And selecting the optimal scale for each pixel point, and overcoming the problems of missed and wrong division of uneven-thickness blood vessels.

Disclosure of Invention

Aiming at the problems in the prior art, the invention provides a fundus blood vessel image segmentation system based on a multi-scale level set.

The technical scheme adopted by the invention is as follows: a fundus blood vessel image segmentation system based on a multi-scale level set comprises an input module, a preprocessing module, a scale construction module, a level set module and an output module;

the input module is used for receiving the fundus blood vessel image input by the user and then transmitting the fundus blood vessel image to the preprocessing module;

the preprocessing module is used for expanding the received fundus blood vessel image, and the expanded result is transmitted to the scale construction module;

the scale construction module is used for finding the optimal blood vessel response scale for the blood vessel image at the bottom of the eye, further constructing a Gaussian curve based on the scale, subtracting the original image from the convolution result of the Gaussian curve and the image to obtain an image M, and transmitting the image M to the level set module;

the level set module is used for carrying out blood vessel segmentation on the image M by adopting a level set algorithm and transmitting a segmentation result to the output module;

and the output module is used for outputting and displaying the segmentation result image.

In another aspect, the present invention provides a fundus blood vessel image segmentation method based on a multi-scale level set, which is implemented by the fundus blood vessel image segmentation system based on a multi-scale level set, and includes the following steps:

step 1: the input module receives an image A of fundus blood vessels to be segmented, wherein the width of the image A is W, and the height of the image A is H;

step 2: the preprocessing module constructs the received fundus blood vessel image A to be segmented to generate a fundus blood vessel image I without an edge background, and the resolution of the fundus blood vessel image I is the same as that of the image A; converting the image A into a gray image J, setting all areas of the image J as omega, directly obtaining an image P if the area to be segmented is full of omega, carrying out mirror image expansion on an effective boundary area in the image to cover a black filling area in the whole image if the area to be segmented is not full of omega, setting the expanded image as P, and executing the step 2.1 if I is equal to P; otherwise, executing step 3;

step 2.1: determining a threshold Q to divide the grayscale image J into a foreground R and a background B, where R and B are binary images, R is 1 and B is 0 in the foreground region, R is 0 and B is 1 in the background region;

step 2.2: constructing a foreground region boundary image R₂Etching the foreground R to obtain R₁，R₂＝R-R₁Wherein the image R₂The image is a binary image;

step 2.3: initializing a variable P to 1, wherein P is an initializing abscissa value of a pixel point in the outer loop;

step 2.4: if p is less than or equal to W and the initialized longitudinal coordinate value q of the pixel point in the outer loop is 1, executing the step 2.5, otherwise, executing the step 3;

step 2.5: if q is less than or equal to H, executing step 2.6, otherwise, making p equal to p +1, and executing step 2.4;

step 2.6: if the current point (p, q) satisfies B ═ 1, step 2.7 is performed, otherwise, let q ═ q +1, step 2.5 is performed;

step 2.7: initializing a variable m which is 1 and is an initial abscissa value of a pixel point in the image inner loop, and initializing

Wherein r is the distance between any two pixels in the initialized image, and the initial coordinates (a, b) are two random positive integers;

step 2.8: if m is less than or equal to W, the variable n is 1, n is the initial value of the vertical coordinate of the pixel point in the inner loop, the step 2.9 is executed, otherwise, the step 2.13 is executed;

step 2.9: if n is less than or equal to H, executing step 2.10, otherwise, making m equal to m +1, and executing step 2.8;

step 2.10: if the point (m, n) satisfies R₂Step 2.11 is performed when n is equal to n +1, otherwise step 2.9 is performed;

step 2.11: calculating a distance dist ((p, q), (m, n)) between the point (p, q) and the point (m, n), wherein dist is a distance metric function, and if r < dist ((p, q), (m, n)), then r ═ dist ((p, q), (m, n)), and let a ═ m, b ═ n;

step 2.12: step 2.9 is skipped to by making n equal to n + 1;

step 2.13: the point (p, q) takes the point (a, b) as the center, R is a symmetrical distance, a corresponding pixel point L in R is obtained, and the value of the pixel point L is assigned to the current point Z; wherein L is_x＝2a-p，L_y2b-q, wherein L_xIs the abscissa L of the pixel point L_yIs the ordinate of the pixel point L;

step 2.14: q is q +1, and the step 2.5 is skipped;

and step 3: for parameter minimum scale sigma_minMaximum scale σ_maxA constant coefficient gamma, a constant coefficient c, a scale step length delta sigma are set, wherein sigma is_min＜σ_max；

And 4, step 4: for the initial scale sigma_iIs set so that σ_i＝σ_min；

And 5: constructing a Hessian matrix of the image I, constructing a Gaussian function through the optimal response scale sigma, and subtracting the image I from the result of convolution operation of the image I and the Gaussian function to obtain an image M, wherein the method comprises the following specific steps:

step 5.1: initialization variable p 1, scale variable k 1, variable v_max＝0, wherein v_maxIs the maximum response value variable of the pixel;

step 5.2: if p is less than or equal to W, the variable q is 1, and step 5.3 is executed, otherwise, step 6 is executed;

step 5.3: if q is less than or equal to H, executing step 5.4, otherwise, making p equal to p +1, and executing step 5.2;

step 5.4: judging the current scale sigma_iWhether or not greater than the maximum scale sigma_maxIf it is greater than σ_maxIf so, let σ_i＝σ_minExecuting step 5.12, otherwise, executing step 5.5;

step 5.5: constructing a Hessian matrix for point (p, q):

wherein ,

D＝∫∫I(p-u,q-v)×A(u,v)dudv

E＝∫∫I(p-u,q-v)×B(u,v)dudv

F＝∫∫I(p-u,q-v)×C(u,v)dudv

a (u, v), B (u, v) and C (u, v) are based on the scale σ_iConstructing a second partial derivative of Gaussian

Wherein D, E, F is the convolution of image I with the corresponding gaussian second order partial derivative A, B, C;

step 5.6: calculating two eigenvalues lambda of the Hessian matrix obtained in step 5.5₁、λ₂, wherein

wherein

Step 5.7: calculated at the current scale sigma_iResponse value of pixel (p, q)

wherein

in the formula ,R_BIs the degree of anisotropy of the linear shape, S is the norm;

step 5.8: if it is

Then order

And κ ═ σ -_i；

Step 5.9: constructing a Gaussian kernel based on the scale k to obtain a Gaussian curve g_κ；

Step 5.10: step 5.9The resulting Gaussian curve g_κPerforming convolution operation with the image I, and subtracting the image I from a convolution result to obtain an image M;

M(p,q)＝(g_κ*I)(p,q)-I(p,q)

step 5.11: updating the scale sigma_i＝σ_i+ Δ σ, Δ σ is the scale step, and step 5.4 is skipped;

step 5.12: updating the coordinate q to be q +1, and jumping to the step 5.3;

step 6, setting an initial zero level set segmentation curve on the image M, and generating a high one-dimensional level set function phi according to the set zero level set⁰；

The set initial zero level set segmentation curve is a closed curve with a specified shape and is a randomly generated closed curve, and the level set function phi is a high one-dimensional symbol distance function generated according to the set initial level set segmentation curve;

step 7, establishing an energy functional E (phi) of level set segmentation, and setting a weighting coefficient α in the level set method₁Weighting factor α₂Curve smoothing term coefficient v and area smoothing term coefficient mu;

the energy functional is as follows:

E(φ)＝E_D+vL(φ)+μP(φ)

wherein ,E_DFor the image data item, L (φ) is a length constraint term, P (φ) is an area constraint term, and H (φ (P, q)) is a Heaviside function of φ (P, q);

E_D(φ)＝α₁∫_Ω|M(p,q)-c₁|²H(φ(p,q))dpdq+α₂∫_Ω|M(p,q)-c₂|²(1-H(φ(p,q)))dpdq

P(φ)＝∫H(φ(p,q))dpdq

wherein c₁Is the mean value of the gray levels in the zero level set segmentation curve, c₂Is the gray average outside the zero level set segmentation curve;

and 8: setting the maximum iteration number J and a convergence threshold epsilon of an energy functional of the image level set;

and step 9: separately updating the image level set function set H (phi)^j)、1-H(φ^j)、C^j、φ^jObtaining an energy functional E of the j iteration^j；

Specifically, the method comprises the following steps:

where at is the step of time in which,

for the time partial derivative of the iteratively updated level set function, the time partial derivative formula is as follows:

step 10: judging whether an iteration termination condition is reached: if the current iteration time J reaches the maximum iteration time J or the difference value of the iteration results of two adjacent times of the energy functional E (phi) is smaller than a convergence threshold epsilon which is set in advance, executing the step 11, otherwise, adding 1 to the iteration time, and returning to the step 9;

step 11: according to the current updated image level set function phi^jAn image, i.e. the segmentation result of the image, is constructed.

Adopt the produced beneficial effect of above-mentioned technical scheme to lie in:

the invention provides a segmentation system and method of fundus blood vessel images based on a multi-scale level set, and improves the existing level set method based on local features. The fundus blood vessels are different in thickness, and the blood vessel segmentation precision is not high due to the fact that the blood vessel segmentation granularity is controlled by a single scale. Based on the characteristic that the blood vessels are different in thickness, the optimal response scale is found for each pixel point, namely the scales of the blood vessels different in thickness are different, and the problems of missing, wrong and the like of the blood vessels with different thicknesses are solved, so that the robustness of the segmentation method is enhanced. By applying the method, an accurate image segmentation result can be obtained.

Drawings

FIG. 1 is a block diagram of a fundus blood vessel image segmentation system;

FIG. 2 is a flowchart of a fundus blood vessel image segmentation method

Fig. 3 shows an original image of fundus blood vessel image input;

FIG. 4 shows the conversion of the original image of the blood vessels in the fundus into a gray scale image;

FIG. 5 is a view after expansion of a fundus blood vessel image;

FIG. 6 shows an enhanced image;

FIG. 7 shows an image after convolution and subtraction of the original image;

FIG. 8 shows an initial level set;

FIG. 9 shows the final segmentation result of the image;

FIG. 10 shows a gold standard for an image;

Detailed Description

The following detailed description of embodiments of the present invention is provided in connection with the accompanying drawings and examples. The following examples are intended to illustrate the invention but are not intended to limit the scope of the invention.

A fundus blood vessel image segmentation system based on a multi-scale level set is shown in figure 1 and comprises an input module, a preprocessing module, a scale construction module, a level set module and an output module;

the input module is used for receiving the fundus blood vessel image input by the user and then transmitting the fundus blood vessel image to the preprocessing module;

the preprocessing module is used for expanding the received fundus blood vessel image, and the expanded result is transmitted to the scale construction module;

the scale construction module is used for finding the optimal blood vessel response scale for the blood vessel image at the bottom of the eye, further constructing a Gaussian curve based on the scale, subtracting the original image from the convolution result of the Gaussian curve and the image to obtain an image M, and transmitting the image M to the level set module;

the level set module is used for carrying out blood vessel segmentation on the image M by adopting a level set algorithm and transmitting a segmentation result to the output module;

and the output module is used for outputting and displaying the segmentation result image.

In another aspect, the present invention provides a fundus blood vessel image segmentation method based on a multi-scale level set, which is implemented by the foregoing fundus blood vessel image segmentation system based on a multi-scale level set, as shown in fig. 2, and includes the following steps:

step 1: the input module receives an image A of fundus blood vessels to be segmented, wherein the image A has the width W and the height H as shown in FIG. 3;

step 2: the preprocessing module constructs the received fundus blood vessel image A to be segmented to generate a fundus blood vessel image I without an edge background, and the resolution of the fundus blood vessel image I is the same as that of the image A; converting the image a into a gray image J, as shown in fig. 4, setting all regions of the image J to be Ω, directly obtaining an image P if the region to be segmented is full of Ω, and performing mirror image expansion on an effective boundary region in the image to cover a black filled region in the whole image if the region to be segmented is not full of Ω, and setting the expanded image to be P as shown in fig. 5, and meanwhile, executing step 2.1 if I is equal to P; otherwise, executing step 3;

step 2.1: determining a threshold Q to divide the grayscale image J into a foreground R and a background B, where R and B are binary images, R is 1 and B is 0 in the foreground region, R is 0 and B is 1 in the background region;

step 2.2: constructing a foreground region boundary image R₂Etching the foreground R to obtain R₁，R₂＝R-R₁Wherein the image R₂The image is a binary image;

step 2.3: initializing a variable P to 1, wherein P is an initializing abscissa value of a pixel point in the outer loop;

step 2.4: if p is less than or equal to W and the initialized longitudinal coordinate value q of the pixel point in the outer loop is 1, executing the step 2.5, otherwise, executing the step 3;

step 2.5: if q is less than or equal to H, executing step 2.6, otherwise, making p equal to p +1, and executing step 2.4;

step 2.6: if the current point (p, q) satisfies B ═ 1, step 2.7 is performed, otherwise, let q ═ q +1, step 2.5 is performed;

step 2.7: initializing a variable m which is 1 and is an initial abscissa value of a pixel point in the image inner loop, and initializing

Wherein r is the distance between any two pixels in the initialized image, and the initial coordinates (a, b) are two random positive integers;

step 2.8: if m is less than or equal to W, the variable n is 1, n is the initial value of the vertical coordinate of the pixel point in the inner loop, the step 2.9 is executed, otherwise, the step 2.13 is executed;

step 2.9: if n is less than or equal to H, executing step 2.10, otherwise, making m equal to m +1, and executing step 2.8;

step 2.10: if the point (m, n) satisfies R₂Step 2.11 is performed when n is equal to n +1, otherwise step 2.9 is performed;

step 2.11: calculating a distance dist ((p, q), (m, n)) between the point (p, q) and the point (m, n), wherein dist is a distance metric function, and if r < dist ((p, q), (m, n)), then r ═ dist ((p, q), (m, n)), and let a ═ m, b ═ n;

step 2.12: step 2.9 is skipped to by making n equal to n + 1;

step 2.13: the point (p, q) takes the point (a, b) as the center, R is a symmetrical distance, a corresponding pixel point L in R is obtained, and the value of the pixel point L is assigned to the current point Z; wherein L is_x＝2a-p，L_y2b-q, wherein L_xIs the abscissa L of the pixel point L_yIs the ordinate of the pixel point L;

step 2.14: q is q +1, and the step 2.5 is skipped;

and step 3: for parameter minimum scale sigma_minMaximum scale σ_maxA constant coefficient gamma, a constant coefficient c, a scale step length delta sigma are set, wherein sigma is_min＜σ_max；

And 4, step 4: for the initial scale sigma_iIs set so that σ_i＝σ_min；

And 5: constructing a Hessian matrix of the image I, constructing a gaussian function through the optimal response scale σ, and subtracting the image I from a result of convolution operation of the image I and the gaussian function to obtain an image M, as shown in fig. 6 and 7, the specific steps are as follows:

step 5.1: initialization variable p 1, scale variable k 1, variable v_max＝0, wherein v_maxIs the maximum response value variable of the pixel;

step 5.2: if p is less than or equal to W, the variable q is 1, and step 5.3 is executed, otherwise, step 6 is executed;

step 5.3: if q is less than or equal to H, executing step 5.4, otherwise, making p equal to p +1, and executing step 5.2;

step 5.4: judging the current scale sigma_iWhether or not greater than the maximum scale sigma_maxIf it is greater than σ_maxIf so, let σ_i＝σ_minExecuting step 5.12, otherwise, executing step 5.5;

step 5.5: constructing a Hessian matrix for point (p, q):

wherein ,

D＝∫∫I(p-u,q-v)×A(u,v)dudv

E＝∫∫I(p-u,q-v)×B(u,v)dudv

F＝∫∫I(p-u,q-v)×C(u,v)dudv

a (u, v), B (u, v) and C (u, v) are based on the scale σ_iConstructing a second-order partial derivative of Gaussian;

wherein D, E, F is the convolution of image I with the corresponding gaussian second order partial derivative A, B, C;

step 5.6: calculating two eigenvalues lambda of the Hessian matrix obtained in step 5.5₁、λ₂, wherein

wherein

Step 5.7: calculated at the current scale sigma_iResponse value of pixel (p, q)

wherein

in the formula ,R_BIs the degree of anisotropy of the linear shape, S is the norm;

step 5.8: if it is

Then order

And κ ═ σ -_i；

Step 5.9: constructing a Gaussian kernel based on the scale k to obtain a Gaussian curve g_κ；

Step 5.10: the Gaussian curve g obtained in the step 5.9_κPerforming convolution operation with the image I, and subtracting the image I from a convolution result to obtain an image M;

M(p,q)＝(g_κ*I)(p,q)-I(p,q)

step 5.11: updating the scale sigma_i＝σ_i+ Δ σ, Δ σ is the scale step, and step 5.4 is skipped;

step 5.12: updating the coordinate q to be q +1, and jumping to the step 5.3;

step 6, setting an initial zero level set segmentation curve on the image M, as shown in FIG. 8, generating a high one-dimensional level set function phi according to the set zero level set⁰；

The set initial zero level set segmentation curve is a closed curve with a specified shape and is a randomly generated closed curve, and the level set function phi is a high one-dimensional symbol distance function generated according to the set initial level set segmentation curve;

step 7, establishing an energy functional E (phi) of level set segmentation, and setting a weighting coefficient α in the level set method₁Weighting factor α₂Curve smoothing term coefficient v and area smoothing term coefficient mu;

the energy functional is as follows:

E(φ)＝E_D+vL(φ)+μP(φ)

wherein ,E_DFor the image data item, L (φ) is a length constraint term, P (φ) is an area constraint term, and H (φ (P, q)) is a Heaviside function of φ (P, q);

E_D(φ)＝α₁∫_Ω|M(p,q)-c₁|²H(φ(p,q))dpdq+α₂∫_Ω|M(p,q)-c₂|²(1-H(φ(p,q)))dpdq

P(φ)＝∫H(φ(p,q))dpdq

wherein c₁Is the mean value of the gray levels in the zero level set segmentation curve, c₂Is the gray average outside the zero level set segmentation curve;

in the embodiment, the larger the image data item coefficient is, the closer the zero level set segmentation curve is to the boundary of the object to be segmented in the image, and the larger the level set segmentation curve length item coefficient is, the smoother the segmentation curve is, so that the weighting coefficient α, which is the control parameter for image level set segmentation, is set according to the image to be segmented₁1.0, weighting factor α₂The coefficient v of the curve smoothing term is 0.001 × 255, and the coefficient μ of the area smoothing term is 2.0.

And 8: setting the maximum iteration number J and a convergence threshold epsilon of an energy functional of the image level set;

in this embodiment, the maximum number of iterations J is set to 50, and the energy functional E (c) of image level set segmentation is set₁,c₂Phi), the convergence threshold is 0.01.

And step 9: separately updating the image level set function set H (phi)^j)、1-H(φ^j)、C^j、φ^jObtaining an energy functional E of the j iteration^j；

Specifically, the method comprises the following steps:

where at is the step of time in which,

for the time partial derivative of the iteratively updated level set function, the time partial derivative formula is as follows:

step 10: judging whether an iteration termination condition is reached: if the current iteration time J reaches the maximum iteration time J or the difference value of the iteration results of two adjacent times of the energy functional E (phi) is smaller than a convergence threshold epsilon which is set in advance, executing the step 11, otherwise, adding 1 to the iteration time, and returning to the step 9;

step 11: according to the current updated image level set function phi^jAs shown in fig. 9, the constructed image, i.e., the segmentation result of the image, fig. 10 is the gold standard of the image.

Finally, it should be noted that: the above embodiments are only used to illustrate the technical solution of the present invention, and not to limit the same; while the invention has been described in detail and with reference to the foregoing embodiments, it will be understood by those skilled in the art that: the technical solutions described in the foregoing embodiments may still be modified, or some or all of the technical features may be equivalently replaced; such modifications and substitutions do not depart from the spirit of the corresponding technical solutions and scope of the present invention as defined in the appended claims.

Claims (2)

1. A fundus blood vessel image segmentation system based on a multi-scale level set is characterized in that: the system comprises an input module, a preprocessing module, a scale construction module, a level set module and an output module;

the input module is used for receiving the fundus blood vessel image input by the user and then transmitting the fundus blood vessel image to the preprocessing module;

the preprocessing module is used for expanding the received fundus blood vessel image, and the expanded result is transmitted to the scale construction module;

the scale construction module is used for finding the optimal blood vessel response scale for the blood vessel image at the bottom of the eye, further constructing a Gaussian curve based on the scale, subtracting the original image from the convolution result of the Gaussian curve and the image to obtain an image M, and transmitting the image M to the level set module;

the level set module is used for carrying out blood vessel segmentation on the image M by adopting a level set algorithm and transmitting a segmentation result to the output module;

and the output module is used for outputting and displaying the segmentation result image.

2. A fundus blood vessel image segmentation method based on a multi-scale level set, which is realized by the fundus blood vessel image segmentation system based on the multi-scale level set of claim 1, and is characterized in that: the method comprises the following steps:

step 1: the input module receives an image A of fundus blood vessels to be segmented, wherein the width of the image A is W, and the height of the image A is H;

step 2: the preprocessing module constructs the received fundus blood vessel image A to be segmented to generate a fundus blood vessel image I without an edge background, and the resolution of the fundus blood vessel image I is the same as that of the image A; converting the image A into a gray image J, setting all areas of the image J as omega, directly obtaining an image P if the area to be segmented is full of omega, carrying out mirror image expansion on an effective boundary area in the image to cover a black filling area in the whole image if the area to be segmented is not full of omega, setting the expanded image as P, and executing the step 2.1 if I is equal to P; otherwise, executing step 3;

step 2.1: determining a threshold Q to divide the grayscale image J into a foreground R and a background B, where R and B are binary images, R is 1 and B is 0 in the foreground region, R is 0 and B is 1 in the background region;

step 2.2: constructing a foreground region boundary image R₂Etching the foreground R to obtain R₁，R₂＝R-R₁Wherein the image R₂The image is a binary image;

step 2.3: initializing a variable P to 1, wherein P is an initializing abscissa value of a pixel point in the outer loop;

step 2.4: if p is less than or equal to W and the initialized longitudinal coordinate value q of the pixel point in the outer loop is 1, executing the step 2.5, otherwise, executing the step 3;

step 2.5: if q is less than or equal to H, executing step 2.6, otherwise, making p equal to p +1, and executing step 2.4;

step 2.6: if the current point (p, q) satisfies B ═ 1, step 2.7 is performed, otherwise, let q ═ q +1, step 2.5 is performed;

step 2.7: initializing a variable m which is 1 and is an initial abscissa value of a pixel point in the image inner loop, and initializing

Wherein r is the distance between any two pixels in the initialized image, and the initial coordinates (a, b) are two random positive integers;

step 2.8: if m is less than or equal to W, the variable n is 1, n is the initial value of the vertical coordinate of the pixel point in the inner loop, the step 2.9 is executed, otherwise, the step 2.13 is executed;

step 2.9: if n is less than or equal to H, executing step 2.10, otherwise, making m equal to m +1, and executing step 2.8;

step 2.10: if the point (m, n) satisfies R₂Step 2.11 is performed when n is equal to n +1, otherwise step 2.9 is performed;

step 2.11: calculating a distance dist ((p, q), (m, n)) between the point (p, q) and the point (m, n), wherein dist is a distance metric function, and if r < dist ((p, q), (m, n)), then r ═ dist ((p, q), (m, n)), and let a ═ m, b ═ n;

step 2.12: step 2.9 is skipped to by making n equal to n + 1;

step 2.13: the point (p, q) takes the point (a, b) as the center, R is a symmetrical distance, a corresponding pixel point L in R is obtained, and the value of the pixel point L is assigned to the current point Z; wherein L is_x＝2a-p，L_y2b-q, wherein L_xIs the abscissa L of the pixel point L_yIs the ordinate of the pixel point L;

step 2.14: q is q +1, and the step 2.5 is skipped;

and step 3: for parameter minimum scale sigma_minMaximum scale σ_maxA constant coefficient gamma, a constant coefficient c, a scale step length delta sigma are set, wherein sigma is_min＜σ_max；

And 4, step 4: for the initial scale sigma_iIs set so that σ_i＝σ_min；

And 5: constructing a Hessian matrix of the image I, constructing a Gaussian function through the optimal response scale sigma, and subtracting the image I from the result of convolution operation of the image I and the Gaussian function to obtain an image M, wherein the method comprises the following specific steps:

step 5.1: initialization variable p 1, scale variable k 1, variable v_max＝0, wherein v_maxIs the maximum response value variable of the pixel;

step 5.2: if p is less than or equal to W, the variable q is 1, and step 5.3 is executed, otherwise, step 6 is executed;

step 5.3: if q is less than or equal to H, executing step 5.4, otherwise, making p equal to p +1, and executing step 5.2;

step 5.4: judging the current scale sigma_iWhether or not greater than the maximum scale sigma_maxIf it is greater than σ_maxIf so, let σ_i＝σ_minExecuting step 5.12, otherwise, executing step 5.5;

step 5.5: constructing a Hessian matrix for point (p, q):

wherein ,

D＝∫∫I(p-u,q-v)×A(u,v)dudv

E＝∫∫I(p-u,q-v)×B(u,v)dudv

F＝∫∫I(p-u,q-v)×C(u,v)dudv

a (u, v), B (u, v) and C (u, v) are based on the scale σ_iConstructing a second partial derivative of Gaussian

Wherein D, E, F is the convolution of image I with the corresponding gaussian second order partial derivative A, B, C;

step 5.6: calculating two eigenvalues lambda of the Hessian matrix obtained in step 5.5₁、λ₂, wherein

wherein

Step 5.7: calculated at the current scale sigma_iResponse value of pixel (p, q)

wherein

in the formula ,R_BIs the degree of anisotropy of the linear shape, S is the norm;

step 5.8: if it is

Then order

And κ ═ σ -_i；

Step 5.9: constructing a Gaussian kernel based on the scale k to obtain a Gaussian curve g_κ；

Step 5.10: the Gaussian curve g obtained in the step 5.9_κPerforming convolution operation with image ISubtracting the image I from the convolution result to obtain an image M;

M(p,q)＝(g_κ*I)(p,q)-I(p,q)

step 5.11: updating the scale sigma_i＝σ_i+ Δ σ, Δ σ is the scale step, and step 5.4 is skipped;

step 5.12: updating the coordinate q to be q +1, and jumping to the step 5.3;

step 6, setting an initial zero level set segmentation curve on the image M, and generating a high one-dimensional level set function phi according to the set zero level set⁰；

The set initial zero level set segmentation curve is a closed curve with a specified shape and is a randomly generated closed curve, and the level set function phi is a high one-dimensional symbol distance function generated according to the set initial level set segmentation curve;

step 7, establishing an energy functional E (phi) of level set segmentation, and setting a weighting coefficient α in the level set method₁Weighting factor α₂Curve smoothing term coefficient v and area smoothing term coefficient mu;

the energy functional is as follows:

E(φ)＝E_D+vL(φ)+μP(φ)

wherein ,E_DFor the image data item, L (φ) is a length constraint term, P (φ) is an area constraint term, and H (φ (P, q)) is a Heaviside function of φ (P, q);

E_D(φ)＝α₁∫_Ω|M(p,q)-c₁|²H(φ(p,q))dpdq+α₂∫_Ω|M(p,q)-c₂|²(1-H(φ(p,q)))dpdq

P(φ)＝∫H(φ(p,q))dpdq

wherein c₁Is the mean value of the gray levels in the zero level set segmentation curve, c₂Is the gray average outside the zero level set segmentation curve;

and 8: setting a maximum iteration number J and a convergence threshold epsilon of an energy functional E (phi) of an image level set;

and step 9: separately updating the image level set function set H (phi)^j)、1-H(φ^j)、C^j、φ^jObtaining an energy functional E of the j iteration^j；

Specifically, the method comprises the following steps:

where at is the step of time in which,

for the time partial derivative of the iteratively updated level set function, the time partial derivative formula is as follows:

step 10: judging whether an iteration termination condition is reached: if the current iteration time J reaches the maximum iteration time J or the difference value of the iteration results of two adjacent times of the energy functional E (phi) is smaller than a convergence threshold epsilon which is set in advance, executing the step 11, otherwise, adding 1 to the iteration time, and returning to the step 9;

step 11: according to the current updated image level set function phi^jAn image, i.e. the segmentation result of the image, is constructed.

Images