# Structurally adaptive mathematical morphology based on nonlinear scale-space decompositions.

INTRODUCTION

Mathematical morphology is a well-known nonlinear image processing methodology based on the application of lattice theory to spatial structures (Serra. 1988; Heijmans and Ronse, 1990). Let E be a subset of the Euclidean Rd or the discrete space [Z.sup.d], considered as the support space of the image, and let T be a set of grey-levels, corresponding to the space of values of the image. It is assumed that T = [bar.R] = R[union]{-[infinity], +[infinity]} a grey-level image is represented by a function.

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

i.e., f [member of] F(E, T) maps each pixel x [member of] E into a grey-level value t [member of] ST: t = f(x). The two basic morphological mappings F(E, T) [right arrow] F(E, T) are the grey-level dilation and the grey-level erosion given respectively by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (2)

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (3)

where f [member of] F(E, T) is the original grey-level image and b [member of] F(E, T), is the fixed structuring function. The further convention for ambiguous expression is considered: f{x - h) + b(h) = -[infinity] when f(x - h) = -[infinity] or b(h) = -[infinity], and that f{x + h) - b(h) = +[infinity] when f(x + h) = +[infinity] or b(h) = The other morphological operators, such as the opening and the closing, are obtained as products of dilation/erosion. For instance, starting from the adjunction pair {[[delta].sub.b], [[epsilon].sub.b]}, the opening and closing of a grey-level image f according to the structuring function b are the mappings F(E, T) [right arrow] F(E, T) given respectively by

[[gamma].sub.b](f)(x) = [[[delta].sub.b]{[[epsilon].sub.b]{f))]{x),

[[psi].sub.b](f)(x) = [[[delta].sub.b]{[[epsilon].sub.b]{f))]{x) (4)

Particularly interesting in theory and in practical applications (Soille, 1999), the flat grey-level dilation and erosion is obtained when the structuring function is flat and becomes a structuring element. More precisely, a flat structuring function of support subspace B is defined as

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII],

where B is a Boolean set, i.e., B [subset or equal to] E or B [member of] P(E), which defines the "shape" of the structuring element. We notice that [B.sup.c] denotes the complement set of B (i. e., B [intersection] [B.sup.c] = [empty set] and B [union] [B.sup.c] = E). The structuring element is defined at the origin o [member of] of] E, then to each point p of E corresponds the translation mapping o to p, and this translation maps B onto [B.sub.p], i.e., [B.sub.p] = {b + p : b [member of] of] B}. Therefore, the flat dilation and erosion of a grey-level image f{x) with respect to the structuring element B are respectively

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (5)

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (6)

where [??] is the reflection of B with respect to the origin. i.e.,[??] = {-b|b [member of] B}.

Standard formulation of morphological operators is translation invariant in the space ("horizontal" direction invariance) and in the grey-level intensity ("vertical" direction invariance), see properties and proofs in (Heijmans and Ronse, 1990), i.e., for a given translated image g(z) = f(x + y) + r, it is obtained that [[delta].sub.b](g)(z) = [delta]b(f)(x + y) + r. But of course, this result is true if and only if the same structuring function b(x) is considered for each point x of the image. A current challenging topic in mathematical morphology is the construction of adaptive operators; or in other words, operators whose structuring functions become dependent on position or on the input image itself. In previous works, the proposed adaptive operators have been based on two main approaches.

--On the one hand, a variability on the support space E: spatially variable shape of structuring functions according to i) the geometric position in the image: perspective-adapted Morphology by Beucher et al. (1987) and by Cuisenaire (2006), ii) the local regularity of image values: morphological amoebas by Lerallut et al. (2005), intrinsic structuring elements by Debayle and Pinoli (2005) morphological bilateral filtering by Angulo (2011), hi) the orientation: curvilinear morphohessian filter by Tankyevych et al. (2008). locally-variant anisotropic morphological filters by Verdu-Monedero et al. (2011).

--On the other hand, a variability on the value space T: variable size of structuring functions according to the local intensity or contrast: viscous watershed by Vachier and Meyer (2005) and viscous openings/closings by Vachier and Meyer (2007), together with a PDE formulation by Maragos and Vachier (2008).

For a recent overview on the state-the-art on adaptive morphology, the interested reader is invited to papers by Maragos and Vachier (2009) and by Verdu-Monedero et al. (2011). Another recent study by Roerdink (2009) is very interesting for understanding the theoretical limitations of input-adaptive morphological operators.

Motivation of the new adaptive framework

The intensity-adaptive operators, or level-varying operators, studied in Vachier and Meyer (2007) and Maragos and Vachier (2008) are typically implemented using a level-set decomposition, followed by the processing of each level set by a flat operator of a size which depends on the index of the particular level set and finally the image reconstitution from the processed level sets. Hence, the adaptavility for each point depends on the absolute intensity value. Two points of the image associated to structures of similar scale but of different intensity are processed differently. This work introduces a new framework of nonlinear adaptive operators: the structurally adaptive mathematical morphology. More precisely, the rationale behind the present approach is to work on a nonlinear multi-scale image decomposition based on a morphological family of operators called levellings (Meyer, 1998; Meyer and Maragos, 2000). Then, the size of the structuring function b{x) at point x is intrinsically adapted to the local scale of the structures. In other words, two regions of the image appearing in a similar scale will undergo a nonlinear processing of similar "size" independently of the absolute intensities. The value of the processed image will be obtained by the reconstitution of the processed scales. We notice also that our approach is not limited to flat morphological operators. The properties of the derived operators are investigated and their practical performances are compared with respect to standard morphological operators using natural image examples. The paper is an extended and improved version of the conference paper (Angulo and Velasco-Forero, 2010).

PAPER ORGANISATION

The rest of the paper is organised as follows. In the first part of Methods Section is presented the structural image model, using a pyramid of levellings, which is the basic ingredient for the present adaptability framework. Then, the second part of Methods Section introduces the structurally adaptive pseudo-dilation and pseudo-erosion and other associated operators, including also the study of the algebraic properties of these new operators. Results Section discusses various examples which illustrate the behavior of the structurally adaptive operators for typical applications of mathematical morphology. Finally, Section on conclusions and perspectives closes the paper.

METHODS

STRUCTURAL IMAGE MODEL USING A PYRAMID OF LEVELLINGS

The classical tool from mathematical morphology for multi-scale image decomposition is based on the granulometric analysis (Serra, 1982; Maragos, 1989 A granulometry is the study of the size distribution of the objects of an image (Serra, 1982). Formally, forthe discrete case, a granulometry is a family of openings T = [{[[gamma].sub.n]}.sub.n [greater than or equal to o]] that depends on a positive parameter n (which expresses a size factor) such as: i) [[gamma].sub.0] (f)(x) = f (x); ii) f(x)[less than or equal to]g(x) [[gamma].sub.n] (f)(x)[less than or equal to][[gamma].sub.n] (g)(x),[for all] [greater than or equal to] 0, [for all],g [[member of] of] F(E, T); iii) [[gamma].sub.n] (f)(x) [less than or equal to] f(x), [[for all].sub.n] [greater than or equal to] 0, [for all] [[member of] of] F(E, and T)and iv) [[gamma].sub.n] verifies the granulometric semi-group absorption law; i.e.,

[[gamma].sub.n] [[gamma].sub.m] = [[gamma].sub.m] [[gamma].sub.n] = [[gamma].sub.max] (n,m), [for all], m [greater than or equal to] 0.

Moreover, a granulometry by closings (or anti-granulometry) can be defined as a family of increasing closings [THETA] = [{([[sigma].sub.n]}.sub.n [greater than or equal to]0]. In practice, the most useful granulometry and anti-granulometry are those associated to morphological flat openings and closings: [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] respectively, where [B.sub.1] = B is a isotropic structuring element of unit size (unitary discrete disk), and [B.sub.n] is the homothetic structuring element of size n, with n = 1,2, .... The granulometric analysis of an image / with respect to a finite T consists in computing the pyramid of opened images and then in calculating the image residue between two successive openings:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

The pyramid of residues of openings [{[[rho].sup.+.sub.i] (f)(x)}.sup.N.sub.i=1] leads to a decomposition of bright structures into N scales. We recall that the subgraph of the opened image by [B.sub.i] is equivalent to the union of the translations of the structuring when it fits the subgraph of the original image (Serra, 1982), or in other terms, image structures in the foreground (i.e., bright ones) that cannot contain [B.sub.i] are removed by the opening. Therefore, the residue [[rho].sup.+.sub.i] extracts the structures removed between [B.sub.i-1] and [B.sub.i]. By duality, the image decomposition into N scales of dark structures is obtained from the pyramid of residues of closings [{[[rho].sup.-.sub.i] (f)(x).sup.N.sub.i=1] where

[[rho].sup.-.sub.i](f)(x) = [gamma][B.sub.i](f)(x) - [gamma][B.sub.i-1](f)(x).

In summary, using morphological openings and closings, it is obtained a pair of nonlinear scale-space image decomposition, i.e.,

f(x) [??] [{{[[rho].sup.+.sub.i](f)(x)}.sup.N.sub.i=1]; {[[rho].sup.-.sub.i](f)(x)}.sup.N.sub.i=1]}

Furthermore, morphological openings/closings deal selectively with the grey-level structures according to their support size but they also "distort" the contours (level sets) of the image objects. In mathematical morphology, the powerful geodesic operators (Serra. 1988; Soille, 1999) (or connected filters) are perfectly able to deal with the objects, without distorting their contours.

For all these reasons, we have considered for the present structurally adaptive operators a multiscale framework which, on the one hand, decomposes simultaneously bright and dark structures (by means of alternate sequential filters); and on the other hand, uses geodesic operators (in particular, the levelling operator) to extract precisely the grey-level objects. Let us consider in detail this methodology.

Nonlinear scale-space based on viscous levellings

The levelling (Meyer, 1998) of the reference image f(x) according to the marker image m(x) (this latter is generally a rough simplification of the reference image), denoted [lambda](f,m)(x), is a morphological geodesic filter which simplifies textures and eliminates small details of the reference image according to the marked structures, but preserving the contours of remaining objects. In practice, the levelling is obtained by iteration of geodesic dilations and erosions until idempotence. More precisely, a straightforward algorithm is as follows (Meyer, 1998):

[lambda] (f, m) (x) = [[LAMBDA].sup.i] (f (x), m (x)) (7)

such that the idempotence has been reached, i.e.,

[[LAMBDA].sup.i] (f,m) (x) = [[LAMBDA].sup.i + 1] (f, m)

where the i-iteration is

[[LAMBDA].sup.i] (f, m) = (f [conjunction] [alpha][[LAMBDA].sup.i-1](f, m))) [disjunction] [beta] [[LAMBDA].sup.i-1] (f, m)), (8)

with the initial step

[[LAMBDA].sup.1] (f, m) = (f [cojunction] [alpha](m)) [disjunction] [beta](m). (9)

Levellings of different nature are obtained according to the pair of transformations ([alpha],[beta]) used in Eq. 8 and 9. More precisely, two of the most interesting ones are (Meyer, 1998)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

Fig. 1 shows a comparison of standard levelling vs. viscous levelling using three different marker images. We observe that even the structures which have been removed from the marker can be reconstructed (by strict geodesic propagation) by standard levelling. On the contrary, viscous levelling produces a reconstruction more faithful to the simplification degree driven by the marker. Consequently, by its skilful properties on structure regularisation, we consider here the viscous levelling as the fundamental operator for our structural decomposition.

[FIGURE 1 OMITTED]

A pyramid of levellings is then obtained by associating to the image a pyramid of markers, i.e.,

f(x) [??] [{[m.sub.i](x)}.sup.N.sub.i=1] [??] [{[lambda](f, [m.sub.i])(x)}.sup.N.sub.i=1] (10)

which leads to a nonlinear scale-space (Meyer and Maragos, 2000). The marker at level i determines the degree of simplification of the corresponding levelled image; the simplification effect increases with i.

The levelling can be considered as the activity supremum between an opening by reconstruction and a closing by reconstruction (Meyer, 1998), leading to a brigh/dark symmetric operator. Let Cf stand for the image complement (or negative) of f [[member of] of] F(E, T), i.e., Cf(x) = -f(x) if T = [bar.R] and Cf (x) = [t.sub.max] - f(x) if T = [0,[t.sub.max]]. We have for the levelling the following property: C[lambda](Cf,Cm) = [lambda](f, m), which means that (f, g) [??] [lambda](f,g) is always a self-dual mapping. In addition, if m derives from f by a self-dual operation, i.e., m = m(f) with Cf(Cf) = m(f) then levelling [??] [lambda](f, m(f)) is self-dual. We notice that this important relation is distinct from that of invariance under complement m(Cf) = m(f) which does not imply self-duality for [??] A (f, m(f)). Hence, to obtain a self-dual and symmetric processing of brigh/dark structures, a self-dual family of marker images is required. Classical examples of self-dual filters are the convolution with Gaussian kernels or the (weighted) median filters. In particular, the nonlinear scale spaces associated to levellings where the markers are Gaussian filters have been studied by Meyer and Maragos (2000). In addition, this kind of marker was also used in (Sofou et al, 2005) to extract the texture layer of an image as the residue of the corresponding levelling.

Nevertheless, we propose in this study to use as multi-scale markers the averaged alternate sequential filters (AASF):

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (11)

where the pair of alternate sequential filters (ASF) (Serra, 1988) are given by the following products of openings and closings

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (12)

The properties of ASF are well known in mathematical morphology (Serra, 1988). In particular, they are idempotent and they follow the semi-group absorption law, i.e., [for all] n [greater than or equal to] m [greater than or equal to] 0,we have

[ASF.sup.I,II.sub.n], [ASF.sup.I,II.sub.m] = [ASF.sup.I,II.sub.n], [ASF.sup.I,II.sub.n] = [ASF.sup.I,II.sub.n = max(n, m)]

and consequently this property is inherited by [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. ASF are very useful for noise cancellation and smoothing since they produces progressive simplification of image structures. This behaviour is related with their property of compatibility with scale modifications (i.e., intuitively that means that the filter of order k works at a scale k as a filter of order 1 at a scale 1). Theoretically, ASF are not self-dual, but our symmetric averaged version [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] leads to an quasi self-dual operator: we have observed that in practice the difference between [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (Cf)(x) and [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (f)(x) is always negligible. The main drawback of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is the computation time, which increases with the size n of the scale. However, using well-known optimized implementations of discrete dilation/erosion for the openings/closings, the final computation time is not nowadays a limiting point. The rationale behind this choice of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] , instead of a Gaussian filter of equivalent scale n, is because we have observed empirically that the viscous levellings associated with alternate sequential filter are more effective for extracting the bright/dark structures of corresponding scale n.

[FIGURE 2 OMITTED]

Structural image model

We have now all the ingredients to introduce the following structural model image:

f(x) [equivalent to] c + s(x) + [N.summation over (i=1)] [t.sub.i](x), (13)

where c is the continuous component, s(x) is the structure component, and [{[t.sub.i](x)}.sup.N.sub.i = l] are the texture components at N scales.

More precisely, given the levelling scale-space [{[lambda] (f, [m.sub.i])(x)}.sup.N.sub.i = l], the components of the present image model can be calculated as follows.

--The continuous component is obtained as the mean value of the last levelled image, i.e.,

c = mean([lambda](f,[m.sub.N])(x)). (14)

We can consider that c varies with the DC level of the image acquisition device. In the operators presented in this paper c does not give any useful information but is still necessary for ongoing research on adaptive self-dual operators (see perspectives section);

--The structure component corresponds to the mean centred last levelled image, i.e.,

s(x) = [lambda]{f,[m.sub.N])(x) - c. (15)

Therefore, the image s(x) represents the most simplified image according to the choice of scale N. In other words, the level-sets of the significant objects having a "support size" bigger than the structuring element [B.sub.N].

--The texture components are obtained as the derivative of the levelling pyramid, i.e.,

[t.sub.1](x) = [lambda](f,[m.sub.i=1])(x) - [lambda](f,[m.sub.i])(x), (16)

where [m.sub.o] (x) = f(x). Each residue image [t.sub.i], represents the brigh/dark objects which were in scale image i - 1 but have been removed in scale image i; that is, level-sets of objects whose "support size" is comprised between both successive scales. The intensity at each point of [t.sub.i] is the relative intensity of the object (which is invariant to monotone increasing intensity transformation on original image). In fact, the residue image between two successive levellings can be considered as a "morphological laplacian pyramid". In natural images, we can suppose that the texture components are sparse images (i.e., only a limited number of pixels are different of 0) and that the complexity of the original image is split into a simpler image (structure component) and a scale-complexity bounded series of images (texture components).

Our additive image decomposition model can be considered as generalisation of the Cartoon + Texture decomposition by Meyer (2002): f(x) = u(x) + v(x) where u(x) is the cartoon component (homogeneous zones of the objects) and v(x) is the texture oscillation, which is usually solved using variational algorithms (Aujol et al, 2006). The terms can be identified as u(x) = c + s(x) and v(x) = [[summation].sup.N.sub.i = 1] (x). The value N of the larger structuring element (scale of structure equivalent) in our model is equivalent to the "scale size" in other Cartoon + Texture decomposition algorithms (Buades et al, 2010).

[FIGURE 3 OMITTED]

Fig. 2 gives an example of three levels of a viscous levelling pyramid and the corresponding residues images, used to define the texture components. We notice that the texture components, with respect to the structure component, are not defined in terms of periodicity or regularity as in the nowadays extended Cartoon+Texture decomposition. That can be observed in the comparative example of our structural decomposition (size of largest structure scale N = 8) and the Cartoon+Texture decomposition recently proposed by Buades et al. (2010) (scale size [sigma] = 3) given in Fig. 3. The Cartoon + Texture decomposition has been obtained using the algorithm available online in (Buades et al, 2009). Therefore, the notion of texture matches here exclusively with the notion of scale. We remark also that, in comparison with a laplacian pyramid associated to a Gaussian scalespace, the main advantages of the derivative levelling pyramid are, on the one hand, the preservation of the contours of extracted structures and on the other hand, at each scale, the size of the structures is bounded by the size of the AASF filter used as marker.

The structural image model introduced above is the key element for the structurally adaptive operators studied in this section. As we will show below, the two basic adaptive operators present a lack of some algebraic properties; and consequently they are not dilations and erosions in stricto sensu. Therefore, we prefer to them structurally adaptive pseudo-dilation and pseudo-erosion.

[FIGURE 4 OMITTED]

Given the image f [[member of] of] F (E, T), structurally decomposed into the scales {m, m + k, m + 2k, * * * ,M} (denoted compactly by {m : k : M}) according to the model in Eq. 13, we define the corresponding structural adaptive dilation as

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (17)

and the structural adaptive erosion as

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (18)

where [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] are the standard dilation and erosion of image g(x) according to the fixed spatially invariant structuring function [b.sub.n](x).

We have considered in this study examples using the two most useful families of isotropic multi-scale structuring functions:

--parabolic shape function of width n, i.e., [b.sub.n](x) = -[[parallel]x[parallel].sup.2]/2n;

--flat disk of radius n, denoted by [B.sub.n].

For the examples given in the paper, the family of disks [B.sub.n] used in the flat operators are implemented digitally by hexagons of radius n pixels. This choice is based on the properties of isotropy of hexagonal grid with respect to the square one and the fact that the hexagons are an approximation to the disk which can be implemented using optimized algorithms. In any case, if the computation time is not critical, the flat disks can be obtained by thresholding the paraboloidal structuring functions. Fig. 4 depicts a comparison of three flat structural adaptive pseudo-dilations for M = 4, 8 and 12 (with m = 2 and k = 2), versus the corresponding standard flat dilations. We notice that the "size" of the operator is determined by the size of the biggest scale M; or in other words, the size of the scale of the structure component. The two other parameters from the structural decomposition: i) the size of the smallest scale m, and ii) the sampling parameter k (or step size between two successive scales) have a lower impact on the adaptive operator, see for instance the example given in Fig. 5. This effect allows to introduce a sampling parameter k > 1 in order to reduce the number of processed discrete scales (reduction of computation time).

[FIGURE 5 OMITTED]

Let us consider the algebraic properties of these two morphological operators. The proofs are not included but are easily obtained from the standard properties of erosion and dilation.

--The structural adaptive pseudo-dilation (pseudo-erosion) is extensive (anti-extensive), i.e., [for all]f [[member of] of] F(E, T),

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

--The structural adaptive pseudo-dilation and pseudo-erosion are dual operators with respect to the grey-level inversion, i.e., [for all]f [member of] F{E, T)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

- However, the fundamental law which links the pair of dilation/erosion, the adjunction property, fails for the introduced structurally adaptive operators. [for all]f, g [member of] F (E, T)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

As recently pointed out by Roerdink (2009), this drawback is common to any adaptive morphology dilation/erosion where the processing at each point depends on the input value.

--Moreover, as a consequence of the lack of adjunction property, other properties are also lost. Typically, the structural adaptive pseudo-dilation (pseudo-erosion) does not commute with the supremum (infimum), i.e., [for all] f, g [member of] F (E, T)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

--Another fundamental property of standard morphology which fails is the increasiness. The structural adaptive pseudo-dilation/pseudo-erosion are generally not increasing operators, i.e., f [less than or equal to] g does not always involve that [[??].sub.m:k-M](f) [less than or equal to] [[??].sub.m:k-M] (g) ? In fact, the property is verified if

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

or in other terms, the ordering must be preserved in the derivative pyramid, which is not always the case.

In mathematical morphology, a pair of dual operators which are extensive and anti-extensive are respectively a thickening and thinning (Serra. 1982; Soille, 1999); therefore, the structurally adaptive operators are just algebraic grey-level thickenings/thinnings. In any case, despite the lack of some properties, the operators [[??].sub.m:k-M] (f) and [[??].sub.m:k-M] (F) can be used to construct structurally adaptive gradient, laplacian, toggle mapping, etc.

Due to the fact that the introduced structurally adaptive operators are not formally a pair of adjunct dilation/erosion, their product [[??].sub.m:k-M] (f) is not an algebraic opening. Hence, we prefer to formulate the structurally adaptive pseudo-opening as

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

and dually, the structurally adaptive pseudo-closing is defined by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

The structural adaptive opening in Eq. 19 (closing in Eq. 20) is anti-extensive (extensive) but it is not increasing. The idempotency, [gamma][gamma](f) = f(f), is the property which guarantees the stability of standard opening/closing. Even if our adaptive operators are not strictly idempotent, their iteration is almost stable when AASF filters are used for the levelling decomposition.

[FIGURE 6 OMITTED]

As we can observe in the example given in Fig. 6. the adaptive opening [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] removes bright objects smaller than the size of the structuring function but. with respect to the standard opening [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], the notion of "smaller than" depends on the scale of decomposition of the object. The product of operators [[??].sub.m:k-M]{f) and [[??].sub.m:k-M](f) can be then used to define other more complex operators.

APPLICATIONS AND RESULTS

Let us illustrate the performance of the present operators by means of four typical applications of mathematical morphology.

[FIGURE 7 OMITTED]

The first case study, depicted in Fig. 7, corresponds to a coronary network acquired by contrast-enhanced radiography. Besides the vessels, the image presents a strong and irregular background (which includes other structures such as the ribs). In order to extract the vessels, the classical top-hat (Serra, 1982) is particularly useful:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (21)

i.e., residue between the original image f and the opening [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII](f). The size of the isotropic opening should be bigger than the diameter of the thickest vessels. For the particular example, we consider a parabolic structuring function of width n = 10. As we can observe, the standard opening produces a rough regularisation by removing all the peaks smaller than [b.sub.10] whereas the structurally adaptive opening removes also the peaks at different scales but preserves better the secondary structures; hence, the comparison of the standard [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (f) against the adaptive top-hat [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (f) shows that the second one presents a better extraction of thin vessels as well as less texture and irregularities.

The second application, given in Fig. 8, deals with an example of image denoising and detail simplification. The mosaic image is quite textured and noisy and the aim is to filter it out in order to obtain a more regular image where the main contours are well preserved. This goal can be achieved by means of the morphological centre operator (Serra, 1988).

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] defined in terms of two products of openings and closings of a chosen size:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (22)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (23)

We compare in particular the denoising effect using a flat structuring element of size N = 4 and a parabolic structuring function of similar width. As we can observe in the images, the adaptive flat or parabolic versions of the operator, [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (f) and [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (f) outperforms clearly the effects of the standard ones [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (f) and [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (f), yielding a more efficient denoising/simplification image, especially the parabolic case, with an excellent contour preservation of remaining structures.

In order to compare quantitatively the behavior of standard vs. adaptive morphological centre, we have evaluated in an experiment the quantitative assessment of noisy images. Fig. 9 summarizes the results of this experiment: starting from an image corrupted with four different levels of Gaussian noise and with four different levels of impulsive noise ("salt and pepper" noise), we have filtered out the eight images using standard parabolic centre [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] and structural adaptive parabolic centre [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], and then, we have computed the PSNR value with respect to the original image. For all the noise levels, the value of PSNR is better for the structurally adaptive case than for the standard one. We note also that the "size" of the filter 4 is not particularly appropriate for low Gaussian noise levels or for high impulsive noise levels, but even in such a case the trade-off between denoising and image preservation is better for the adaptive version of the morphological center.

[FIGURE 8 OMITTED]

The next problem concerns an aerial image of an airport, including objects of different size, but not well contrasted, see Fig. 10. The morphological laplacian (or nonlinear laplace filter as introduced by van Vliet et al. (1989)) is defined in terms of dilation and the erosion as

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (24)

A simple but efficient image contrast enhancement is obtained by subtracting the laplacian image from the original one:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (25)

[FIGURE 10 OMITTED]

We observe that the structurally adaptive laplacian [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] produces a more balanced contrast enhancement than the standard one [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. In particular, the small objects are enhanced in a small neighborhood and the large in a proportional neighborhood. In the standard case, all the structures are equally enhanced. By the way, the laplacian-based enhancement is an example of structurally adaptive operator which probably is better to be implemented using a full scale decomposition, with m = 1 and k = 1.

The last case study deals with an edge detector which can be implemented using only a pair of morphological erosion and dilation, and which is useful for complex images, like the one given in Fig. 11. In fact, the zero crossings of the laplacian [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] as m the classical Marr-Hildreth model. correspond to the edges of image f. In order to select only the most prominent edges, the zero crossing detector is "multiplied" by the binary image obtained with a threshold by hysteresis of the morphological gradient (Serra, 1982; Soille, 1999) given by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (26)

Besides the size parameter of the laplacian and the gradient, the single parameter of this edge detection model is the threshold value for the gradient. We observe again in this example that the structurally adaptive gradient [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], using a parabolic structuring function, is more appropriate than the standard counterpart. The comparison by simple visual assessment of edge detection is not conclusive, but it seems that in the adaptive case the edges of small objects are better detected.

CONCLUSIONS AND PERSPECTIVES

We have discussed an additive image model associated to a nonlinear multi-scale image decomposition using a family of viscous levellings. Then, working on this decomposition, we have introduced structurally adaptive morphological operators, where each component of the image is processed with a structuring function of size intrinsically adapted to the corresponding local scale.

We have shown that, in practical applications, the proposed operators perform better than standard ones for object extraction, image denoising, image enhancement, etc. This positive behaviour is justified by the fact that, using geodesic operators for image decomposition, the different objects are almost individually processed and the "interferences" between adjacent objects are notably reduced.

In addition, we have proved that another advantage of the adaptability is the fact that the choice of the "size" for the erosions/dilations and the openings/closings is much less critical.

However, we have also demonstrated that the underlying algebraic structure of the structurally adaptive operators is less rich than the standard ones.

[FIGURE 11 OMITTED]

As we have discussed, in the present image model, the structure and texture terms are positive/negative signals and consequently, in ongoing work, we will formulate adaptive self-dual erosion and opening. Finally, sparsity properties of our image model suggest us to consider it as a starting point to explore the notion of sparse mathematical morphology in future research.

doi: 10.5566/ias.v30.p111-122

REFERENCES

Angulo J, Velasco-Forero S (2010). Structurally Adaptive Mathematical Morphology on Nonlinear Scale-Space Representations. In: Proc IEEE Int Conf Image Proc (ICIP'10) 121-4.

Angulo J (2011). Morphological bilateral filtering and spatially-variant adaptive structuring functions. In: Soille P, Pesaresi M, Ouzounis G, eds. Proc Int Symp Math Morpho (ISMM'11). Lect Notes Comput Sci 6671:212-23.

Aujol JF, Gilboa G, Chan T, Osher S (2006). Structure-Texture Image Decomposition - Modeling, Algorithms, and Parameter Selection. Int J Comput Vision 67(1):111-36.

Buades A, Le T, Morel J.-M, Vese L.A (2009). Cartoon+Texture Image Decomposition. Image Processing On Line (IPOL) Journal, www.ipol.im.

Buades A, Le T, Morel J.-M, Vese L.A (2010). Fast cartoon + texture image filters. IEEE T Image Proc 19(18):1978-86.

Beucher S, Blosseville J.M, Lenoir F (1987). Traffic Spatial Measurements Using Video Image Processing. In: Casasent DP, Hall EL, eds. Proc SPIE 848:648-55.

Cuisenaire O (2006). Locally adaptable mathematical morphology using distance transformations. Pattern Recogn39(3):405-16.

Debayle J, Pinoli J (2005). Spatially adaptive morphological image filtering using intrinsic structuring elements. Image Anal Sterol 24:145-58.

Lerallut R, Decenciere E, Meyer F (2005). Image Filtering Using Morphological Amoebas. In: Proc Int Symp Math Morpho (ISMM'05) 13-22.

Maragos P (1989). Pattern Spectrum and Multiscale Shape Representation. IEEE T Pattern Anal 11:701-16.

Maragos P, Vachier C (2008). A PDE Formulation for Viscous Morphological Operators with Extensions to Intensity-Adaptive Operators. In: Proc IEEE Int Conf Image Proc (ICIP'08) 2200-3.

Maragos P, Vachier C (2009). Overview of adaptive morphology: trends and perspectives. In: Proc IEEE Int Conf Image Proc (ICIP'09) 2241-4.

Meyer F (1998). The levelings. In: Proc Int Symp Math Morpho (ISMM'98) 199-206.

Meyer F, Maragos P (2000). Nonlinear scale-space representation with morphological levelings. J Vis Commun Image R 11:245-65.

Meyer Y (2002). Oscillating Patterns in Image Processing and Nonlinear Evolution Equations. University Lecture Series Vol. 22, Providence: American Mathematical Society.

Heijmans H, Ronse C (1990). The Algebraic Basis of Mathematical Morphology. Part I: Dilations and Erosions. Comput Vision Graph 50(3):245-95.

Roerdink J (2009). Adaptivity and group invariance in mathematical morphology. In: Proc IEEE Int Conf Image Proc (ICIP'09) 2253-6.

Serra J (1982). Image Analysis and Mathematical Morphology, Vol I. London: Academic Press.

Serra J, ed. (1988). Image Analysis and Mathematical Morphology Vol.2: Theoretical Advances. London: Academic Press.

Sofou A, Evangelopoulos G, Maragos P (2005). Coupled geometric and texture PDE-based segmentation. In: Proc IEEE Int Conf Image Proc (ICIP'05) II:650-3.

Soille P (1999). Morphological Image Analysis. Berlin: Springer-Verlag.

Tankyevych O, Talbot H, Dokladal P (2008). Curvilinear morpho-hessian filter. In: Proc 5th IEEE Int Symp Biomed Imag: From Nano to Macro (ISBI'08) 1011-4.

Vachier C, Meyer F (2005). The Viscous Watershed Transform. J Math Imaging Vis 22:251-67.

Vachier C, Meyer F (2007). News fromviscousland. In: Proc Int Symp Math Morphol (ISMM'07) 189-200.

van Vliet L.J, Young I.T, Beckers G.L (1989). A Nonlinear Operator as Edge Detector in Noisy Images. Comput Vision Graph 45:167-95.

Verdu-Monedero R, Angulo J, Serra J (2011). Anisotropic Morphological Filters With Spatially-Variant Structuring Elements Based on Image-Dependent Gradient Fields. IEEE T Image Process 20(12):200-12.

JESUS ANGULO AND SANTIAGO VELASCO-FORERO

CMM-Centre de Morphologie Mathematique, Mathematiques et Systemes, MINES Paristech, France e-mail: jesus.angulo@mines-paristech.fr, santiago.velasco@mines-paristech.fr

(Accepted June 20, 2011)
```Figure 9. Evaluation of image denoising using
morphological center: first row and fourth row,
original images; second row and fifth row, parabolic
centre (paraboloid of width n=4); third and sixth row,
structural adaptive parabolic centre (M = 4, m = 1,
k = 1). Original image was corrupted with Gaussian
noise in first row: first column s = 10, second column
s = 20, third column s = 40, fourth column s = 80;
and with impulsive noise ("salt and pepper" noise
which consists in randomly replacing a percentage of
pixels with white pixels and the same percentage with
black pixels) in fourth row: first column Pr = 2.5%,
second column Pr = 5%, third column Pr = 7.5%,
fourth column Pr = 12.5% . Table (g) provides the
corresponding values of PSNR.

[sigma]   [sigma]  [sigma]  [sigma]
= 10      = 20     = 40     = 80

[[zeta]b.sub.4]           24.65      24.18    20.58    15.36
[[zeta]b.sub.1:1:4]       25.18      24.30    22.45    18.40
2.5%       5%       7.5%    12.5%
[[zeta]b.sub.4]           26.08      19.02    14.93    12.22
[[zeta]b.sub.1:1:4]       27.59      24.59    24.45    23.34
```
COPYRIGHT 2011 Slovenian Society for Stereology and Quantitative Image Analysis
No portion of this article can be reproduced without the express written permission from the copyright holder.