An Overview of the Topological Gradient Approach in Image Processing: Advantages and Inconveniences

Abstract

Image analysis by topological gradient approach is a technique based upon the historic application of the topological asymptotic expansion to crack localization problem from boundary measurements. This paper aims at reviewing this methodology through various applications in image processing; in particular image restoration with edge detection, classification and segmentation problems for both grey level and color images is presented in this work. The numerical experiments show the efficiency of the topological gradient approach for modelling and solving different image analysis problems. However, the topological gradient approach presents a major drawback: the identified edges are not connected and then the results obtained particularly for the segmentation problem can be degraded. To overcome this inconvenience, we propose an alternative solution by combining the topological gradient approach with the watershed technique. The numerical results obtained using the coupled method are very interesting.

1. Introduction

The goal of topological optimization is to find the optimal decomposition of a given domain in two parts: the optimal design and its complementary. Similarly in image processing, the goal is to split an image in several parts, in particular, in image restoration the detection of edges makes this operation straightforward.

The goal of this work is to show how it is possible to solve some image analysis problems using topological optimization tools. Let us recall that the basic and first idea is based on the work of Amstutz et al. [1] in which the authors propose a solution for crack detection problem using the topological gradient approach. Then, this approach was extended to different problems in image processing like restoration, classification, segmentation, inpainting, and tomography application [2–6]. This paper presents an overview of the applications of the topological gradient approach to major image processing problems with an illustration of both inconveniences and advantages of this technique in image analysis.

Let 𝑣 be a given noisy image defined in a domain Ω⊂ℝ² and 𝑢 the restored image. We recall that a classical way to restore the image 𝑢 is to solve the following partial differential equations (PDE) problems()𝜕−div𝑐∇𝑢+𝑢=𝑣inΩ,_𝑛𝑢=0on𝜕Ω, (1.1) where 𝑐 is a small positive constant and 𝑛 denotes the outward unit normal to 𝜕Ω. This method is well known to give poor results: it blurs important structures like edges and in order to improve this method, nonlinear isotropic and anisotropic methods were introduced. In our topological optimization approach, 𝑐 takes only two values: 𝑐₀ in the smooth part of the image and a small value 𝜌 on edges. Then classical nonlinear diffusive approaches could be seen as a relaxation of our method that increases the instability of the restoration process. The basic principle of the topological gradient approach is presented briefly in Section 2. The other sections are presented as follows The restoration problem with edge detection is presented in Section 3. Section 4 deals with the classification problem. In Section 5, we present the application of the topological gradient approach to tomography problem. Section 6 is devoted to the segmentation problem; we present in this section the major drawback of the topological gradient approach and propose an alternative to this inconvenience. We end this paper with some concluding remarks. Several numerical tests are also presented for all the applications studied.

2. Topological Gradient Approach: Basic Principle

Let Ω be an open bounded domain of ℝ², and let, 𝑗(Ω)=𝐽(𝑢_Ω) be a cost function to be minimized, where 𝑢_Ω is the solution to a PDE problem defined in Ω. For a small 𝜌≥0, let Ω_𝜌=Ω⧵𝜎_𝜌 be the perturbed domain by the insertion of a crack 𝜎_𝜌=𝑥₀+𝜌𝜎(𝑛), where 𝑥₀∈Ω, 𝜎(𝑛) is a straight crack and 𝑛 a unit vector normal to the crack. The topological sensitivity theory provides an asymptotic expansion of 𝑗 when 𝜌 tends to zero. It takes the general form𝑗Ω_𝜌−𝑗(Ω)=𝑓𝜌𝐺𝑥₀+∘𝑓𝜌, (2.1) where 𝑓(𝜌) is an explicit positive function going to zero with 𝜌 and 𝐺(𝑥₀) is called the topological gradient at point 𝑥₀. It is clear that to minimize the criterion 𝑗, we have to insert small holes at points where 𝐺 is negative. Using this gradient type information, it is possible to build fast algorithms. In most applications, a satisfying approximation of the optimal solution is reached at the first iteration of the optimization process. A topological sensitivity framework allowing to obtain such an expansion for general cost functions has been proposed in the work of Masmoudi [1, 7].

We recall that a classical way to restore an image 𝑢 from its noisy version 𝑣 defined in a domain ℝ² is to solve the following PDE problem()𝜕−div𝑐∇𝑢+𝑢=𝑣inΩ,_𝑛𝑢=0on𝜕Ω, (2.2) where 𝑐 is a small positive constant, called the conductivity, 𝜕_𝑛 denotes the normal derivative and 𝑛 is the outward unit normal to 𝜕Ω. The idea of the topological gradient approach is to consider only two values for the conductivity 𝑐: 𝑐₀ in the smooth part of the image and 0 (or rather undefined) in its complementary (the edges). For this reason, classical nonlinear diffusive approaches, where c takes all the values of ]0,𝑐₀[, could be seen as a relaxation of our method. By enlarging the set of admissible solutions, relaxation increases the instability of the restoration process, and this could explain why our method is so efficient with respect to the algorithm cost.

3. Application to Image Restoration

The goal of this section is to use the topological gradient approach as a tool for detecting edges for image restoration [5, 6]. In fact, an image can be viewed as a piecewise smooth function and edges can be considered as a set of singularities. The basic idea consists then in applying the topological gradient approach used for crack detection [1]. More precisely, let Ω_𝜌 be the perturbed domain obtained by the insertion of a crack 𝜎_𝜌. The perturbed solution 𝑢_𝜌∈𝐻¹(Ω_𝜌) satisfies−div𝑐∇𝑢_𝜌+𝑢_𝜌=𝑣inΩ_𝜌,𝜕_𝑛𝑢_𝜌=0on𝜕Ω_𝜌. (3.1) The variational formulation of problem (3.1) is given byfind𝑢_𝜌∈𝐻¹Ω_𝜌suchthat𝑎_𝜌𝑢_𝜌,𝑤=𝑙_𝜌(𝑤)∀𝑤∈𝐻¹Ω_𝜌, (3.2) where 𝑎_𝜌 is the bilinear form, defined on 𝐻¹(Ω_𝜌)×𝐻¹(Ω_𝜌) by𝑎_𝜌(𝑢,𝑤)=_Ω𝜌(𝑐∇𝑢∇𝑤+𝑢𝑤)𝑑𝑥, (3.3) and 𝑙_𝜌 is the linear form defined on 𝐿²(Ω_𝜌) by𝑙_𝜌(𝑤)=_Ω𝜌𝑣𝑤𝑑𝑥. (3.4)

Edge detection is equivalent to look for a subdomain of Ω where the energy is small. So, our goal is to minimize the energy norm outside edges𝑗𝜌=𝐽_𝜌𝑢_𝜌=_Ω𝜌∇𝑢_𝜌². (3.5) The variation of the cost function 𝑗(𝜌)=𝐽_𝜌(𝑢_𝜌) induced by the insertion of the small crack 𝜎_𝜌 is given by the topological gradient theory [1]:𝑗𝜌−𝑗(0)=𝜌²𝐺𝑥₀,𝑛+∘𝜌², (3.6) with 𝐺𝑥₀,𝑛=−𝑐𝜋∇𝑢₀𝑥₀⋅𝑛∇𝑣₀𝑥₀⋅𝑛−𝜋∇𝑢₀𝑥₀⋅𝑛², (3.7) where 𝑢₀ is the solution of the direct problem (3.1) such that 𝜌=0 and 𝑣₀ is the solution of the following adjoint problems−div𝑐∇𝑣₀+𝑣₀=−𝜕_𝑢𝐽(𝑢)𝜕inΩ,_𝑛𝑣₀=0on𝜕Ω. (3.8) We note here that the solutions 𝑢₀ and 𝑣₀ are calculated in the initial domain without cracks, in the sense that the crack size is equal to zero. Then, the topological gradient can be rewritten in the following ways𝐺(𝑥,𝑛)=⟨𝑀(𝑥)𝑛,𝑛⟩, (3.9) where 𝑀(𝑥) is the symmetric matrix given by 𝑀(𝑥)=−𝑐𝜋∇𝑢₀(𝑥)∇𝑣₀(𝑥)^𝑇+∇𝑣₀(𝑥)∇𝑢₀(𝑥)^𝑇2−𝜋∇𝑢₀(𝑥)∇𝑢₀(𝑥)^𝑇. (3.10) We can deduce that 𝐺(𝑥,𝑛) is minimal when the normal 𝑛 is the eigenvector associated to the smallest (i.e., most negative) eigenvalue of the matrix 𝑀(𝑥). In the following, this eigenvalue will be considered as the topological gradient. Our algorithm consists in inserting small heterogeneities in regions where the topological gradient is smaller than a given threshold 𝛼<0. These regions are the edges 𝜔_𝜌 of the image. The algorithm is as follows.

Restoration Algorithm
(i)Initialization: 𝑐=𝑐₀.(ii)Calculation of 𝑢₀ and 𝑣₀ the solutions of the direct (3.1) and adjoint (3.8) problems respectively.(iii)Computation of the 2×2 matrix 𝑀 and its lowest eigenvalue 𝜆_min at each point of the domain.(iv)Set 𝑐₁=𝜌if𝑥∈Ω,𝜆_min𝑐<𝛼<0,𝜌>0₀elsewhere. (3.11)(v)Calculation of 𝑢₁ solution to problem (3.1) with 𝑐=𝑐₁.

For numerical experiments, we considered, in a first step, grey level images (see Figure 1). We note that the perturbed image is obtained by adding a gaussian noise (𝜎=20) to the original image, which grey level takes its values in the interval [0, 255]. We present in Figure 1 images restored by different approaches. More precisely, the restored image using the topological gradient approach is presented and compared to both linear and nonlinear diffusion methods [8]. Notice that the convergence of the nonlinear approach is achieved after 53 iterations. Our method can be interpreted as a linear isotropic diffusion method defined in the complementary of edges. In the anisotropic case, not only the amount of diffusion is adapted locally to the data but also the direction of smoothing. Like in the anisotropic nonlinear method, we smooth along image edges while inhibiting smoothing across edges. This can be achieved by replacing the diffusivity scalar by a diffusivity tensor. To show the efficiency of our method, we test an anisotropic algorithm with the topological gradient as follows

(a)

(b)

(c)

(d)

Figure 1 

(a) Initial 512×512 Lena image, (b) restored image using the linear diffusion method, (c) restored image using a nonlinear diffusion method, and (d) restored image using the topological gradient method.

Anisotropic Topological Gradient Algorithm
(i)Initialization: 𝑐=𝑐₀.(ii)Calculation of 𝑢₀ and 𝑣₀: solutions of the direct (3.1) and adjoint (3.8) problems. (iii)Computation of the 2×2 matrix 𝑀 and its lowest eigenvalue 𝜆_min at each point of the domain. (iv)Set 𝑐=𝐴𝑐_𝜏𝑛𝐴^⋆, with 𝐴=(𝜏,𝑛) and 𝑐_𝜏𝑛=(^𝑐00_0𝜀), where 𝜏 and 𝑛 represent the tangent and normal directions.(v)Calculation of 𝑢₁ solution to problem (3.1).

In Figure 2, we show that there is no significant difference between the results obtained using the isotropic or the anisotropic approach with the topological gradient method. Our two algorithms are tested on several images. All numerical tests show that the results obtained using the two algorithms are very similar. We show in Figure 2 a zoom of the 512×512 Lena image and its restoration by the nonlinear method presented in [8], and the isotropic and the anisotropic topological gradient algorithms.

(a)

(b)

(c)

(d)

Figure 2 

(a) zoom of the noisy Lena image (𝜎=20), (b) zoom of the restored Lena image using the nonlinear diffusion method, (c) zoom of the restored Lena image using topological gradient approach with the isotropic algorithm, and (d) zoom of the restored Lena image using topological gradient approach with the anisotropic algorithm.

We consider now color images such that the RGB model is used. Then, the image 𝐼 is defined by a 3D vector 𝐼(𝑥)∈ℝ³ which represents the intensity of the three colors: red, green, blue, and since each component can be treated as a grey level image, then the topological gradient approach can easily be used, provided the cost function to be minimized has no coupling terms between the several channels. For the restoration problem, one usually set the cost function to the square norm of the gradient of the image. If we consider the most simple case, that is, the standard 𝐿² norm of the gradient, then the cost function can be easily written as a sum of the three independent terms, one for each channel. Then, the color image problem can be decomposed into three channel-by-channel independent problems. We mention here that the three channels can be supposed to be coupled, and a tensor gradient of the image can be considered as it is proposed in the work of Di Zenzo [9].

Some other coupling techniques have been defined in order to process multichannel images. We can cite for instance the Beltrami diffusion-type process, in which images are considered as manifolds. A nonlinear structure tensor is defined and locally adjusted according to a gradient-type measure and allows one to control the diffusion process [10–12].

For the numerical tests, we first consider in Figure 3 the natural extension of the topological gradient method to multi-channel images, using the standard gradient norm and we compare this approach with the Ambrosio-Tortorelli version of the Mumford-Shah functional.

(a)

(b)

(c)

(d)

Figure 3 

(a) original image, (b) noisy image (SNR=10.59), (c) restored image by the classical topological gradient method (SNR=19.44, CPU=61.03), and (d) restored image using Ambrosio-Tortorelli algorithm (SNR=18.47, CPU=76.77).

Then to improve the restoration process, we propose in Figure 4 to present numerical results using the Di Zenzo approach coupled with the topological gradient approach. Figure 4 shows the restoration algorithm applied to the Di Zenzo gradient using a 500×380 color image. The noisy image is obtained with an additive gaussian noise, such that the signal-to-noise ratio (SNR) is equal to 10.54. The third image corresponds to the restored image by our algorithm. The SNR of the restored image is 19.23. The identified edges given in the fourth image, correspond to the pixels where the most negative eigenvalue of the matrix 𝑀 is smaller than a given threshold.

(a)

(b)

(c)

(d)

Figure 4 

(a) the initial 500×380 image, (b) noisy image, (c) the restored image, and (d) edges of the restored image.

We must note here that, for the case corresponding to the classical topological gradient method, the solution 𝑢 of (2.2) can be written as 𝑢=(𝑢¹,𝑢²,𝑢³), where 𝑢^𝑘 are solutions of the following problems:−div𝑐∇𝑢^𝑘+𝑢^𝑘=𝑣^𝑘𝜕inΩ,_𝑛𝑢^𝑘=0on𝜕Ω, (3.12) for 𝑘=1,2,3, which implies that the cost function 𝑗 can also be decomposed into three terms, one for each channel, as follows:𝑗𝜌=³_𝑘=1𝑗^𝑘𝜌, (3.13) with𝑗^𝑘𝜌=_Ω𝜌∇𝑢^𝑘_𝜌²𝑑𝑥, (3.14) and the topological gradient 𝐺(𝑥) can be rewritten in the following way:𝐺(𝑥)=⟨𝑀(𝑥)𝑛,𝑛⟩, (3.15) where 𝑀(𝑥) is the symmetric 2×2 matrix given by𝑀(𝑥)=³_𝑘=1−𝜋𝑐∇𝑢^𝑘₀(𝑥)∇𝑣^𝑘₀(𝑥)^𝑇+∇𝑣^𝑘₀(𝑥)∇𝑢^𝑘₀(𝑥)^𝑇2−𝜋∇𝑢^𝑘₀(𝑥)∇𝑢^𝑘₀(𝑥)^𝑇, (3.16) with 𝑢^𝑘₀ and 𝑣^𝑘₀ being the solutions of the direct and adjoint problems with respect to each channel 𝑘.

On the other hand, for the case corresponding to the topological gradient method using the Di Zenzo gradient, a multispectral tensor associated to the image vector field is considered and the largest eigenvalue of the tensor corresponds to the norm of the gradient, called the Di Zenzo gradient [9]. The authors of [13] have proved the following result.

Lemma 3.1. For 𝑢=(𝑢¹,𝑢²,𝑢³) in 𝐿²(Ω,ℝ³), the Di Zenzo gradient is given by ‖∇𝑢‖_DZ=∇𝑢¹²+∇𝑢²²+∇𝑢³²1+1−4𝑓(∇𝑢)2, (3.17) where 𝑓 is the function defined by 𝑓(∇𝑢)=det²∇𝑢¹,∇𝑢²+det²∇𝑢¹,∇𝑢³+det²∇𝑢²,∇𝑢³∇𝑢¹²+∇𝑢²²+∇𝑢³²², (3.18) with det²∇𝑢^𝑠,∇𝑢^𝑡=𝜕𝑢^𝑠𝜕𝑥₁𝜕𝑢^𝑡𝜕𝑥₂−𝜕𝑢^𝑡𝜕𝑥₁𝜕𝑢^𝑠𝜕𝑥₂²for𝑠,𝑡=1,2,3, (3.19) and |⋅| is the standard norm in ℝ². As in the previous section, the topological gradient can be written as 𝐺(𝑥,𝑛)=⟨𝑀(𝑥)𝑛,𝑛⟩, (3.20) where 𝑀(𝑥) is the 2×2 symmetric matrix given by 𝑀(𝑥)=³_𝑘=1−𝜋𝑐∇𝑢^𝑘₀(𝑥)∇𝑣^𝑘₀(𝑥)^𝑇+∇𝑣^𝑘₀(𝑥)∇𝑢^𝑘₀(𝑥)^𝑇2−𝜋𝐻∇𝑢₀(𝑥)∇𝑢^𝑘₀(𝑥)∇𝑢^𝑘₀(𝑥)^𝑇, (3.21) where 𝐻 is given by 𝐻∇𝑢₀=1+1−4𝑓∇𝑢₀2, (3.22) and 𝑓 is defined by (3.18).

4. Application to Image Classification

This section is concerned with the problem of classifying an original image using 𝑛 predefined classes 𝐶_𝑖, 1≤𝑖≤𝑛, by choosing the color intensity as a classifier. More precisely, let 𝑢₀ be the original image defined on an open set Ω of ℝ². Our goal is to find 𝑛 color 𝐶_𝑖, 1≤𝑖≤𝑛, and a partition Ω_𝑖, 1≤𝑖≤𝑛, of Ω such that 𝑢₀ is close to 𝐶_𝑖 in Ω_𝑖. The classified image 𝑢 will then be defined by𝑢(𝑥)=𝐶_𝑖∀𝑥∈Ω_𝑖, (4.1) where {Ω_𝑖}_{𝑖=1,…,𝑛} are defined byΩ_𝑖=𝑥∈Ω;𝑥belongstothe𝑖thclass. (4.2) We recall that the variational approach consists in defining a cost function measuring the root mean square difference between the original image and the classified image𝐽(Ω_𝑖)_{𝑖=1,…,𝑛}=^𝑛_𝑖=1Ω𝑖𝑢₀(𝑥)−𝐶_𝑖²𝑑𝑥. (4.3) The minimization of 𝐽 is easy, because for each point 𝑥∈Ω, we only have to find 𝑖_𝑥=argmin{|𝑢₀(𝑥)−𝐶_𝑖|;𝑖=1,…,𝑛} and add 𝑥 to subset Ω_𝑖𝑥. This can be called the closest class algorithm because each pixel of the original image is assigned in the classified image to its closest class. Inspired by the work of Samson et al. [14] in which the authors propose a classification model coupled with a restoration process, we propose to use the topological gradient approach applied to image restoration problem for both supervised and unsupervised classification problems [3, 4]. Our idea is to still consider the PDE problem given by (2.2) and instead of considering 𝑐=0 (or 𝑐=𝜌 from a numerical point of view) on the edge set and 𝑐=𝑐₀ elsewhere, we consider 𝑐=𝜌 on the contours of the image and 𝑐=1/𝜌 elsewhere. As 𝜌 is supposed to be small, then 𝑢 and 𝑢₀ are almost the same. Otherwise, 𝑐=1/𝜌, which implies that the PDE defined in (2.2) is nearly equivalent to ▵𝑢=0 with Neumann boundary condition (normal derivative) equal to zero, and this will provide a really smooth image. The classification problem is then as follows.

Supervised Classification Algorithm
(i)Initialization: 𝑐=𝑐₀.(ii)Calculation of 𝑢₀ and 𝑣₀ the solutions of the direct (3.1) and adjoint (3.8) problems respectively.(iii) Computation of the 2×2 matrix 𝑀 and its lowest eigenvalue 𝜆_min at each point of the domain. (iv)Set 𝑐₁=𝜌if𝑥∈Ω,𝜆_min𝑐<𝛼<0,𝜌>0₀𝜌elsewhere. (4.4)(v)Calculation of 𝑢₁ solution to problem (3.1) with 𝑐=𝑐₁.(vi) Application of the k-mean classification algorithm to 𝑢₁.

As in the restoration process, the classification algorithm is applied to color images by decomposing the image in the RGB color space and then we deal separately with the three component images. The topological gradient of the sum can be considered as the sum of the three topological gradients with respect to each channel. Figure 5 shows the classification of a 500×380 color image, using 30 color levels (whereas the original image has 256³ different colors), without regularization in a first time, then with regularization. We should mention here that the control of the regularization process can be made by considering progressively different values of 𝑐∶𝑐=𝑐₀/𝜌,  𝑐=𝑐₀/𝜌²,  𝑐=𝑐₀/𝜌³, and so on.

(a)

(b)

Figure 5 

From (a) The unregularized classified image and and (b)the regularized classified image.

On the other hand, the unsupervised classification corresponds to a classification problem in which the classes are not given. In this case, it is possible to determine them in an optimal way, still by using the topological gradient method. The idea is to study the impact of changing the value of a class 𝐶_𝑖∶=𝐶_𝑖+1 or 𝐶_𝑖−1 on the cost function which implies that for each class, we add (or substract, depending on which variation provides the most negative variation) one to the value of the class. This algorithm has been applied in the previous numerical tests in order to determine the optimal value of the classes: in Figure 5, the values of   3,2, and 5 classes with respect to the channels RGB have been determined with this algorithm. One may notice that this algorithm needs at least that the number of classes is known. But in unsupervised classification, it is not the case. The idea is then to add another term in the cost function defined in (4.3), measuring the number of classes. This can be seen as a regularization term, because we usually do not want to get too many classes. We propose here to rewrite the cost function given by (4.3) as follows:𝐽Ω_{𝑖𝑖=1,…,𝑛}=^𝑛_𝑖=1Ω𝑖𝑢₀(𝑥)−𝐶_𝑖²𝑑𝑥+𝛼𝑛, (4.5) where 𝛼 is a positive regularization coefficient. It is then possible to bound 𝑛 as the two terms of the cost function are positive. The unsupervised classification algorithm is as follows

Unsupervised Classification Algorithm
(i)Initialization with 𝑛=1. (ii)Determination of the optimal value.(iii) While 𝑛 is smaller than its upper bound: (iv) set 𝑛∶=𝑛+1 and determine the optimal values of the new 𝑛 classes.(v)Stop if the value of the cost function at the optimum (with respect to 𝐶_𝑖) is larger than the previous optimal value with 𝑛−1 classes.

This algorithm has been applied in Figure 5 in order to find the optimal number of classes and their optimal values.

5. Application to Tomography Problem

In this section, we use the topological gradient as a tool for the reconstruction problem in tomography. We recall that a standard approach for regularizing the ill-posed problem of tomographic imaging consists in the following optimization problem:min_𝑓𝐴𝑓−𝑔²+𝜆𝜙𝑓, (5.1) where 𝐴 denotes a system matrix defining the discrete Radon transform, 𝑔 is the measured projection, 𝜙 is a regularization functional, and 𝜆 represents a parameter which controls the tradeoff between a good fit to the data and the smoothness of the solution. Inspired by the work of. Chan et al. [15], we propose to consider the following minimization problem:min_𝑓Ω𝐴𝑓−𝑔²𝑑𝑥+𝑐_Ω∇𝑓^𝛼𝑑𝑥, (5.2) where 𝑐 is a positive constant and 𝛼 is equal to 1 on the edges and to 2 elsewhere. Using the same idea of the topological gradient approach, we consider, for a given regular function 𝑔, the following direct problem:−div𝑐∇𝑓_𝜌+𝐴^#𝐴𝑓_𝜌=𝐴^#𝑔inΩ_𝜌,𝜕_𝑛𝑓_𝜌=0on𝜕Ω_𝜌, (5.3) where 𝑛 denotes the outward unit normal to 𝜕Ω_𝜌 and 𝐴 and 𝐴^# are respectively the discrete Radon transform and the dual discrete Radon transform.

Let 𝑣 be the solution of the adjoint problem−div(𝑐∇𝑣)+𝐴^#𝐴𝑣=−𝜕_𝑓𝐽𝑓𝜕inΩ,_𝑛𝑣=0on𝜕Ω. (5.4) By considering the same algorithm applied for image restoration, we compute the reconstruction solution by solving the following proble:−div𝑐₁𝑥₀∇𝑓+𝐴^#𝐴𝑓=𝐴^#𝜕𝑔inΩ,_𝑛𝑓=0on𝜕Ω, (5.5) with𝑐₁𝑥₀=1∇𝑓𝑥₀,𝑥₀∈𝑥∈Ω,𝜆_min<𝛼₀,𝑐<0₀elsewhere, (5.6) where 𝑐₀ is a given constant. In order to obtain the optimal value of the coefficient 𝑐₀, we have studied the behavior of the SNR of the reconstructed image as a function of 𝑐₀. In fact, we have experimentally confirmed that the quality of the reconstructed image is first increasing rapidly with 𝑐₀ and then decreasing slightly, as we can see it in Figure 6. The optimal value is 𝑐₀=30, and it has then been considered for our numerical tests. We have tested our method on the well known Shepp-Logan head phantom, for which projections data have been computed using a discrete Radon transform. Figure 7(a) shows the original image. A gaussian noise is added to the data with a signal to noise ratio SNR=22, the sinogram is represented in Figure 7(b). Figure 7(c) shows the reconstruction result (PSNR=26.18). We also present in this paper a numerical comparison using several methods. Table 1 shows the PSNR (peak signal to noise ratio, in dB), SNR (signal-to-noise ratio, in dB), SSIM (structural SIMilarity), and MSE (mean squared error) of the reconstructed images for several methods: FBP (filtered backprojection), FBP + Hamming (with a Hamming filter on the FBP image), and the topological gradient approach using either the TV (total variation) norm or both TV and 𝐿² norms. The sinogram has a noise level such that the SNR is 24.5. We can see in Table 1 that the best results for all indicators are obtained with the hybrid 𝐿¹/𝐿² norm within the topological gradient approach. The peak SNR is higher of more than 3 decibels, the structural similarity index is larger than 94% and the mean squared error is smaller than 1%. Note that both topological gradient methods provide much better results than the filtered backprojection, which is not surprising.

Table 1 

Comparison between several reconstruction methods, quantified by several ratios and similarity indexes, for a noisy sinogram with SNR=24.5.

Figure 6 

Evolution of the PSNR (in dB) of the reconstructed image by the topological gradient method, as a function of the level of the regularization coefficient 𝑐0.

(a) Original image

(b) Noised sinogram (SNR=22)

(c) Reconstructed image (PSNR=26.18)

Figure 7 

Reconstruction of the Shepp-Logan head phantom using the topological gradient method.

6. Application to Image Segmentation

We propose in this part a new algorithm for the segmentation problem based on the topological gradient approach. Through the different and various applications which were presented, we can easily conclude that the topological gradient approach is a very interesting method in image analysis. First, the topological gradient method can easily be applied to many problems in image analysis. Second, we should mention that we obtain excellent numerical results and the computing time is very interesting and confirms the theoretical complexity 𝑂(𝑛⋅log(𝑛)) of our algorithms, where 𝑛 is the size of the image. However, the topological gradient approach presents a major drawback: the identified edges are not connected and then the results obtained for the segmentation problem can be degraded, particularly for complex images. So, the main idea in this part is to take advantage of the topological gradient efficiency, to detect the main contours with an interesting computational cost (the topological gradient algorithms require only three system resolutions) and to find an alternative solution for the discontinuous edges detected by the topological gradient approach. The approach proposed in this paper is based on mathematical morphology [16, 17] and more precisely on the watershed transformation [18]. A major advantage of the method proposed is that it attenuates the oversegmentation usually observed with the watershed technique. This attenuation is due to the fact that the topological gradient provides a global analysis of the image then the almost unwanted contours due to the noise added to a given image are effectively reduced by our approach. The basic idea observed using this coupled method is the use of the topological gradient instead of the morphological gradient classically used with watersheds. The property of the topological gradient to be less sensitive than the morphological or classical gradients to noise and small variations of the image, will contribute significantly to give the best identification of the edges and then to attenuate the oversegmentation usually obtained using watershed techniques. We illustrate in this work the efficiency of the combined method proposed for the segmentation problem. Other numerical tests and comparisons with other methods can be found in [19]. In this part, we consider two real images which have various difficulties (curves, circles, straight lines, etc.). The first image is a 302×280 fruit-basket gray level image and the second one is a 399×246 gray level image which represents a road scene. Figure 8 shows the segmentation result of the two previous images, according to the coupled method for image segmentation based on a watershed algorithm [18] and the topological gradient approach previously presented in this paper. The results shown in Figures 8(c) and 8(d) are obtained with 𝑐=1, we can improve the attenuation of the identified regions number after the segmentation process by choosing larger values of the coefficient 𝑐, but this will depend on what we want to segment and what we want to obtain at the end of the segmentation process: if one would obtain a maximum of details then it suffices to take 𝑐=1, otherwise, one should consider larger values of 𝑐 which gives a considerable attenuation of the segmented regions. Figures 8(e) and 8(f) show the numerical results of the segmentation process with 𝑐=50. We can clearly observe the attenuation of the obtained regions such that the main edges of the two images are conserved. One should remark that we can either define the watershed of the function 𝑢 or of its gradient: the difference between the two definitions is that in the first case we obtain the influence zones of the processed image, while the second case gives the image edges. It is well-known that in both cases, the watershed gives an oversegmentation, and to avoid this drawback, the alternative to combine it with the topological gradient approach seems to be a good choice to overcome this inconvenience as shown in Figure 8.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 8 

Segmentation using a watershed transformation combined with the topological gradient: (a) is a 302×280 fruit-basket gray level image, (b) is a 399×246 road scene image, (c) is the segmented fruit-basket image with 𝑐=1, (𝑒) is the same segmented image such that 𝑐=50, (d) is the segmented road scene image with 𝑐=1, and (f) represents the same segmented image such that 𝑐=50.

Finally, due to a large number of approaches used for the segmentation problem, it clearly appears that it is important to compare our experimental results with methods already proposed in the literature. Particularly, we propose to compare our method with an active contour model based on the level set approach, which is well known to be an efficient method, extensively used in many applications over the last decade. The basic idea in active contour models is to evolve a curve subject to constraints from a given image, in order to detect different objects in that image. To achieve this goal, we start with a curve around the object to be detected, the curve moves toward its interior normal and has to stop on the boundary of the object. This was the first idea of classical snakes and active contour models proposed by [20]. Our numerical tests are based on a level set method proposed by [21] and given by the following evolution equation: 𝜕𝑢𝜕𝑡=𝑔(|∇𝐼|)div∇𝑢||||+()∇𝑢+𝛼∇𝑢∇𝑔,∇𝑢,in0,∞×Ω, (6.1) with boundary and initial conditions given by 𝜕𝑢/𝜕𝑛=0 on (0,∞)×𝜕Ω and Φ(0,𝑥,𝑦)=Φ₀(𝑥,𝑦) in Ω, 𝛼≥0 is a positive given coefficient, 𝑔(|∇𝐼|) is an edge detector function defined in our numerical tests by 𝑔(𝑠)=1/(1+𝑠²), and 𝜙₀ represents the initial level set function. We refer the reader to [8] for more details about geodesic active contours and level set methods. Figure 9 shows the numerical results obtained. We should mention that in order to obtain the final contour illustrated in Figure 9(e), the algorithm requires 500 iterations and a computational time of 78 seconds. We present in Figure 9 some solutions as time evolves: during the evolution, the initial curve is shrinking and stopping as soon as it is close to an object boundary and spilling in order to detect the others objects.

(a)

(b)

(c)

(d)

(e)

(f)

Figure 9 

Segmentation results of a synthetic image using an active contour model: different iterations are displayed from (b) to (e), (f) is the segmented image using our new approach.

7. Conclusion

We presented in this paper many applications in image processing of the topological gradient approach with an illustration of both advantages and inconveniences of this technique. We proved that the topological gradient approach provides an excellent frame for solving different image processing problems, same if an alternative way to overcome its major drawback to give not connected contours is to combine it with other approaches. It has been applied to image restoration, edge detection, image classification, and image segmentation for both grey level and color images. It has also been applied to tomography problem, and our main goal consists in applying it to other real life problems mainly that in all the applications proposed, we obtained good results and the computing time is very short. We propose in a forthcoming paper to perform the application of the topological gradient approach to image segmentation, eventually by combining it with other approaches and particularly using marker criteria. In fact, as the edges are detected in regions where the topological gradient is the most negative, then it suffices to extract some points which belong to the edge set and then use these points as a selected marker set. We also intend to extend this work to color image segmentation and three-dimensional segmentation.