Edge Extraction with an Anisotropic Vector Field using Divergence Map

Giuliani Donatella
International Journal of Image Processing (IJIP) Volume (6), Issue: (4) : 2012 255
Edge Extraction With an Anisotropic Vector Field
Using Divergence Map
Giuliani Donatella giulianidonatella@libero.it
Scientifc-Didactic Polo of Rimini
University of Bologna
Via Angherà 22, Rimini, Italy
Abstract
The aim of this work is edge detection by a deformable contour procedure, using
an external force field derived from an anisotropic flow, with different external and
initial conditions. By evaluating the divergence of the force field, we have
generated a divergence map associated with it in order to analyze the field
convergence. As we know, divergence measures the intensity of convergence or
divergence of a vector field at a given point, so by means level curves of the
divergence contour map, we have automatically selected an initial contour for the
deformation process. The initial curve must include areas from which the vector
field diverges pushing it towards edges. Furthermore the divergence map brings
out the presence of curves pointing to the most significant geometric parts of
boundaries corresponding to high curvature values, in this way it will result better
defined geometrical shape of the extracted object.
Keywords: Edge Extraction, Active Contours, Anisotropic Flow, GGVF
1. INTRODUCTION
Image segmentation and boundary extraction are diffusely handled topics of
research in image processing. These problems have been dealt with in various
forms [1], [2] [3], [4],[5],[6],[7] nevertheless the initial conditions and the ability to
change the topology of the evolving curve are the main limitations for most of
these methods. Edge detection may be realized by a deformable contour process
with active contours or surfaces embedded within an image domain that move
under the influence of internal and external forces [1],[6]. Internal forces, which
are defined within the curve or surface itself, are designed to smooth the model
during deformation. External forces, which are computed from image data, are
used to move the model toward researched features in the image [8],[9],[10]. In
this paper the external force field is derived from an anisotropic flow, in a more
general framework respect to the GGVF [10],[11], in order to come to an overall
view of the field flow, analyzed through its divergence values. A careful overview
of divergence map, showing the divergence of the external force field, gives us the
opportunity to circumscribe areas within which we can place an initial contour that
may be subsequently deformed.

Giuliani Donatella
2. PRELIMINARIES
In the traditional parametric models an active contour or snake, expressed
explicitly by parametric equations [ ]1,0,))(),(()( ∈= ssysxsx
r
, is defined within a
given image I(x,y) and subjected to modifications under the action of forces, until
the evolving curve fits well into the final contour [1],[7]. The final shape of the
contour to be extracted will be such as to minimize an energy functional
associated with it, so given:
( ) ( ) ( )∫ += dssxEsxEsxE ExtInt )()()(
rrr
(1)
the first term is the internal energy that expresses a priori knowledge of the model
in relation to the degree of flexibility of an active contour:
( ) ( ) ( )∫ 







⋅+⋅= ds
ds
xd
s
ds
xd
ssxEInt
2
2
22
2
1
)(
rr
r
βα (2)
the term α(s) controls the contour tension, while β(s) regularises its rigidity. The
second term ( ))(sxEExt
r
represents the external energy related to a potential
energy function ( ))(sxP
r
deriving from the image I(x,y), whose local minima
correspond to edges of the features to be extracted. By using a variational
approach [12], the contour that minimizes the total energy must satisfy the Euler-
Lagrange equation:
02
2
2
2
=∇+





+





− P
ds
xd
ds
d
ds
xd
ds
d
rr
βα (3)
where∇ is the gradient operator. Through energy minimization we derive a
resolution of a static problem. We could construct a deformable model able to
create a geometrical shape that evolves over time, by introducing a time-variable
parametric equation )),(),,((),( tsytsxtsx =
r
. To this end, indicating with µ(s) and γ(s)
the density of mass and the damping coefficient respectively, instead of equation
(3) we obtain:
P
s
x
ss
x
st
x
t
x
−∇=





∂
∂
∂
∂
+





∂
∂
∂
∂
−
∂
∂
+
∂
∂
2
2
2
2
2
2 rrrr
βαγµ (4)
The equilibrium is reached when internal and external forces are equal, this
implies achieving a steady state in which time derivatives will become null.
When the solution ),( tsx
r
of equation (4) stabilizes, we carry out the result of
equation (3). Neglecting the inertial term and thus the second order derivatives,
considering dumping, elasticity and rigidity as constant functions, the numerical
solution of equation (4), can be reduced to that one of the dynamic equation:

Giuliani Donatella
( ) ( ) ( ) ( )
( ) ( )




=
+
∂
∂
−
∂
∂
=
∂
∂
sxsx
xF
s
tsx
s
tsx
t
tsx
Ext
0
4
4
2
2
0,
,,,
rr
rrrrr
βαγ (5)
where ( )sx0
r
is an initial contour. In order to reduce the sensitivity of this model to
the initial conditions, edge extraction may be realized using a different class of
external forces, the GGVF force field or Generalized Gradient Vector Flow
[10],[11],[13],[14] obtained by solving a diffusion problem. In the GGVF contour
generation, the external force field will be referred as:
( )xvF GGVF
Ext
rrr
= (6)
Using the calculus of variations [12], the GGVF force field ),( yxv
r
can be found by
solving the following diffusion equation:
( ) ( ) ( )fvfhvfgvt ∇−⋅∇−∇⋅∇=
rrr 2
(7)
where 2
∇ is the Laplacian operator, ),( yxf∇ is the gradient of the edge map
f(x,y), derived from the gradient of brightness function I(x,y), that can be computed
as ( )yxIyxf ,(),( ∇= or using any other edge detector. The edge map gradient
are vectors directed towards boundaries to detect with norms significantly different
from zero in proximity of them. In equation (7) ( )fg ∇ and ( )fh ∇ are space
varying weighting functions, being dependent on absolute gradient of edge map,
generally not uniform. The function ( )fg ∇ will be monotonically non-increasing,
since the vector field ),( yxv
r
will be weakly variable far from edges to be extracted
where image intensities are uniform. On the other hand, ( )fh ∇ should be
monotonically non-decreasing, therefore, when f∇ is large, the vector field
),( yxv
r
should have a trend nearly equal to f∇ .
The main shortcomings to be overcome in the generation of deformable contours
are: 1) the initialization problem, i.e. the excessive sensitivity to shape and initial
position of an initial curve 2) the weak convergence of models towards edges,
especially in regions with highly variable concavities 3) the capture range, i.e. the
size of area inside which an active contour can be initialized to be able achieving
the desired boundary. We would like to explore a method that tries to overcome
one of the main drawbacks for most of the approaches that were introduced until
now: the initialization problem. The convergence of arbitrary initial contours for the
GGVF leads only partially to the expected results, giving seemingly
incomprehensible outcomes, as we can see in Fig.1a, Fig1b and Fig.2a, Fig2b for
two test images. This paper is concerned with the problem of identifying
automatically [15],[16],[17] an initial contour using the divergence of a force field
derived from an anisotropic vector flow.

Giuliani Donatella
FIGURE 1a: Initialization problem FIGURE 1b:
Initialization problem
FIGURE 2a: Initialization problem FIGURE 2b: Initialization problem
3. THE DIVERGENCE MAP OF AN ANISOTROPIC VECTOR
FIELD
In this section we suggest a more general analysis of the diffusion process for an
external force field. To this aim, we could see the solution of the GGVF equation
(7) as the steady state of the following generalized parabolic equation:
( )( ) ( )
( ) ( )
)8(
,0,, 0


=
+∇⋅∇=
yxvyxv
vFvfgdivvt
rr
rrrr
where div is the divergence operator of the diffusion component and ( )vF
rr
is a
term generating an external force for the diffusion process of ( )yxv ,
r
, g(.) is the
conduction coefficient that must be null or tend to zero at boundaries. So, in
adiabatic condition hypothesis the vector flow is realized inside or outside the
region. If the scalar-valued diffusivity function g(.) is monotonically non-increasing
the diffusion has been stopping across the edges to be extracted. We obtain the
GGVF field [11] if ( ) 0=∇⋅∇∇ vfg
r
( ) ( ) ( )fvfhvF ∇−⋅∇−=
rrr
( ) k
f
efg
∇
−
=∇ ( ) ( )fgfh ∇−=∇ 1

Giuliani Donatella
with initial conditions fv ∇=0
r
and k as a constant positive value. If we consider an
anisotropic flow for the vector field ( )yxv ,
r
expressed as follows:
( )( )
( )
)9(
0,,


∇=
∇⋅∇=
fyxv
vfgdivvt
r
rr
where ( ) 0=vF
rr
, we carry out results rather similar to those of the GGVF field,
given that the term ( ) ( ) ( )fvfhvF ∇−⋅∇−=
rrr
tends to zero whether near edges, for the
chosen initial conditions, or far from them, since the function ( )fh ∇ is getting
irrelevant, consequently the overall contribution of ( )vF
rr
is not significant. But now
considering the problem through the parabolic equation (9), we could interpret the
process as a field flow from boundaries toward inside or outside without crossing
edges and pointing to them because initially equal to the edge map gradient
(Fig.3b). We could note from Fig.3a that the vector field can capture object
boundaries from either sides. As we know, divergence is a measure of a field
convergence or divergence at a given point by means signed scalar values. So if
we compute the divergence of ),( yxv
r
, we would obtain negative values
corresponding to object boundaries towards which the vector field converges,
sinks of the flow, whereas positive values would define regions from which the
field flow springs out, the flow sources.
FIGURE 3a: The vector field ),( yxv
r
FIGURE 3b: The vector field
),( yxv
r

Giuliani Donatella
As a consequence, the analysis of divergence of the external force field ),( yxv
r
may be very useful to feature detection, since it allows to demarcate regions from
which the flow has been originated. Furthermore in correspondence of long and
narrow boundary indentations, the field convergence towards opposite very near
sides gives rise to zones of high positive divergence (Fig.3a,3b), hence
divergence values may well highlight them. An efficient numerical scheme to
resolve the anisotropic diffusion equation (9) may be that one proposed by
Perona-Malik [18] and realized recurring to the 4-neighbors discretization of the
Laplacian operator:
[ ]t
WW
t
EE
t
SS
t
NN
tt
jijijijijiji
vgvgvgvgvv ,,,,,,
1 rvrrrr
∇+∇+∇+∇+=+
λ
jijijiW
jijijiE
jijijiS
jijijiN
vvv
vvv
vvv
vvv
,1,,
,1,1,
,,1,
,,1,
rrr
rrr
rrr
rrr
−=∇
−=∇
−=∇
−=∇
−
++
+
−
(10)
with λ=
4
1 , and values of the edge-stopping function g(.) evaluated through a
gradient approximation of f∇ along horizontal and vertical directions, such as:
( )
( )
( )
( )jijiW
jijiE
jijiS
jijiN
ffgg
ffgg
ffgg
ffgg
,1,
,1,
,,1
,,1
−=
−=
−=
−=
−
+
+
−
The anisotropic vector field so obtained, from now on AVF, is partially shown in
Fig.4a and Fig.4b for different numbers of time iterations, it results very similar to
that one evaluated with the GGVF approach already seen. At first we could note
that the size and the shape of capture range vary with the number of iterations
performed in equation (9).

Giuliani Donatella
FIGURE 4a: AVF field from equation (9) FIGURE 4b: AVF field from equation
(9)
after 20 iterations after 140 iterations
If now we consider equation (8) with an external force ( )vF
rr
significantly different
from zero, the convergence of ),( yxv
r
to the steady state will result faster than the
GGVF solution. By using the diffusion equation:
( )( ) ( )
( )


∇=
+∇⋅∇=
fyxv
vFvfgdivvt
0,,
)11(
r
rrrr
with, for example, an the external force ( ) ( ) ( )fvfgvF ∇−⋅∇=
rrr
we will be able to
speed up the convergence process, because ( )vF
rr
is approximately null near
edges but increases as moving away from them. This implies an extension of
capture range for an equal number of time iterations. In Fig.5a and Fig.5b we
could compare the AVF field derived from equation (11) with GGVF field obtained
through equation (7). With the same number of iterations it has been realized an
increase of capture range, the region inside of which must be placed any initial
contour. Consequently the evaluation of the AVF vector field will be less time
consuming respect to the GGVF because we could choose a reduced number of
iterations to achieve similar results. In Fig.5c and Fig.5d we have compared the
results obtained by AVF field with those of GGVF using the same number of
iterations (in the case in point 80) and the same initial polygonal. Recurring to 120
iterations to research the solution of equation (7), the GGVF field has been able to
detect edges completely.
FIGURE 5a: AVF field from equation (11) FIGURE 5b: GGVF field from
equation (7)

Giuliani Donatella
FIGURE 5c: AVF field from equation (11) FIGURE 5d: GGVF field from equation (7)
with the same initial contour after 80 iterations
As an alternative approach we could evaluate the vector field ),( yxv
r
using
diffusion equation (11) but recurring to different external forces and initial
conditions, so we have:
( )( ) ( )
( )


=
+∇⋅∇=
00,,
)12(
yxv
vFvfgdivvt
r
rrrr
where now ( )vF
rr
is an external attractive force directed towards edges, for
example of the type ( ) ( ) ( )fvfhvF ∇−⋅∇=
rrr
that stops acting when the vector field,
initially null everywhere, tends to the edge map gradient. Again, by a comparison
between the AVF field generated with equation (12) and the GGVF field after the
same number of iterations, we could verify an acceleration of the convergence
process (Fig.6a,6b). Therefore the area, within which the field is able to act on a
deformable initial contour pushing it towards edges, is larger with the AVF field of
equation (12) than GGVF. As a consequence, we need a smaller number of
iterations to obtain a vector field as effective as the GGVF.

Giuliani Donatella
FIGURE 6a: AVF field from equation (12) FIGURE 6b: GGVF field from equation (7)
After these considerations, we may suggest a method to generate automatically
an initial curve, that could be able to adequately reconstruct the final contour.
Indeed the effects due to initialization problem and capture range extension are
reduced using the divergence map of the normalized force field, defined as:
( ) ( )
( ) 







=
tyxv
tyxv
divtyxI D
,,
,,
,, r
r
In the following figures it has been shown divergence maps corresponding to
vector fields resulted as numerical solutions of equation (9) at different times t, i.e.
for different numbers of iterations, referred to the first test image ( Fig.7a,7b,7c,7d)
and to the second one (Fig.8a,8b,8c,8d) respectively:
FIGURE 7a: Divergence Map of AVF FIGURE 7b: Divergence Map of AVF
from equation (9) after 0 iterations from equation (9) after 30 iterations

Giuliani Donatella
FIGURE 7c: Divergence Map of AVF FIGURE 7d: Divergence Map of
AVF
FIGURE 8a: Divergence Map of AVF FIGURE 8b: Divergence Map of AVF
FIGURE 8c: Divergence Map of AVF FIGURE 8b: Divergence Map of AVF
By an accurate analysis of the shown figures, we may point out that each
divergence map is characterized by a gray background with divergence values
near zero, dark curves with negative divergence corresponding to edges, towards
which the vector field converges, and a system of light curves with positive values,
defining regions from which the vector field comes out. These regions are varying
with the number of iterations after which the flow process has been interrupted.

Giuliani Donatella
The zones between opposite values of divergence define the capture range size.
Moreover, by increasing the iteration number, the sides of areas that delimit parts
of image inside which the field is null, collapse each other, especially in
conjunction with deep and narrow concavities (see Fig.8b and Fig.8c). Thereby
these curves, corresponding to positive divergence values, remain as traces of
those edge sections with high curvatures and concur in forming the skeleton of the
figure, as we can see in Fig.7d after 240 iterations for the solution of equation (9)
or in Fig.9 only after 120 iterations for the solution of equation (11). In Fig.9 and
Fig.10 we may compare divergence maps for AVF field evaluated by equation
(11) and GGVF field respectively.
FIGURE 9: Divergence Maps of AVF field FIGURE 10: Divergence Maps of GGVF
from equation (11) for different iteration numbers for different iteration numbers
As a consequence shape and position of any initial contour also depends on the
number of iterations performed in the numerical resolution. We have to take into
account boundaries of areas with positive divergence, because within them the
vector field is null, so the sections of an evolving curve remain entrapped in their
interiors, as we can see in Fig.11a, where an arbitrary initial contour is
superimposed to the divergence map of our image.

Giuliani Donatella
FIGURE 11a: Divergence Map of GGVF field with an arbitrary initial contour after 80
iterations
If we now review the results shown in Fig.5c and Fig.5d realized using AVF field
and GGVF field respectively, after 80 iterations and we superimpose the initial
polygonal with the corresponding divergence maps, we are able to understand the
anomalous shape of the extracted edges of Fig.5d, apparently incomprehensible.
Some portions of the initial contour fall within the areas where the GGVF field is
null, thereby those traits remain trapped inside them (see Fig.11c). Instead after
120 iterations, the capture range of the GGVF field is become more extended, so
the evolving curve converges to right results
FIGURE 11b: Divergence Map of AVF field FIGURE 11c: Divergence Map of GGVF
field

Giuliani Donatella
FIGURE 11d: Divergence Map and active contour with GGVF field after 120 iterations
We would remember that divergence will be zero either when the field is null or
when the inward and the outward flux are equal. In the divergence map relative to
the second test image, composed by two disjoint boundaries, we may well see
two closed curves corresponding to positive divergence values. Within the first
curve the field is zero, instead inside the second, divergence values are near zero
because the flux is null, whereas the vector field is rather uniform and changes its
direction converging to the inner circle, as we can see in Fig.12a. Then the traits
of an initial contour that fall inside this region are pushed towards the inner circle
rather than the outer edge (Fig.12b,112c).
FIGURE 12a: AVF field after 80 iterations

Giuliani Donatella
FIGURE 12b: Divergence Map of AVF field FIGURE 12c: Divergence Map of AVF field
with an arbitrary initial contour with an arbitrary initial contour
Thereby it results much easier to understand the anomalous behaviour in the
convergence of each polygonal sections of Fig.11a, Fig.11b, Fig.11c or Fig.12b,
Fig.12c using maps of divergence associated with the GGVF or AVF field. The
traits of an initial contour are pushed along the field lines that, according to the
Divergence Theorem, start at sources, i.e. areas of positive divergence, and end
at sinks with negative divergence. In addition we would like to make some
observations about the capture range extension. Looking closely at the
divergence map, we could note that the area enclosed between the external or
internal curve with positive divergence, and the curve of negative divergence,
corresponding to the object edge, delimits the capture range. As we may well see
in Fig.7a, 7b, 7c, 7d and in Fig.8a, 8b, 8c, 8d, shape and size of this area are
varying with the number of iterations chosen. So we may conclude that
divergence map, at any given time, is strictly related to the initialization problem
and the capture range, because it allow us to clearly distinguish between sources
and sinks for our image, enabling us to properly place an initial curve, as shown in
Fig.13. Accordingly we suggest to generate a contour map of the field divergence
in order to visualize it through level sets, i.e. the set of points in the image domain
where divergence is constant. In this way we could automatically extract an
appropriate initial contour using level curves corresponding to high divergence
values.

Giuliani Donatella
FIGURE13: Divergence Map with initial contours converging towards edges to detect
Given an image I(x,y), it has been carried out ),( yxv
r
resolving numerically
equation (11). After evaluating its divergence map CRDID →⊂ 2
: in the gray
colour space C = [0;255], it has been generated a contour map concerning to it
(Fig.14, Fig.15) and subsequently it has been selected as initial contour the level
curve of intensity 153 for the two test images. In this way the chosen curve
certainly encloses source areas from which the field diverges, pushing the active
contour towards object edges. As we can see from Fig.16 to Fig.19, edges have
been extracted correctly, establishing a quite good effectiveness of the proposed
approach.
FIGURE 14: Level curves of Divergence Map FIGURE15: Level curves of Divergence Map
FIGURE.16: Initial contour automatically selected FIGURE.17: Initial contour
automatically selected

Giuliani Donatella
FIGURE. 18: Initial contour automatically selected FIGURE.19: Initial contour
automatically selected
We have used real images for the validation of the proposed framework and
promising experimental results were obtained with Gray Matter images derived
from a data set of patients affected by Alzheimer disease (Fig.20) and subjects of
a control group. This model doesn’t require the user to input manually an initial
curve close to the edges, that may be fastidious to do for an application dealing
with a large number of images.
FIGURE 20: Initial contour automatically selected for a GM image after 10 Iterations
4. CONCLUSIONS
In this work, we have analyzed edge extraction by a deformable contour
procedure, using an external force field derived from an anisotropic vector flow,
with different diffusion equations and initial conditions. We have suggested a new

Giuliani Donatella
approach to overcome some limitations of traditional parametric snakes, specially
regarding to the initialization problem. The external force field, resolution of an
anisotropic diffusion problem, should be considered in relation to its divergence
map, because the source areas, from which the vector field diverges pushing the
deformable initial contour towards the edges to detect, are well defined and
recognizable in it. We have observed that only if source areas are completely
included inside the initial contour, the edge detection will turn out well. Moreover
we have proposed as automatic initialization procedure the use of level curves of
divergence contour map, in such a way it doesn’t be required a specification of an
initial curve by the user.
ACKNOWLEDGMENT
I sincerely would like to thank Prof. Naldi Giovanni (Department of Mathematics
“Federigo Enriques”, University of Milan – Italy) for the helpful contribution given
for this study.
REFERENCES
[1] M. Kass, A. Witkin, D. Terzopoulos, “ Snakes: active contour models” , Int.
Jour. Comp. Vision, Vol. I, N.4, 1988, pp 321-331
[2] T. Mc Inerney, D. Terzopoulos, “T-snakes:Topology adaptive snakes” ,
Medical Image Analysis,4,2000,pp. 73-91.
[3] J. Liang,T.Mc Inerney, D. Terzopoulos, “United Snakes” , Medical Image
Analysis,10,2006pp. 215-233
[4] S.J. Osher, J.A. Sethian, “Fronts propagating with curvature dependent
speed: algorithms based on Hamilton-Jacobi formulations”, J. of
Computational Physics, 79, 1988, pp 12-49
[5] S.J. Osher, R.P. Fedkiw, “Level set methods: an overview and some recent
results”, J. of Computational Physics, 169, 2001, pp 463-502
[6] C. Xu, L.J. Prince, “Snakes ,shapes, and gradient vector flow”, IEEE Trans.
on Image Processing, 7, 1998
[7] C. Xu, A. Yezzi, and J. Prince, “On the Relationship between Parametric and
Geometric Active Contours,” Proc. Asilomar Conf. Signals, Systems, and
Computers, pp. 483-489, 2000.
[8] B. Li, T.A Scott, “Active contour external force using vector field convolution
for image segmentation”,. IEEE Transactions on Image Processing, 16,
2007, pp.2096-2106
[9] A. Kovács, T. Szirányi, “Improved force field for vector field convolution
method”,
ICIP IEEE 2011, pp. 2909-2912.

Giuliani Donatella
[10] C. Xu, L.J. Prince, “Gradient Vector Flow: a New External Force for Snakes”,
IEEE Proc. Conf. on Comp. Vis. Patt. Recogn., CVPR’97
[11] C. Xu, L.J. Prince, “Generalized gradient vector flow external forces for
active contours”, Signal Proc., an Intern. Journal, 71, 1998
[12] R. Courant,D. Hilbert, “Methods of Mathematical Physics”, Vol I, Ed.
Interscience, New York,1953
[13] J. Cheng, W.F. Say, “Dynamic Directional Gradient Vector Flow for Snakes”,
IEEE Trans. on Image Processing, 15, 2006
[14] N. Paragios, O. Mellina-Gottardo, V. Ramesh, “Gradient Vector Flow Fast
Geometric Active Contours”, IEEE Trans. on Pattern Analysis and Machine
Intelligence, 26, 2004
[15] C. Li, J. Li, M.D. Fox, “Segmentation of Edge Preserving Gradient Vector
Flow: an Approach Toward Automatically Initializing and Splitting of
Snakes”, IEEE Proc. Conf. on Comp. Vis. Patt. Recogn., CVPR’05
[16] C. Tauber, H. Batiata, A. Ayache, “A Robust Active Initialization and Gradient
Vector Flow for Ultrasound Image Segmentation”, MVA2005 IAPR
Conference on Machine Vision Applications, 2005
[17] F. Yabin, L. Caixia, Z. Bingsen, P. Zhenkuan, “An Improved Algorithm of
Contour Initialization in Active Contour Model”,Image and Graphics, ICIG
2007
[18] P. Perona, J. Malik, “Scale-space and edge detection using anisotropic
diffusion” IEEE Trans. on Pattern Analysis and Machine Intelligence, 12,
1990

Edge Extraction with an Anisotropic Vector Field using Divergence Map

Recommended

More Related Content

What's hot (19)

Similar to Edge Extraction with an Anisotropic Vector Field using Divergence Map (20)

Recently uploaded (20)

Edge Extraction with an Anisotropic Vector Field using Divergence Map