Coherence-Enhancing Diffusion Filtering.
The completion of interrupted lines or the enhancement of flow-like structures is a challenging task in computer vision, human vision, and image processing. We address this problem by presenting a multiscale method in which a nonlinear diffusion filter is steered by the so-called interest operator (second-moment matrix, structure tensor). An m-dimensional formulation of this method is analysed with respect to its well-posedness and scale-space properties. An efficient scheme is presented which uses a stabilization by a semi-implicit additive operator splitting (AOS), and the scale-space behaviour of this method is illustrated by applying it to both 2-D and 3-D images.
Introduction
Oriented flow-like structures arise in many computer vision and image processing problems: Within the
field of texture analysis they appear for instance in the automatic grading of fabrics or wood surfaces,
but they are also of importance for fingerprint analysis in forensic applications. They are present in
many scientific imaging problems ranging from fluid dynamics to meteorology, and last but not least
in medical imaging, for instance in the analysis of trabecular structures in bones. Interestingly, related
tasks such as gap completion or the completion of interrupted lines also play a role in human vision
tasks such as perceptual grouping. Moreover, the success of paintings by Munch or van Gogh suggests
that emphasizing flow-like structures may create effects which fascinate many people.
Analysing flow-like patterns is an active research topic for certainly more than one decade; see for
instance Kass and Witkin [28]. Besides methods such as Gabor filters or steerable filters, many of
these approaches are equivalent to the so-called structure tensor (interest operator, second moment
[17, 44]. In the meantime this field has even entered textbooks [21], and it can be regarded as
well-established.
Much less attention has been payed to the question how to enhance flow-like patterns. Poor quality
of fingerprint or trabecular bone images is not unusual. In those cases it would be desirable to have
a tool which improves the quality of flow-like structures without destroying for instance semantically
important singularities like the minutiae in fingerprints.
For problems like the grading of fabrics or applications to fluid dynamics it is also useful to have a
multiscale simplification of the original image by embedding it into a scale-space in order to obtain
a subsequently coarser, more global impression of the main flow-like structures. Of course, such a
scale-space should take into account the coherence of the structures by smoothing mainly along their
preferred orientation instead of perpendicular to it.
Since flow-like structures can also be present in higher dimensional data sets, e.g. 3-D images of
trabecular bones, it should be possible to generalize such a method to arbitrary dimensions.
The preceding problems will be addressed in this paper by presenting an m-dimensional scale-space for
the enhancement of coherent structures. The underlying concept can be motivated from ideas of Perona
Current address: Department of Computer Science, University of Copenhagen, Universitetsparken 1, DK-2100
Copenhagen, Denmark.
and Malik for improving edge detection by creating a feedback loop of an edge detector and a nonlinear
diffusion process [42]. In a similar way we embed a classical method for describing flow-like structures
- namely the structure tensor approach - into a nonlinear diffusion process. This turns a method for
analysing coherent pattern into a technique for enhancing and simplifying them. In contrast to most
nonlinear diffusion filters, however, we use an approach where the process is steered by a diffusion
tensor instead of a scalar-valued diffusivity. This enables true anisotropic behaviour by adapting the
diffusion process not only to the location, but allowing also different smoothing in different directions.
We shall see that this filter belongs to a class of nonlinear diffusion methods for which many well-posedness
and scale-space properties can be proved. This is also in contrast to the Perona-Malik filter,
whose continuous formulation is only locally well-posed [29]. Details of the numerical implementation
will be discussed as well. In particular, an efficient novel stabilization technique based on a semi-implicit
additive operator splitting (AOS) is presented. It extends previous encouraging experiments with AOS
schemes for nonlinear diffusion filtering with a scalar-valued diffusivity to the fully anisotropic case
with a diffusion tensor. All theoretical and numerical discussions hold in the m-dimensional case. As
examples, results from 2-D and 3-D implementations are presented.
The outline of the papers is as follows: Section 2 surveys the underlying structure tensor method for
describing coherence in images. This method is used in Section 3 for constructing a nonlinear diffusion
process which smoothes along coherent flow-like structures. This process is called coherence-enhancing
diffusion (CED). Section 4 gives a detailed analysis of the theoretical properties of a more general class
of diffusion filters comprising CED. In Section 5 numerical and algorithmical questions are addressed;
in particular, a more efficient alternative to the explicit scheme is presented. Section 6 applies CED
to 2-D and 3-D images from different application areas. The paper is concluded with a summary in
Section 7. Some results in this paper have been presented earlier at conferences [50, 56].
Related work. Since the work of Perona and Malik [42] numerous nonlinear diffusion filters have
been proposed; see e.g. [23, 53] for an overview. Nevertheless, most of them use a (spatially varying)
scalar diffusivity, not a diffusion tensor. Thus, they act inhomogeneously (nonuniformly) on the image,
but - in our terminology - they remain isotropic.
True anisotropic diffusion filtering is studied in the reaction-diffusion model of Cottet and Germain
[15]. Its diffusion tensor uses the eigenvectors v 1 k ru oe , its eigenvalues are given by
(j ? 0):
This choice is similar to our method it that sense that it diffuses mainly along strongly anisotropic
structures. However, there are two important differences: Firstly, we observe that this diffusion tensor
cannot be treated within a scale-space framework using uniformly positive definite diffusion tensors.
Moreover, the Cottet-Germain model uses an additional reaction term which leads to nontrivial
steady-states and qualifies it as a pure restoration method without scale-space ambitions. Secondly,
the eigendirections of D are adapted to ru oe , not to the eigendirections of the structure tensor. We
shall see that the introduction of an integration scale in the structure tensor is an essential model
feature in order to obtain reliable orientation estimates for flow-like structures. Recently Cottet and
El Ayyadi have proposed a modified image restoration model which replaces the spatial regularization
by a temporal regularization [14].
Nitzberg and Shiota [39] pioneered shape-adapted Gaussian smoothing, where the shape of an anisotropic
Gaussian convolution kernel is a function of the structure tensor. Later on similar proposals
have been made in [33, 61] and supplemented with scale-selection methods [1]. It should be noted that
a common feature of all the preceding shape-adapted Gaussian smoothing methods is the combination
of isotropic smoothing inside a region with anisotropic smoothing along edges. Space-variant Gaussian
blurring is, however, not equivalent to an inhomogeneous diffusion process, and it does not preserve the
average grey value. Anisotropic diffusion filters with isotropic diffusion within regions and anisotropic
diffusion along edges can be found in [51] and the references therein. These methods are different from
the present approach: CED is basically a pure 1-D diffusion, where a minimal amount of isotropic
smoothing is added only for regularization purposes.
Other anisotropic partial differential equations for smoothing images rely on morphological methods
such as the mean-curvature motion (geometric heat equation, Euclidean shortening flow) [30, 4]
with - being the direction perpendicular to ru.
Since mean-curvature motion propagates isophotes in inner normal direction with curvature-dependent
speed, we should not expect such completely local methods to be capable of closing interrupted line-like
structures. The same holds for the affine invariant version of this process [3, 46]. All these methods do
not take into account semilocal information as can be gained from the structure tensor. Recent results
by Carmona and Zhong [12] confirm the importance of semilocal estimates of the smoothing direction
also for equations of mean-curvature type, if one is interested in specific goals such as enhancement
of line-like structures.
Three-dimensional nonlinear diffusion filters have been investigated first by Gerig et al. [18] in the
isotropic case, and by Rambaux and Gar-con [43] in the edge-enhancing anisotropic case. A generalization
of coherence-enhancing anisotropic diffusion to higher dimensions was first proposed in [56]. A
recent three-dimensional PDE-based filter by Krissian et al. [32] and a 3-D reaction-diffusion process
by Payot et al. [41] may be related to these anisotropic diffusion techniques.
The discussion above shows that the distinctive features of our approach is a semi-local analysis by
means of the structure tensor combined with 1-D diffusion along one of its eigenvectors.
2 The structure tensor
First we review a reliable tool for analysing coherent flow-like structures. Consider an m-dimensional
image
domain\Omega := (0; a 1 ) \Theta \Delta \Delta \Delta \Theta (0; am ), and let an image u(x) be represented by a bounded mapping
A very simple structure descriptor is given by ru oe , the gradient of a Gaussian-smoothed version of
u:
The standard deviation oe denotes the noise scale, since it makes the edge detector ignorant of details
smaller than O(oe). Convolution on a finite
domain\Omega can be defined as convolution between K oe and
~
u on IR m , where ~ u denotes an extension of u by mirroring. This results in applying no-flux boundary
conditions.
Although ru oe is useful for detecting edges, it is unsuited for finding parallel structures, as we can see
from

Figure

1. The left image shows an original fingerprint. Figure 1(b),(c) illustrates the gradient
orientation using grey values: vertical gradients are depicted in black, and horizontal ones in white.
We observe that for small oe high fluctuations remain, while larger oe lead to entirely useless results.
This is due to the fact that for larger oe neighbouring gradients with the same orientation, but opposite
sign cancel one another. Gradient smoothing averages directions instead of orientations 2 . To make the
structure descriptor invariant under sign changes, we may replace ru oe by its tensor product
oe\Omega ru oe := ru oe ru T
This matrix is symmetric and positive semidefinite, and its eigenvectors are parallel and orthogonal
to ru oe , respectively. The corresponding eigenvalues jru oe j 2 and 0 describe just the contrast in the
eigendirections. Now that we have replaced directions by orientations, we can average the orientations
by applying a componentwise convolution with a Gaussian
This matrix is named structure tensor, interest operator or second-moment matrix. It is useful for
many different tasks, for instance for analysing flow-like textures [44], corners and T-junctions [17, 39],
shape cues [33] and spatio-temporal image sequences [27]. Equivalent approaches have been discovered
In our terminology, gradients with opposite sign share the same orientation, but point in opposite directions.

Figure

1: Local orientation in a fingerprint image. From Top Left To Bottom
Right: (a) Original
Orientation of smoothed gra-
dient, Orientation of smoothed gradient, 5. (d) Structure tensor
independently in [7, 28]. A book by J-ahne [27] gives a nice overview of these methods and clarifies
their mutual relations.
It is not hard to verify that the symmetric matrix J ae is positive semidefinite. Let its eigenvalues
be ordered such that
and let fw 1 ; :::; wm g denote the corresponding orthonormal set of eigenvectors. Since the eigenvalues
integrate the variation of the grey values within a neighbourhood of size O(ae), they describe the
average contrast in the eigendirections w . Thus, the integration scale ae should reflect the
characteristic size of the texture. Usually, it is large in comparison to the noise scale oe. The eigenvector
wm corresponds to the smallest eigenvalue -m . It is the orientation with the lowest fluctuations, the
so-called coherence orientation 3 .

Figure

1(d) depicts the coherence orientation for the fingerprint image. We observe that it is exactly
the desired average orientation of the lines. It should be noted how well the singularity corresponds
to the singularity in the original fingerprint image.
The eigenvalues of J ae provide useful information on the coherence of a structure, i.e. the actual amount
3 For the case that we have p identical minimal eigenvalues, one may regard every orientation within the span of
as a coherence orientation. For natural images, this situation happens almost never and can be neglected.
of anisotropy. As a measure for the coherence, one can define
It becomes large for strongly differing eigenvalues, and it tends to zero for isotropic structures.
Coherence-enhancing anisotropic diffusion in m dimensions
Now that we know how to analyse coherent structures, we draw our attention to the question of how to
enhance them. This can be done by embedding the structure tensor analysis into a nonlinear diffusion
filter.
The principle of nonlinear diffusion filtering is as follows. One calculates a processed version u(x; t) of
f(x) with a scale parameter t - 0 as the solution of a diffusion equation with f as initial condition
and reflecting boundary conditions:
(D ru)
on\Omega \Theta (0; 1); (8)
@\Omega \Theta (0; 1): (10)
Hereby, n denotes the outer normal and h:; :i the usual inner product. If one wants to adapt the
diffusion process to the image itself one should choose the symmetric positive definite diffusion tensor
as a function of the local image structure. In the case of enhancing flow-like
patterns, it is natural to adapt the diffusion tensor D to the structure tensor J ae (ru oe ).
How should this function look like? For enhancing coherence in images with flow-like structures, we
need a smoothing process which acts mainly along the flow direction wm (with the notations from
Section 2) and the smoothing should increase with the strength of its orientation given by the coherence
-. This may be achieved in the following way:
We construct D such that it has the same eigenvectors as J ae and its eigenvalues are given by
-m :=
ae ff if -=0,
serves as a threshold parameter: For - AE C we get -m - 1, and - C leads to -m - ff.
The exponential function and the small positive parameter ff 2 (0; 1) were introduced mainly for two
theoretical reasons: First, this guarantees that the smoothness of the structure tensor carries over to
the diffusion tensor. The second reason is that the process never stops: Even if the structure becomes
remains some small linear diffusion with diffusivity ff ? 0. Thus, ff serves as
a regularization parameter which keeps the diffusion tensor uniformly positive definite. In the next
section we shall see that these are useful requirements in order to establish many theoretical properties
for CED.
4 A general well-posedness and scale-space framework
Coherence-enhancing anisotropic diffusion filtering can be regarded as a special case within a more
general framework for nonlinear diffusion filtering, for which many well-posedness and scale-space
properties can be found. This framework shall be analysed now.
We study a diffusion filter with the following properties:
Assume that f
Let a:=inf
\Omega f , and consider the problem
\Omega \Theta (0; 1);
@\Omega \Theta (0; 1);
where the diffusion tensor satisfies the following properties:
Smoothness:
positive definiteness:
If w
there exists a positive
lower bound -(K) for the eigenvalues of D(J ae (w)).
Evidently, coherence-enhancing anisotropic diffusion satisfies the preceding requirements. Under the
assumptions (P c ) the following theorem, which generalizes and extends results from [13, 51], can be
proved.
Theorem 1 (Well-posedness properties, smoothness, extremum principle)
The problem (P c ) has a unique solution u(x; t) in the distributional sense, which is smooth for t ? 0:
-\Omega \Theta (0; 1)):
This solution depends continuously on f with respect to
(\Omega\Gamma , and it fulfils the extremum principle
a - u(x;
on\Omega \Theta (0; 1): (13)
Proof: See Appendix A.1.
Remarks:
(a) We observe a strong smoothing effect which is characteristic for many diffusion processes: boundedness
of the initial image is sufficient to obtain an infinitely often differentiable solution for arbitrary
small positive times. Additional requirements - for instance that f should be uniformly
continuous in order to apply the theory of viscosity solutions - are not necessary in our case.
(b) The continuous dependence of the solution on the initial image has significant practical impact
as it ensures stability with respect to perturbations of the original image. This is of importance
when considering stereo image pairs, spatio-temporal image sequences or slices from medical CT
or MRI sequences, since we know that similar images remain similar after filtering.
(c) The extremum principle offers the practical advantage that, if we start for instance with an
image within the range [0; 255], we will never obtain results with grey value such as 257.
(d) Coherence-enhancing anisotropic diffusion in not the only filter which is covered by the theorem
above. The well-posedness theory (as well as the subsequent scale-space framework) comprises
the regularized Perona-Malik filter by Catt'e, Lions, Morel and Coll [13] as well as the edge-enhancing
anisotropic diffusion filter described in [51]. It is also possible to establish a semidiscrete
and fully discrete well-posedness and scale-space theory which proceeds in the same way;
see [54] for more details.
4.1 Scale-space properties
Scale-space representations embed an image f into a family fT t f j t - 0g of gradually smoother,
simplified versions of it. Long before this concept became popular in the western world by the work
of Witkin [59] and Koenderink [31], it had been discovered and axiomatically justified in Japan by
Iijima [26, 55]. Most axiomatic scale-space representations are devoted to the linear diffusion scale-space
(Gaussian scale-space), but one can also create nonlinear scale-spaces which can be classified in
a unique way [3, 8, 40]. All these approaches are uncommitted, since their unique classification does
not give us the liberty to incorporate a-priori knowledge.
Anisotropic nonlinear diffusion filtering, however, does offer the possibility to incorporate knowledge
into its evolution. So let us now study some of its scale-space properties.
It is evident that the filter class (P c ) satisfies typical scale-space properties such as the semi-group
property, invariance under Euclidean image transformations, grey-level shifts, or contrast reversion.
Moreover, since it can be written in divergence form and uses reflecting boundary conditions, it is also
not hard to verify that the average grey level
Z
is not affected by nonlinear diffusion filtering [54, pp. 63-64]:j\Omega j
Z
Average grey level invariance is a property in which diffusion scale-spaces differ from morphological
scale-spaces. In general, the evolution PDEs of the latter ones are not of divergence form and do not
preserve the mean grey value. A constant average grey level is useful for scale-space based segmentation
algorithms such as the hyperstack [37, 49]. It is also a desirable quality for applications in medical
imaging where grey values measure physical qualities of the depicted object.
4.1.1 Causality in terms of nonenhancement of local extrema
Let us now turn to the question in which sense an evolution equation of type (P c ) can be considered
as a smoothing, information-reducing image transformation.
Koenderink [31] required that a scale-space evolution should not create new level curves when increasing
the scale parameter. If this is satisfied, iso-intensity linking through the scales is possible and a
structure at a coarse scale can (in principle) be traced back to the original image (causality). For this
reason, he required that at spatial extrema with nonvanishing determinant of the Hessian isophotes
in scale-space are upwards convex. He showed that this constraint can be written as
A sufficient condition for the causality equation (16) to hold is requiring that local extrema with
positive or negative definite Hessians are not enhanced: an extremum in - at scale ' satisfies @ t u
if - is a minimum, and @ t is a maximum. This implication is easily seen: In the first case,
for instance, the eigenvalues j 1 ,.,j m of Hess(u) are positive. Thus,
which gives just the causality requirement (16).
Nonenhancement of local extrema has first been used by Babaud et al. [5] in the context of linear
diffusion filtering. However, it is also satisfied by nonlinear diffusion scale-spaces, as we shall see now. 4
4 As in the linear diffusion case, nonenhancement of local extrema generally does not imply that their number is
nonincreasing for dimensions
Theorem 2 (Nonenhancement of local extrema).
Let u be the unique solution of
2\Omega is a local extremum
of u(:; ') with nonvanishing Hessian. Then,
is a local maximum, (18)
is a local minimum. (19)
Proof: See Appendix A.2.
Nonenhancement of local extrema is just one possibility to end up with Koenderink's causality re-
quirement. Another way to establish causality is via the extremum principle (13) following Hummel's
reasoning; see [25] for more details.
4.1.2 Lyapunov functionals and behaviour for t !1
Interestingly, causality in terms of nonenhancement of local extrema can be regarded as a first representative
of a much larger class of smoothing properties of nonlinear diffusion filtering. Other representatives
result from studying the Lyapunov functionals of these processes. Lyapunov functionals
are energy-like expressions which decrease during the evolution and which have a lower bound. They
can be useful for finding attractors of a process for t ! 1. The theorem below establishes a class of
Lyapunov functionals and shows that processes of type (P c ) converge to a flat steady-state which is
given by the average grey value of the original image. This is a desirable property for scale-spaces,
since such a steady-state can be regarded as the coarsest, most global representation of the original
image.
Theorem 3 (Lyapunov functionals and behaviour for t !1).
Suppose that u is the solution of (P c ) and let a, b, - be defined as in (P c ) and (14). Then the following
properties are valid:
(a) (Lyapunov functionals)
For all r 2 C 2 [a; b] which are convex on [a; b] the function
Z
is a Lyapunov functional:
(i) It is bounded from below by the value of \Phi for a constant image with the same average grey
level:
continuous in 0, and decreasing for all t ? 0.
(b) (Convergence)
converges to a constant image with the same average grey value:
lim
Proof: See Appendix A.3.
What are interesting representatives of this family of Lyapunov functionals? Considering the Lyapunov
functions associated with respectively, the preceding
theorem gives the following corollary.
Corollary 1 (Special Lyapunov functionals).
Let u be the solution of (P c ) and a and - be defined as in (P c ) and (14). Then the following functions
are decreasing for t 2 [0; 1):
(a) ku(t)k L p for all p - 2.
(b) M 2n [u(t)] :=j\Omega j
Z
(c) H [u(t)] :=
Z
Corollary 1 offers multiple possibilities of how to interpret nonlinear anisotropic diffusion filtering as
a smoothing transformation.
As a special case of (a) it follows that the energy ku(t)k 2
is reduced by diffusion. Using Parseval's
equality we know that a decreasing energy is also equivalent to a decreasing sum of the squared Fourier
coefficients.
Part (b) gives a probabilistic interpretation of anisotropic diffusion filtering. Consider the intensity
in an image f as a random variable Z f with distribution F f (z), i.e. F f (z) is the probability that an
arbitrary grey value Z f of f does not exceed z. By the average grey level invariance, - is equal to the
expected value
Z
z dF u(t) (z); (23)
and it follows that M 2n [u(t)] is just the even central moment
Z
The second central moment (the variance) characterizes the spread of the intensity about its mean. It
is a common tool for constructing measures for the relative smoothness of the intensity distribution.
The fourth moment is frequently used to describe the relative flatness of the grey value distribution.
Higher moments are more difficult to interpret, although they do provide important information for
tasks like texture discrimination [22, pp. 414-415]. All decreasing even moments demonstrate that the
image becomes smoother during diffusion filtering. Hence, local enhancement effects, which object to
increase central moments, are overcompensated by smoothing in other areas.
If we choose another probabilistic model of images, then part (c) characterizes the information-theoretical
side of our scale-space. Provided the initial image f is strictly positive on \Omega\Gamma we may
regard it also as a two-dimensional density. 5 Then,
Z\Omega
is called the entropy of u(t), a measure of uncertainty and missing information [11]. Since anisotropic
diffusion filters increase the entropy, the corresponding scale-space embeds the genuine image f into
a family of subsequently likelier versions of it which contain less information. Moreover, for
the process reaches the state with the lowest possible information, namely a constant image. This
information-reducing property indicates that anisotropic diffusion might be generally useful in the
context of image compression. In particular, it helps to explain the success of nonlinear diffusion
filtering as a preprocessing step for subsampling as observed in [16]. The interpretation of the entropy
in terms of Lyapunov functionals carries also over to generalized entropies, see [48] for more details.
From all the previous considerations, we recognize that anisotropic diffusion does really simplify the
original image in a steady way.
5 Numerical aspects
To approximate CED numerically, we replace the derivatives by finite differences. Since continuous
CED has the structure
5 Without loss of generality we omit the normalization.
its simplest discretization is given by the finite difference scheme
\Deltat
In this notation, U describes a vector containing the values at each pixel. The upper index denotes
the time level and L ij is a central difference approximation to the operator @ x
Since (27) can be rewritten as
I \Deltat
we observe that U k+1 can be calculated explicitly from U k without any matrix inversions. For this
reason it is called an explicit scheme.
Unfortunately, such explicit schemes require very small time steps \Deltat in order to be stable. Therefore,
it is desirable to replace (28) by an implicit scheme which has the same first-order Taylor expansion
in \Deltat, but better stability properties. One possibility is the AOS-stabilized scheme
ll
I \Deltat
This method achieves a stabilization through the nonnegative matrices
a semi-implicit discretization of the diffusion caused by the l-th diagonal entry of the diffusion tensor.
The typically used step size \Deltat = 2 is about one order of magnitude larger than the ones for which
an explicit scheme of type (27) is stable in 2-D or 3-D. For the standard approximations with central
derivatives within a (3 \Theta 3)-stencil, the matrix inversions in (29) come down to solving diagonally
dominant tridiagonal systems of linear equations. This can be performed in linear complexity with
a modified Gaussian algorithm (Thomas algorithm). It should be noted that (29) uses an additive
operator splitting (AOS) instead of the usual multiplicative operator splitting from the mathematical
literature such as [34]. This guarantees that all axes are treated in the same way. AOS schemes have
been introduced in [57] as an efficient and reliable method for isotropic nonlinear diffusion filtering,
which perform well on parallel computer architectures [58]. Equation (29) shows one way how to
extend them to anisotropic processes with a diffusion tensor.
The final CED algorithm is as follows:
1. Calculation of the structure tensor in each pixel. This requires one convolution with a Gaussian
K oe , and
convolutions with K ae . The convolutions were implemented in
the spatial domain exploiting the separability and symmetry of the Gaussian.
2. Principal axis transformation of the structure tensor in each pixel. In our case this was done by
a cyclic Jacobi algorithm [47].
3. Calculation of the diffusion tensor in each pixel. Here a significant speed-up is possible by creating
a look-up table for the function (12) at the beginning of the programme run.
4. Calculation of
I \Deltat
5. Calculation of
l :=
ll
by means of the Thomas algorithm [57].
6. Calculation of
U
This algorithm is computationally less complex than it looks at first glance: typical execution times for
one iteration on a HP 9000/889 are 0.3 seconds for a 256 \Theta 256 image and 9 seconds for a 64 \Theta 64 \Theta 64
image.

Figure

2: Anisotropic equations applied to the fingerprint image. (a) Top Left:
Original image. (b) Top Right: Mean curvature motion,
Left: Anisotropic diffusion similar to the Cottet-Germain model,
Coherence-enhancing anisotropic diffusion,
6 Examples

Figure

2 illustrates the importance of choosing the appropriate smoothing direction in anisotropic
processes. In Figure 2(b) mean curvature motion is applied to the fingerprint image. We observe
that - although mean curvature motion is a completely anisotropic technique with many merits in
computer vision - it cannot be used for closing interrupted line-like structures: topologically connected
structures remain connected, but gaps become larger, since both ends move in opposite directions.
MCM does not exploit any semilocal information.

Figure

2(c) shows the effect of coherence-enhancing anisotropic diffusion where ff and the integration
scale ae are set to 0. This model is close in spirit of the Cottet-Germain filter [15]. For such a model
it is difficult to find a useful noise scale oe for smoothing along flow-like structures: for oe ! 0, ru
becomes eigenvector of D and the diffusion is halted, since the corresponding eigenvalue - 1 is zero. If
oe becomes larger, the effect illustrated in Figure 1(c) dominates: opposite gradients cancel one another
leading to more or less random smoothing directions.
Coherence-enhancing diffusion filtering with a nonvanishing integration scale is depicted in Figure 2(d).
We observe that diffusion along the coherence orientation w 2 is well-suited for closing interrupted
lines in flow-like textures. Due to its reduced diffusivity at noncoherent structures, the location of
the semantically important singularity in the fingerprint remains the same. This is an important
prerequisite that any image processing method has to satisfy if it is to be applied to fingerprint

Figure

3: Scale-space behaviour of coherence-enhancing diffusion
(a) Top Left: Original fabric
analysis. In this and all subsequent CED images the parameters
used.

Figure

3 depicts the scale-space behaviour of coherence-enhancing anisotropic diffusion applied to a
fabric image arising in computer aided quality control. The temporal behaviour of this diffusion filter
seems to be appropriate for visualizing coherent fibre agglomerations (stripes) at different scales, a
difficult problem for the automatic grading of nonwovens.
Let us now investigate the impact of coherence-enhancing diffusion on images, which are not typical
texture images, but still reveal a flow-like character. To this end, we shall process expressionistic
paintings by Vincent van Gogh.
Fig. 4 shows the restoration properties of coherence-enhancing anisotropic diffusion when being applied
to a selfportrait of the artist [19]. We observe that the diffusion filter can close interrupted lines and
enhance the flow-like character which is typical for van Gogh paintings.
The next painting we are concerned with is called "Road with Cypress and Star" [20]. It is depicted
in Fig. 5. In order to demonstrate the influence of the integration scale ae, all filter parameters are
fixed except for ae. In Fig. 5(b) we observe that a value for ae which is too small does not lead to the
visually dominant coherence orientation and, thus, the filtered structures reveal a lot of undesired
fluctuations. Increasing the value for ae improves the image significantly (Fig. 5(c)). Interestingly, a
further increasing of ae does hardly alter this result (Fig. 5(d)), which indicates that this van Gogh
painting possesses a uniform "texture scale" reflecting the characteristic painting style of the artist.

Figure

6 illustrates the potential of CED for medical applications. It depicts a human bone. Its internal

Figure

4: Image restoration using coherence-enhancing anisotropic diffusion. (a)
Left: "Selfportrait" by van Gogh (Saint-R'emy, 1889; Paris, Muse'e d'Orsay),
structure has a distinctive texture through the presence of tiny elongated bony structural elements,
the trabeculae. There is evidence that the trabecular formation is for a great deal determined by the
external load [60, 35, 6]. For this reason the trabecular structure constitute an important clinical
parameter in orthopedics. Examples are the control of recovery after surgical procedures, such as
the placement or removal of metal implants, quantifying the rate of progression of rheumatism and
osteoporosis, the determination of left-right deviations of symmetry in the load or establishing optimal
load corrections for physiotherapy. The high resolution of contemporary CT using slipring technology
and MR scanners now enables in vivo analysis of detailed trabecular structure.
While in vivo analysis of the trabecular bone by means of a structure tensor is investigated in [24, 36],
the present paper points out ways to enhance the coherence information in medical images, in order
to ease such an analysis.
From

Figure

6(b),(c) we observe that CED is indeed capable of closing interrupted lines. Parallel flow-
like structures are enhanced and a subsequent coherence analysis becomes much easier. Increasing
the time t gives a coarser representation of the coherence. The entire evolution creates a task-driven
scale-space which is dedicated to flow-like patterns.
In medical applications physicians often insist in only very small amounts of smoothing, in order
to reduce some noise without affecting clinically relevant structures too much. Simplifications as in

Figure

6(c) would be considered as misleading. Figure 7 shows the results for two iterations of three-dimensional
CED filtering. The data set reveals a size of 256 \Theta 256 \Theta 128, and the voxel dimensions are
0:25 \Theta 0:25 \Theta 0:5 mm 3 . It depicts a CT scan of a foot area by means of a Philips SR8000 slipring CT
at 140 kV. Even for 3-D data sets, such a slight preprocessing using only a few iterations can easily
be performed within the acquisition time. The price of a suitable hardware would be a small fraction
of the price for a CT scanner.
It should be noted that in some 3-D applications it might also be desirable to smooth along two
directions instead of one as in CED. Diffusion methods of this type have been studied by Rambaux
and Gar-con [43]. They are 3-D versions of the edge-enhancing anisotropic diffusion filter from [51].
7 Summary and conclusions
In the present paper we have treated the problem of enhancing flow-like patterns. For such tasks a
reliable measurement of local orientation is needed. Our experiments demonstrate that the structure

Figure

5: Impact of the integration scale on coherence-enhancing anisotropic
diffusion 8). (a) Top Left: "Road with Cypress and Star"
by van Gogh (Auvers-sur-Oise, 1890; Otterlo, Rijksmuseum Kr-oller-M-uller),
Filtered with ae = 1. (c) Bottom Left:

Figure

(a) Top: High resolution slipring CT scan of a femural bone, showing the trabecular
formation. Slice thickness 1 mm, field of view 6
Left: Filtered by coherence-enhancing anisotropic diffusion, oe =0:5, ae =6, t =16. (c) Bottom
Right: Dito with t =128.
tensor satisfies this requirement. Unlike many other applications, we do not restrict its application to
pure image analysis, we use it as a tool for steering a scale-space evolution.
To this end, we evolve the original image by means of a nonlinear anisotropic diffusion equation. Its
diffusion tensor reflects the local image structure by using the same set of eigenvectors as the structure
tensor. The eigenvalues are chosen in such a way that diffusion acts mainly along the direction with
the highest coherence, and becomes stronger when the coherence increases. The resulting coherence-enhancing
diffusion (CED) process gives a theoretically well-founded scale-space representation: proofs
have been presented which show that its unique solution is stable under perturbations of the initial
image, satisfies Koenderink's causality requirement, and creates a large family of Lyapunov functionals
which ensure that the process is image simplifying and converges to a constant steady-state.
This theoretical framework is valid in any dimension. We have also presented an m-dimensional numerical
scheme where a stabilization based on additive operator splitting (AOS) allows time steps
which are one order of magnitude larger than in the explicit case.
Finally, the use of CED was illustrated by applying it to 2-D and 3-D data sets. The results indicate
a variety of possible application areas ranging from computer aided quality control over fingerprint
enhancement to medical imaging. This, however, gives only one aspect of the filter concept. Coherence-enhancing
diffusion allows also generalizations to vector-valued images. Details in this direction can
be found in [52].

Acknowledgments

. The author thanks Robert Maas and Stiliyan Kalitzin for useful comments on
a draft version of this paper, and Wiro Niessen for providing the 3-D trabecular bone data set.

Figure

7: Three-dimensional CED applied to a CT data set of size 256 \Theta 256 \Theta 128.
sections depicting slice no. 180 in x, no. 128 in y, and no.
100 in z, respectively. Right Column: Filtered, oe =0:5, ae =6, t=4.



--R

Enhancement of fingerprint images using shape-adapted scale-space operators

Axioms and fundamental equations in image processing
Image selective smoothing and edge detection by nonlinear diffusion.
Uniqueness of the Gaussian kernel for scale space filtering
The orientation of apatite crystals in bone
Optimal orientation detection of linear symmetry

Op'erateurs maximaux monotones et semi-groupes de contractions dans les espaces de Hilbert

Maximum entropy in action
Adaptive smoothing respecting feature directions
Image selective smoothing and edge detection by nonlinear diffusion
Nonlinear PDE operators with memory terms for image processing
Image processing through reaction combined with nonlinear diffusion

A fast operator for detection and precise location of distinct points
Nonlinear anisotropic filtering of MRI data

Road with cypress and star
Signal processing for computer vision



Representations based on zero-crossings in scale space
Basic theory on normalization of pattern (in case of typical one-dimensional pattern)
Lecture Notes in Comp.
Computer Vision
Maximum and comparison principles for anisotropic diffusion
Toward a computational theory of shape: An overview
The structure of images
Directional anisotropic diffusion applied to segmentation of vessels in 3D images

Splitting and alternating direction methods
Die Architektur der Spongiosa

Nonlinear multiscale representations for image segmentation

Nonlinear image filtering with edge and corner enhancement
Classification and uniqueness of invariant geometric flows
An adaptive and constrained model for 3D X-ray vascular reconstruction
Scale space and edge detection using anisotropic diffusion
Nonlinear anisotropic diffusion filtering of 3D images
Computing oriented texture fields
Generic evolution of edges on families of diffused greyvalue surfaces


ter Haar Romeny
Probabilistic multiscale image segmentation

Theoretical foundations of anisotropic diffusion in image processing

A review of nonlinear diffusion filtering
Anisotropic diffusion in image processing
On the history of Gaussian scale-space axiomatics





Structure adaptive anisotropic filtering
--TR

--CTR
S. Tabik , E. M. Garzn , I. Garca , J. J. Fernndez, High performance noise reduction for biomedical multidimensional data, Digital Signal Processing, v.17 n.4, p.724-736, July, 2007
Jose-Jesus Fernandez , Sam Li, Anisotropic Nonlinear Filtering of Cellular Structures in Cryoelectron Tomography, Computing in Science and Engineering, v.7 n.5, p.54-61, September 2005
J.-O. Lachaud , B. Taton, Deformable model with a complexity independent from image resolution, Computer Vision and Image Understanding, v.99 n.3, p.453-475, September 2005
J. Martn-Herrero, Hybrid object labelling in digital images, Machine Vision and Applications, v.18 n.1, p.1-15, January 2007
Guy Gilboa , Nir Sochen , Yehoshua Y. Zeevi, Image Sharpening by Flows Based on Triple Well Potentials, Journal of Mathematical Imaging and Vision, v.20 n.1-2, p.121-131, January-March 2004
Leah Bar , Nahum Kiryati , Nir Sochen, Image Deblurring in the Presence of Impulsive Noise, International Journal of Computer Vision, v.70 n.3, p.279-298, December  2006
scale-space;nonlinear diffusion;image enhancement;texture