Ghulam-Sakhi Shokouh
email: ghulam-sakhi.shokouh@mines-ales.fr

Baptiste Magnier
email: baptiste.magnier@mines-ales.fr

Binbin Xu
email: binbin.xu@mines-ales.fr

Philippe Montesinos
email: philippe.montesinos@mines-ales.fr

Ridge Detection by Image Filtering Techniques: A Review and an Objective Analysis

Keywords: ridge detection, salient feature analysis, image filtering

are the useful geometric features due to their wide varieties of applications, mainly in image analysis problems such as object detection, image segmentation, scene understanding, etc. Many characterizations have contributed to formalize the ridge notion. The signification of each characterization rely however on its actual application. The objective analysis of ridge characterized as thin and complex image structure is thus essentially important, for choosing which parameter's values correspond to the suitable configuration to obtain accurate results and optimal performance. In this article an extensive analysis followed by a supervised and objective comparison of different filtering-based ridge detection techniques is led. Furthermore, the optimal parameter configuration of each filtering techniques aimed for image salient feature analysis tool have been objectively investigated, where each chosen filter's parameters corresponds to the width of the desired ridge or valley. At last, the comparative evaluations and analysis results are reported on both synthetic images, distorted with various types of noises and real images.

INTRODUCTION

Roughly, ridges (resp., valley) as one of the many images' salient features, are set of curves whose points are local maxima (resp., minima) in an image, as shown in Fig. 1. The precise detection, localization and extraction of this salient features along with their accurate characterization of its geometric structure are important image processing tasks, related to its wide range of application. Exhaustive researches have been accomplished on this significant image features. Ridges have shown to be the most eminent and useful structure amongst for image analysis and various related applications. Digital images comprise varying types of salient features, such as edges, blobs, corners, textures, whereas ridges (also called crest lines or roof edges) represent a special type of contours, as shown in Fig. 2. Classical edge detection techniques are optimized to extract step or ramp edges [START_REF] Canny | A computational approach to edge detection[END_REF]; nonetheless, they fail to detect ridges or crest lines in images. A step/ramp edge extraction functions return two edges at both sides of the crest line, because narrow ridges or valleys on the image surface are composed of two locally parallel step or ramp edges. Roof edges are defined as thin nets inside the image; describing roads or rivers in satellite images, blood vessels in medical images or plant roots. Therefore, finding these dense and thin structures is an important task in image processing.

Concurrently, there are many inescapable challenges in image processing task, such as noise and artifacts that necessitate to be researched concerning the specific application. Generally, the acceptable ridge/valley detection involves cumbersome and manual tuning in order to overcome the issues for specific application. In consequence, to contribute the research process in the domain of ridges detection and extraction techniques, an extensive evaluation of the different state-of-the-art filtering techniques and approaches in the scope of its most useful application, is crucially necessary. This article is an extension and/or improved version of [START_REF] Shokouh | An objective comparison of ridge/valley detectors by image filtering[END_REF] aimed for objective and extensive analysis of state-of-the-art filtering techniques for ridge (resp., valley) detection and extraction.

Regarding a curve in a gray-level image I, ridges correspond to directional maxima, valleys correspond to direction minima, resp. [START_REF] Haralick | Ridges and valleys on digital images[END_REF], as illustrated in Fig. with and being the partial derivatives of , respectively, along the and the axis. Ridges and valleys are given by the points where the values of are maxima (resp., minima) in the orthogonal direction of the curve at . Consequently, to detect and extract ridges and valleys in a signal, the Laplacian is used as a high pass filter to enhance its high frequency component. As a result, to extract peaks in a one-dimensional (1D) signal , the (continuous) Laplacian operator is simply the second derivative of : = . The discrete Laplacian is an approximation to the continuous Laplacian. The common approximation to the second derivative of a discrete signal is:

(

∂ = = ∂ ∂   =   ∂   , , 1 0 , , , 1 0 x
∂ = ∂ ∂   = =   ∂   , , , ( ) ( ) ( ) ( ) ( )
, x y ( ) S t Δ S ( ) ΔS t ( ) ∂ ∂ 2 2 S t t ( ) S t ( ) ( ) ( ) ( ) { } - -- + Δ = ∈ … 2 2
with 1,2,3 ;

S t S t h S t h S t h h ital images, the Laplacian operator becomes = +

, tied to the horizontal and the vertical second derivatives (called and directions, resp.). This 2D operator allows approximating the second derivative of the image, in order to detect edges by zero-crossing [START_REF] Marr | Theory of edge detection[END_REF] and highlight ridges or valleys. Unfortunately, it will also increase noise. So, a good practice is, before applying the Laplacian, to smooth the signal by convolving the signal with Gaussian as preprocessing step, which consequently will not only reduce the noise but also avoid producing ringing artifacts that yields misclassified ridges and valleys.

In the next section, first, the mathematical definition of ridge has been detailed. Then, an extensive theoretical review of the most commonly used ridge detection and extraction filtering techniques is performed. Additionally, extensive analysis of ridge characterizations, driving parameters, and its usefulness with different filtering techniques have been examined. Section 3 is devoted to objective experimental evaluations of the ridge detection and extraction filtering techniques on both synthetic and real images. Furthermore, the evaluation of the ridge filtering techniques is explored on the context of common noises in the images. At last, the main conclusions have been inferred in Section 4.

RIDGE EXTRACTION IN IMAGES

Initially, a discrete definition for ridge appears in [START_REF] Crowley | A representation for shape based on peaks and ridges in the difference of low-pass transform[END_REF], where the underlying function is the image convolved with a difference of low-pass (DOLP) transform. Considering two different low-pass filters and (i.e., two supports of different widths) both positioned over the center coefficient at the point (0, 0), ridges, valleys, and blobs may be extracted efficiently with the DOLP transform: these features are high- lighted by applying two different low-pass filters to the same image and then subsequently subtracting these two filtered images. Note that the difference of the filters may be applied before convolving the image with the obtained DOLP filter. Afterwards, crest lines are extracted when the support of the low-pass filter is smaller than the support of the low-pass filter and inversely regarding valleys. The final step of the ridge extraction consists in the suppression of the local nonmaxima of the magnitude of the ridge/valley intensity, finally the image is thresholded in four directions to obtain thin nets of pixels. Regarding DOLP filter, the selected pixels correspond to points being local maxima in one of the 4 orientations (modulo 180°) associated with the 8-neighborhood of the pixels. Figure 5 illustrates this process. Even though the results obtained with square shapes are acceptable (as illustrated in Figs. 3b and3c), the DOLP filter formed by subtracting circularly low-pass filters is preferable. Nevertheless, for their isotropy and circular symmetry properties, the sampled Gaussian filter represents a good achievement. Indeed, the Difference of Gaussians (

Δ ∂ ∂ 2 2 x ∂ ∂
) remains effective in ridge detection and is an approximation of the Laplacian of Gaussian ( ) when the ratio of the size filters is roughly equal to 1.6 1 L 2 L DoG LoG [START_REF] Marr | Theory of edge detection[END_REF]. Usually called Mexican hat or Sombrero filter, the 2D equation of the is given by:

(

where represents the pixel coordinates and is the standard deviation of the Gaussian. A discrete is presented in Fig. 4d and a ridge extraction example in Fig. 4c after a non-maxima suppression in the 4 orientations associated with the 8-neighborhood of the pixels (same process than with DOLP filter, see Fig. 5). Also, in Haralick's approach [START_REF] Haralick | Ridges and valleys on digital images[END_REF], the image function is approximated by a cubic polynomial which, sometimes, may distort the detection.

DOLP transform and allow extracting roughly ridges and valleys, but suffer when the desired objects are too thin, thus the detection is disturbed by noise or undesirable artifacts. Besides, the angle selectivity may be improved by applying other operators, as presented in the following section.

Hessian Matrix

In image filtering, the second order derivative may be used to determine the location of the ridges.

LoG

( ) Indeed, bright or dark ridges correspond to, respectively, a maximum or minimum of the image intensity in the direction orthogonal to them and a constant image intensity in the direction parallel to them. Considering a grey level image and its partial derivatives:

+ -σ   + = -   πσ σ   2
• , the 2nd image derivative along the axis (see Fig. 6b);

•

, the 2nd image derivative along the axis (see Fig. 6c);

•

, the crossing derivative of (see Fig. 6d), the Hessian matrix is often computed in image analysis:

(3) Image derivatives can be calculated by convolving the image with the ±[-1 0 2 0 -1] or the ±[-1 0 1] masks in the and/or directions. Note that ±[-1 0 2 0 -1] ±[-1 0 1] [-1 0 1] in the discrete domain.

The matrix is symmetric, diagonalizing provides the local normal to the ridge or the valley (that is given by the eigenvector related with the highest eigenvalue) and its sharpness (that is related to the values of these eigenvalues) [START_REF] Eberly | Ridges for image analysis[END_REF][START_REF] Steger | An unbiased detector of curvilinear structures[END_REF]. Theoretically, eigenvalues are computed by:

I = ∂ ∂ 2 2 / xx I I x x = ∂ ∂ 2 2 / yy I I y y = ∂ ∂ ∂ 2 / xy I I x y I * ( ) ( ) ( ) ( ) ( )     = =         * * * * *
, k k [START_REF] Crowley | A representation for shape based on peaks and ridges in the difference of low-pass transform[END_REF] they are visible in Figs. 6e, 6f. Then, eigenvectors, tied to the direction perpendicular to the ridge/valley, are given by: [START_REF] Deriche | Recursively implementing the Gaussian and its derivatives[END_REF] The two eigenvalues and correspond to the two main curvatures of the local surface. Besides, there exists several functions indicating the local image contrast [START_REF] Tremblais | Algorithms for the extraction of various diameter vessels[END_REF]:

• , corresponding to the main eigenvalue [START_REF] Steger | An unbiased detector of curvilinear structures[END_REF],

• , see [START_REF] Lindeberg | Edge detection and ridge detection with automatic scale selection[END_REF],

• , see [START_REF] Lindeberg | Edge detection and ridge detection with automatic scale selection[END_REF],

• , see [START_REF] Tremblais | Algorithms for the extraction of various diameter vessels[END_REF]. These functions are plotted in Fig. 8 as a function of and .

Eventually, a pixel is labeled as a ridge/valley pixel if is maximum in the direction. It is

( ) ( ) ( )
( )  = +   - + +     = +   + + +   1 1
 θ =   -   * * 12 1 1 1 . k 1 k 2 k { } ∈ , 1,2,3,4 i i D = 1 1 D k = + 2 2 2 1 2 D k k ( ) = + 2 2 2 3 1 2 D k k = -⋅ + 4 1 2 1 2 D k k k k 1 k 2 k { } ∈ , 1,2,3,4 i i D θ
Fig. 6. Image derivatives and eigenvalue images using Hessian matrix (Eq. ( 3)) or Weingarten W (Eq. ( 6)). selected after non-maximum suppression [START_REF] Canny | A computational approach to edge detection[END_REF] where the values of the magnitude are linearly interpolated between closest points in the 8-neighborhood. Figure 9 illustrates a ridge detection after non-maxima suppression in the direction.

In practice, regarding real images, due to the luminance variation, acquisition and/or compression noise, the detection of pure ridges/valleys is almost impossible. So, in order to more reliably extract the ridges, the convolution of the image with a low-pass filter is considered, as detailed in Section 2.3.

Weingarten

Weingarten map represents the differential of the Gauss map [START_REF] Do | Differential Geometry of Curves and Surfaces: Revised and Updated Second Edition[END_REF]. This expression can be computed directly from the first (i.e., and ) and second derivatives in the and directions of the images. The linear invariants of the Weingarten map are the intrinsic curvatures of the surface: the eigenvalues are the principal curvatures, the trace is the mean curvature, and the determinant is the Gaussian curvature: [START_REF] Do | Differential Geometry of Curves and Surfaces: Revised and Updated Second Edition[END_REF] The eigenvalues and eigenvectors of are extracted with the same procedure as in Eqs. ( 4) and [START_REF] Deriche | Recursively implementing the Gaussian and its derivatives[END_REF], regarding coefficients of the matrix . The same procedure applies for the non-maxima suppression in the direction. In [START_REF] Armande | Thin nets extraction using a multi-scale approach[END_REF], ridges or valleys are extracted θ

∂ = ∂ x I I x ∂ = ∂ y I I y x y ( )   + -   = ⋅ ⋅       - +     + + 2 3 2 2 2 2 1 1 ( , )
.

Low Pass Filters for Ridge Detection

A low pass filter is the basis for most smoothing methods, as it is needed to apply smoothing process over an image prior to the feature extraction. The optimization criteria, based on the Canny theory, are: (i) detection efficiency, (ii) location accuracy of the detected contour, and (iii) uniqueness condition of filter response to its output for an input signal [START_REF] Canny | A computational approach to edge detection[END_REF]. Based on this theory, several low pass filters have been proposed in the literature. In the following, three lowpass filters and their second derivatives are discussed for ridge and valley detection.

2.3.1. Ziou filter Z. In [START_REF] Ziou | Optimal line detector[END_REF], the author described an optimal line detector allowing an economic temporal complexity because it represents a second order recursive filter. Considering , the equation of the 1D low pass filter is given by [START_REF] Eberly | Ridges for image analysis[END_REF] where the filter parameter represents a positive constant. The second derivative of is obtained by derivation as a function of , two times: [START_REF] Freeman | The design and use of steerable filters[END_REF] Note that the same procedure is available to obtain the 1st derivative of the filter Z, as for the following presented filters.

D ∈ R t Z ( ) ( ) -⋅ = + 2 1 1 , z s t z z Z t s t e s z s Z t ( ) ( ) -⋅ = -1 . z s t z z t
s t e 2.3.2. Gouton filter R. Gouton et al. [START_REF] Gouton | Ridge-line optimal detector[END_REF] described a third order recursive filter. This ridge/valley-line detector function is able to modify its shape as a function of its parameter : [START_REF] Getreuer | A survey of Gaussian convolution algorithms[END_REF] with: , , and Thus, the second derivative of has the following form:

r s ( ) ( ) ( ) ( ) - = + + sin cos , r s t r r R t K s t D s t E e ( ) ( ) - - = + 2 2 2 2 2 2 r r r r s s t A s s t = 4 1 4 r K s = 2 4 2 4 r r s A D s + = 3 . r r r As s E s R ( 10
)
the more the parameter decreases, the more enhances fine ridges/valleys. Furthermore, when decreases, the shape of is nearly a Gaussian, as shown in Fig. 10.

Gaussian filter.

Gaussian kernels are regularly used for their effectiveness in edge detection [START_REF] Canny | A computational approach to edge detection[END_REF], the 1D equation is:

(11) ( ) ( ) ( ( ) ( )) - = - - + cos sin 1 , r r s t r r r r t s t s s t s e r s r r s R ( ) -σ σ = πσ
D 1 D 2 D 3 |k 1 | |k 2 | |k 1 | |k 2 | |k 2 | |k 1 | (a) (b) (c) D 1 function D 2 function D 3 function 1 D 2 D 3 D ( ) 1 2
, k k with representing the standard deviation of the Gaussian. Thus, the second derivative of the Gaussian is given by: (

The two dimensional Gaussian is built by combining with " ," a product of convolution, horizontally and vertically: [START_REF] Iverson | Logical/linear operators for image curves[END_REF] where represents transpose. Using the Gaussian , the strategy is the same as to compute the second derivative on an image, with and , as an example for an image derivative in (see Fig. 11d). Moreover,

σ σ G ( ) ( ) - σ σ σ ∂ -σ = = ∂ π σ 2 2 2 2 2 2 2 5 . 2 t G t g t t e t σ G * σ G ( ) ( ) + -σ σ σ σ = * = πσ 2 2 2 2 2 1 , , , 2 x y x y G G x y e G σ G σ G σ g y
Section 2.4.4 is dedicated to the strategies of the twodimensional filters implementation.

Furthermore, these filters in Eqs. ( 7), [START_REF] Geusebroek | Fast anisotropic Gauss filtering[END_REF], and ( 12) are useful to smooth the image before extracting edges by computing matrix presented in Eq. ( 3). Additionally, it is also possible to use the Weingarten (cf. Eq. ( 6)) with the Gaussian, as in [START_REF] Armande | Thin nets extraction using a multi-scale approach[END_REF].

2.3.4. Parameters. The three above-mentioned filters are suitable for ridge and valley detection. Considering one filter, it is adjustable by tuning only one parameter which is the same for the low pass and the derivative filter. Accordingly, parameters , , and are chosen as a function of the width of the ridge or of the valley. Consequently, these parameters are selected by increasing the width of the filter as robust as possible in order to extract suitably the feature. Here, the main idea is to compare equivalently the 3 filters , , and as a function of the feature width. Indeed, the objective is to tune each filter for a specific width by selecting an appropriate parameter ( , , or ). Thus, in the discrete domain, parameters and are decreasing, and parameter is increasing until the filter coefficients cross 0 and the shape filter contains the width of the feature at the same time. Figure 10 illustrates the selected filters computed with different parameters as a function of the width of the feature. In addition, Table 1 references the optimum parameters for each filter as a function of the features size from 1 to 15 pixels. Finally, the parameter of the Gaussian has the same properties regarding oriented filters widths.

Oriented Filters

One common task in image processing and computer vision is applying the same filter on different angles in order to detect directional responses as steerable filter [START_REF] Freeman | The design and use of steerable filters[END_REF][START_REF] Jacob | Design of steerable filters for feature detection using canny-like criteria[END_REF], anisotropic Gaussian kernel [START_REF] Lopez-Molina | Unsupervised ridge detection using second order anisotropic Gaussian kernels[END_REF], and logical linear filter [START_REF] Iverson | Logical/linear operators for image curves[END_REF].

Steerable filter. Two-dimensional Gaussian kernels

are very useful for their properties of isotropy, steerability or decomposability properties as separability (see Section 2.4.4). Freeman and Adelson proposed an efficient architecture to design oriented filters of arbitrary orientations from linear combinations of basis filters [START_REF] Freeman | The design and use of steerable filters[END_REF]. Thus, applying filter steered in different directions, and then computing the filter responses gives significant description of the orientation for the considered pixel; also, it allows determining analytically the filter output as a function of orientation. The initial step to extract ridges or valleys in images is to estimate their orientation by using even steered filters. Conse-

z s r s σ z s r s σ σ σ G
quently, the steerable filter of second order ( ) considering in Eq. ( 13) oriented at the angle is: [START_REF] Jacob | Design of steerable filters for feature detection using canny-like criteria[END_REF] This allows computing an even filter at a specific orientation, as illustrated in Fig. 11e. At the end, the calculation of the ridges or the valleys corresponds to the filter energy in the direction of the maximum response of the template.

Jacob and Unser [START_REF] Jacob | Design of steerable filters for feature detection using canny-like criteria[END_REF] extend the idea of the steerable filter of order 2 ( ) with operators having a better orientation selectivity. Indeed, they proposed higher order functions, issued of higher order derivatives of the Gaussian (2nd and 4th: , , , ,

), resulting in more elongated templates, as visible in Fig. 11f. Regarding ridge detection, this filter is specified so as to provide the best compromise in terms of signal-to-noise ratio, false detection, and localization (as illustrated in Fig. 7d). Thus, the even steerable filter of 4th order ( ) is formulated as [START_REF] Laligant | Noise estimation from digital step-model signal[END_REF] with • , • ,

• ,

• ,

• . Such parameters control the template not to produce undesirable oscillations and side-lobes along which is contrary to the 3rd Canny criterion: unicity (see [START_REF] Jacob | Design of steerable filters for feature detection using canny-like criteria[END_REF]). This 2D template, presented in Fig. 11f, can be steered in different orientations , as detailed in [START_REF] Jacob | Design of steerable filters for feature detection using canny-like criteria[END_REF], to extract ridges and valleys.

SF σ G θ ( ) ( ) ( ) ( ) θ σ σ σ ∂ = θ ∂ ∂ ∂ + θ θ + θ ∂ ∂ ∂ 2 2 2 2 2 2 2 2
cos cos sin sin .

SF x x y y

G G G 2 SF σ G yy G xx G yyyy G xxyy G xxxx G 4 SF ( ) = α + α + α + α + α 4 1 2 3 4 5
, , 1. The negative of filters detect ridges.

yy xx yyyy xxyy xxxx SF x y G G G G G α = - σ 1 0.392 α = σ
(a [START_REF] Perona | Steerable-scalable kernels for edge detection and junction analysis[END_REF]. As an example, as it corresponds to a thin filter, the orientation selectivity becomes more reliable with the anisotropic Gaussian derivatives than the isotropic Gaussian filter. Moreover, the isotropic Gaussian often makes the ridge extraction difficult as crossing lines in images; also, parallel lines could be blurred into one line due to isotropic smoothing, especially if the smoothing parameter is too large (i.e., parameter in Eqs. ( 11) and ( 12)). The origin is that anisotropic property is more efficient at level of straight lines. It corresponds to a narrow filter which is oriented in different directions to extract the edges when it is steered in the edge direction. The anisotropic Gaussian filter can thus take advantage of this property and overtake drawbacks of the isotropic filter. Thereafter, it is necessary to filter the image with a set of kernels oriented in a variety of directions, as such, leading to the characterization of the partial derivatives in different orientations. The most evident option to produce a single output from that information is to retain the result produced by the oriented kernel with the maximum absolute value. An anisotropic Gaussian filter in two dimensions is not built with the combination of isotropic kernels, as Eq. ( 15), its direct equation is given by [START_REF] Geusebroek | Fast anisotropic Gauss filtering[END_REF]: [START_REF] Laligant | Regularization preserving localization of close edges[END_REF] Here, represent the two parameters of the anisotropic Gaussian, i.e., the standard deviations. When , the kernel reduces to an isotropic Gaussian kernel or . To extract ridges, the second-order anisotropic Gaussian kernel (SOAGK) can be applied [START_REF] Lopez-Molina | Unsupervised ridge detection using second order anisotropic Gaussian kernels[END_REF]. Considering the vertical anisotropic Gaussian directed at :

) (b) (c) (d) (e) (f) [1 0 2 0 1] Z u z T R u r T SF 2 = 30q G V u g V T T SF 4 = 0q
thus, the second derivative of in the direction is calculated by: [START_REF] Lindeberg | Edge detection and ridge detection with automatic scale selection[END_REF] The choice of enables to build a narrow filter smoothing mostly in the direction while high-σ To this end, the anisotropic parameter produces a smoothing along the ridge/valley, which helps to extract easily elongated features, even disturbed by noise. On the contrary, kernels having parameters highlight undesirable features as noise which are interpreted as small, nonelongated ridges [START_REF] Lopez-Molina | Unsupervised ridge detection using second order anisotropic Gaussian kernels[END_REF].

°360 /Δθ °360 /Δθ ( ) ( ) σ σ θ   θ+ θ - θ+ θ - +     σ σ   = πσ σ
v σ ,σ u σ = σ v u & σ v G σ u G θ = 0 ( )   - +     σ σ   σ σ θ= = πσ σ 2 2 2 2 1 2 , , 0 1 , , 2 u u y x u x y e v v v & σ σ θ= v & , , 0 u x ( ) ( ) σ σ θ= σ σ θ=   - +     σ σ   ∂ = ∂ -σ = πσ σ
2.4.3. Logical linear filter. Similar to the SOAGK, Iverson and Zucker proposed a hybrid filter by combining directional linear filters and a linear-logical (L/L) operator which helps to reduce the false positive pixels of ridges/valleys [START_REF] Iverson | Logical/linear operators for image curves[END_REF]. This technique allows selecting any inflection points within the 1D signal region , with (see Fig. 14). First, it depends on the Gaussian (see Eq. ( 11)) and its derivatives of the first and third order and by computing the four parameters: [START_REF] Lopez-Molina | Unsupervised ridge detection using second order anisotropic Gaussian kernels[END_REF] thereby, they can be applied to a signal, as shown in Fig. 14.

Concretely, the L/L operator can be utilized on different edge types as ridge (positive contrast lines), valleys (negative) and edges (ramp or step). In this study, only and are focused and evaluated. These denote functions and , respectively, and combine linear operators in Eq. ( 18) by using the logical operator such that: [START_REF] Magnier | Edge detection evaluation: A new normalized figure of merit[END_REF] where the logical operator is represented by, for two hypotheses :

In this way, and contribute to extract convex and concave points, as shown in Fig. 14. Next, to extract ridges or valleys and their tied directions, the 2D operator is expressed as the Cartesian product of orthogonal, 1D L/L operators or , and a tangential operator T(t). Moreover, this 2D operator is oriented and uses strategies of the logical operator with the tangential operator T(t) to (a) discriminate between locally continuous and discontinuous curves along their tangent direction in the image and (b) align the line termination with the line ending (illustrated in

x σ σ ≈ v v / 1 [ ] - + e e , t t > e 0 σ G σ ' G ( ) σ 3 G ( ) ( ) ( ) ( ) ( ) ( ) σ σ σ σ  = + = -   = + = -   3 3 (3) (3)
£ - = = - - - (3) (3) (3) (3)
' ' , ' ' , 13, for more details refer to [START_REF] Iverson | Logical/linear operators for image curves[END_REF]). To sum up, the L/L operator is similar to the SOAGK, with the parameter tied to the normal operator ( and ) and for the tangential operator (see Fig. 12a).

P l r l r N l r l r n n n n n n n n 6 6 £ £ £ £ £ £ £ ( ) , a b + > ∧ >   > ∧ ≤   ≤ ∧ >   + ≤ ∧ ≤  if 0 0, if 0 0, if 0 0, if 0 0. a b a b b a b a b a a b a b a b £ 6 P 6 N 6 P 6 N £ Fig.

Implementation and complexity.

Presented filters , , and may be implemented with different

σ u 6 P 6 N σ v T Z R σ G
strategies. Filters and have been designed to be implemented recursively. Also known as infinite impulse response (IIR) filters, they represent filters where the output sample is a linear combination of some number of previous inputs and outputs. Even though there are other ways to apply them (i.e., via fast Fourier transform), the recursive implementation Z R

3 l n ( ) 3 r n = e 2 6 P 6 N
strategy is compared here. First, they correspond to separable filters, so they can be written as a product of two 1D filters, which is equivalent to a typical 2D convolution operation providing a reduction of the operator computational cost. Secondly, to reduce even more the number of operations per pixel, each 1D , , or filter may be also implemented recursively, Table 2 reports the order of these 3 detailed filters. To reduce the number of operations per pixel, an -order recursive filter is obtained by calculating its Z-transform. Thus, the two-sided sequence of a filter is the superposition of a causal filter and anti-causal filter : = + , for . To minimize the computational complexity, the authors of [START_REF] Van Vliet | Recursive Gaussian derivative filters[END_REF] proposed to decompose series interconnection into a product of the causal and anti-causal parts, leading to a 3rd-order Gaussian filter, a 4th-order first derivative filter and a 5th-order second derivative filter (many fast approximations of the Gaussian have been proposed, some of them are detailed in [START_REF] Getreuer | A survey of Gaussian convolution algorithms[END_REF]). Now, the first and second derivatives of the original image can be computed easily by applying the [-1 0 1] mask one or two times, respectively, to the smoothed image (i.e., smoothed image obtained by applying the , , or filters both in and directions). Besides, the derivatives of an image are computed by combining in the two directions and the different 1-dimensional filters presented in the beginning of Section 2.3 (which are implemented with the different strategies above). As an example, the second derivative in the direction of an image with the filter can be obtained by applying the low pass filter in the y direction then the second derivative of , called in the direction of the filtered image.

Subsequently, Table 3 specifies the required number of image computations as a function of the segmentation technique ( , , , , , or L/L) and Fig. 15 roughly schematizes the complexity. The Hessian Matrix needs the second derivatives of the image , , and , using , , or filters. Obviously, is more computationally complex than the or because it needs more image derivatives. Regarding the steerable filters, an operation of filter rotation with an angle is necessary (with 360°/ total rotations, where is the angular step); and 5 derivative images are calculated for the steerable filter of order 4 (see Eq. ( 15)). On the other hand, the number of basis filters is large extract features with the SOAGK, and the basis filters are nonseparable, requiring high computational loads. In [START_REF] Geusebroek | Fast anisotropic Gauss filtering[END_REF], the anisotropic Gaussian is decomposed into two Gaussian 1D filters by considering 360°/ steps of rotation, allowing reducing the operation number per pixel (to approximate the SOAGK, the difference of anisotropic Gaussian by differentiating the whole image array

Z R σ G M F - F + F ( ) F n ( ) - F n ( ) + F n { } = 1,..., n M Z R σ G x y x y x Z Z z Z x LoG * W 2 SF 4 SF * xx I yy I xy I Z R G W LoG * θ Δθ Δθ
Δθ with two different standard deviations in Eq. ( 16) is calculated, see [START_REF] Geusebroek | Fast anisotropic Gauss filtering[END_REF]). Also, the L/L filter contains several steps of interpolation for the normal operator ( and ) and for the tangential operator which are directed in different directions in the image. Moreover, the L/L uses other strategies such as the endline or the stabilizer to qualify the segmentation; these steps add more filter complexity.

EXPERIMENTAL RESULTS AND EVALUATION

Experiments are performed on synthetic and real images, showing qualitative and quantitative results. A first result presented in Fig. 16 illustrates the advantage to use sharp and narrow filter to extract thin and close objects, as filters and .

The aim here is to extract branches inside the dragonfly wings; as this image does not contain any noticeable noise, the Hessian matrix with finite filters like [1 0 -2 0 1] gave interesting results for these thin objects, but created many undesirable edge points around certain valleys (similar segmentation also by

). Elsewhere, with the Gaussian and brings similar but less complete result. Segmentation obtained with and , , and are worse with a lot of missing edge points, as with . However, the valley extraction using is perfectible. On the other hand, the result using with filter is quite perfect (Fig. 16b); this justifies the need to use low pass filter. Among all the ridge/valley detectors, exponential ( or ) filters do not delocalize contour points [START_REF] Laligant | Regularization preserving localization of close edges[END_REF], whereas they are sensitive to noise. Techniques using Gaussian filters are less sensitive to noise, but suffer from rounding bends and junctions like the oriented filters , , and the SOAGK. The more the 2D filter is elongated, the more the segmentation remains robust against noise. In the following sections, quantitative results are reported with different types and levels of noise in synthetic images. Then, evaluations will involve real images.

Error Quantification and Evaluation Procedure

Evaluations are reported using synthetic images where the true positions of the edges are known. Let be the reference contour map corresponding to the ground truth and the detected contour map of an Recursive order 2, see [START_REF] Ziou | Optimal line detector[END_REF] 3, see [START_REF] Gouton | Ridge-line optimal detector[END_REF] 4, see [START_REF] Deriche | Recursively implementing the Gaussian and its derivatives[END_REF] or 5, see [START_REF] Van Vliet | Recursive Gaussian derivative filters[END_REF] image . Comparing pixel by pixel and , a basic evaluation is composed of statistics:

σ u 6 P 6 N T z r H 4 SF H σ g 1 D * 2 D 3 D 4 D 2 SF W * z z r 2 SF 4 SF t G c D
• true positive ( ), common points of both and ;

• false positive ( ), spurious detected edges of ;

• false negative ( ), missing boundary points of ;

• true negative ( ), common non-edge points, where denotes the cardinality of a set. Several edge detection evaluation techniques involving only statistics have been developed (see [START_REF] Magnier | A review of supervised edge detection evaluation methods and an objective comparison of filtering gradi-ent computations using hysteresis thresholds[END_REF]). It is clearly proved that poorly located or missing pixels should be penalized according to the distance from the position where they should be localized. Also, as demonstrated in [START_REF] Magnier | A review of supervised edge detection evaluation methods and an objective comparison of filtering gradi-ent computations using hysteresis thresholds[END_REF], the evaluation of FP and FN should not be symmetrical, because such a penalty could alter the visibility of the outlines of the desired objects in an objective evaluation (see [START_REF] Magnier | A review of supervised edge detection evaluation methods and an objective comparison of filtering gradi-ent computations using hysteresis thresholds[END_REF]): some measures calculate a large error for a single FP at a sufficiently large distance, while many desired contours are missing, but unfortunately, they are not penalized enough. Thus, described

I t G c D TP t G c D FP c D FN c D TN ⋅
in [START_REF] Magnier | Edge detection evaluation: A new normalized figure of merit[END_REF], the normalized edge detection evaluation measure is, for or :

where represent two scale parameters [START_REF] Magnier | Edge detection evaluation: A new normalized figure of merit[END_REF], denotes the cardinality of a set, and is the minimal Euclidian distance between a pixel and a set [START_REF] Magnier | A review of supervised edge detection evaluation methods and an objective comparison of filtering gradi-ent computations using hysteresis thresholds[END_REF]. Therefore, the measure calculates a standardized dissimilarity score; the closer the evaluation score is to 1, the more the segmentation is qualified as suitable. On the contrary, a score closes to 0 corresponds to a poor detection of contours.

The aim here is to get the best contour map in a supervised way. For that, the contours are extracted after a suppression of the local non-maxima, then a threshold by hysteresis is applied to obtain a binary segmentation [START_REF] Canny | A computational approach to edge detection[END_REF]. Theoretically, to be objectively compared, the ideal contour map of a measure must be a at which the supervised evaluation gets the highest score [START_REF] Magnier | Edge detection evaluation: A new normalized figure of merit[END_REF][START_REF] Magnier | A review of supervised edge detection evaluation methods and an objective comparison of filtering gradi-ent computations using hysteresis thresholds[END_REF]. For each better segmentation tied to , the FP and percentage of TP relative to the total number of edge pixels of are also displayed (TP/ ). In addition, the last evaluation measure concerns the angle tied to the ridge/valley, . Considering , the set of contour chains in (i.e., at least 2 pixels per chain), the angle evaluation is computed as follows:

where represents a contour pixel belonging to , a 3 × 3 window centered on , is the direction tied to , and is the number of contour pixels in minus the central pixel. This evaluation linearly ranges from 0 for identical angles of and to 1 for angles that differ. Note that and angles belong to [0°; 180°[ and when one direction approximates 0° and the other Complexity schema, depending on the recursive filters order, the number of calculated images (Table 3) and the filter rotation.

N > 0 FN > 0 FP ( ) ( ) ( ) ∈ ∈
= + δ + + κ   2 2 1 , 1 1 , 1 c t t c t c c p D G t p G D FP N G D D d p FN G d p ( ) δ κ ∈ 2 , ]0,1] ⋅ ( ) A d p p A N c D N t G t G θ c D C c D ( ) ∈ ∈ ω   °-θ -θ θ = -   °      90 1 , 1 / , 90 k c c D k c p d D k D p C d E C c C k d ω p θ k d k d k c ω θ p θ k d θ k d θ p
H with Z H with R H with G V W with G V SF 2 SF 4 SOAGK L/L

Complexity

Ridge/valley detection technique direction approximates 180°, the evaluation remains close to 0. Also, from proper binary confusion matrix, the precision ( ) and recall ( ) evaluations are computed, given the overall quality expressed in terms of the -measure: with allows an equal penalization between and . These scores are presented throughout the remainder of this study, according to different images and noise types.

Synthetic Images Corrupted by Poisson Noise

The first image in Fig. 17 contains ridges of width 1 pixel and is corrupt by Poisson noise. This noise distribution typically models shot noise in a sensor in which the time between photon arrivals is governed by Poisson statistics and appears for example in medical imaging system. The Poisson noise density models the time statistics between photon arrivals, where the defined expected number of occurrence of integer , the noise probability corresponds to observed number of instead of in pixel [START_REF] Zhang | PCAbased spatially adaptive denoising of CFA images for single-sensor digital cameras[END_REF] . Given an integer, the maximum probability is obtained for and the variance of the distribution is also ; at a pixel , the

P R F P R TP TP P R TP FP TP FN α = 0.5 FN FP Λ n λ k λ = λ t λ x
equation is given by = . Poisson noise appears from quantum effects of photons. The Poisson noise estimate is the square root of the number of detected photons, then the Cameras cab be calibrated to know how many photons a certain pixel value is [START_REF] Akiyama | Pseudo four-channel image denoising for noisy CFA raw data[END_REF].

As shown in Fig. 17, except SOAGK and L/L, all the other filters are robust to Poisson noise at this scale.

performed exceptionally well, with True Positive value (and ) almost close to 1. It only started to drop from SNR of 5 dB. Filters with , , , and ranked in the second place, their performances are still sufficiently well (TP/ over 0.9). But seemed to be less robust to Poisson noise, whose TP dropped sharply from SNR of 6 dB. In comparison, SOAGK showed its relatively poor performance to resist the noise-starting with TP 0.85, ending with TP 0.3 at SNR of 3 dB. L/L failed completely in this task. It detected barely any true positive ridges. This poor performance of oriented filters is caused by the small size of these filters where small-scale orientation deforms the kernels. The angular score is the best for and filters (using ), because they correspond to sharp filters, especially suitable for thin ridges.

Synthetic Images Corrupted by Speckle Noise

The second image in Fig. 18 contains ridges of width 3 pixels and is disturbed by a speckle noise. This multiplicative noise appears with the image acquisi-( ) tion due to the level of noise in the sensor of a CCD or CMOS camera, increasing in proportion to luminosity [START_REF] Laligant | Noise estimation from digital step-model signal[END_REF]. This noise model can be formulated as: = I + σηI, where represents the observed image, is the noise-free image, and is a normalized Gaussian noise distribution centered at 0 of standard deviation . Compared to Poisson noise, which is correlated to the original image, speckle noise adds some independent noise to the images that could corrupt more the image's geometric structure. For filters providing quite good results in previous situation (Fig. 17), they are less efficient in case of speckle noise. Instead of starting with TP in the range of [0.93, 1], they are now under the threshold of 0.85 corresponding a performance drop of 10% at SNR of 9 dB. And the robustness to noise level decreased much more. This decreasing behavior in the interval [0.85, 0.55] is similar to speckle noise's granular effect property. At SNR of 3 dB, the TP are around 0.55, meaning only about 55% ridge pixels are perfectly detected. Unlike the total failure with Poisson noise, L/L filter worked correctly. However, its performance is still the worst compared to other filters. The main reason could be due to the fact that L/L is by definition a 1D filter with additional processing as the endline or the stabilizer which are noise sensitive. This makes it be much less robust to structure-correlated noise. As the structurecorrelated noise could destroy the two-dimensional visual structures transformed in 1D filter space and cause thus the failure of detection. On the other hand, among the techniques using non-oriented filters, it is noticeable that with and with obtain best scores. Additionally, the extracted ridges are more continuous and less disturbed by undesirable false positive pixels. with allows a better quality of

λ x P t -λ λ ! t e t 2 SF / t TP G α F Z R σ G 4 SF α F 4 SF E Z R H
J J I η σ H 1 D W 1 D W 1 D
detected ridges than with other non-oriented filters.

Finally, the angular score obtained by the SOAGK is less penalized because it corresponds to an elongated kernel applied on close-right structures, so the oriented filter is generally the same along these structures (same remark for the Gaussian noise).

Synthetic Images Corrupted by Gaussian Noise

The last experiment with synthetic images in Fig. 19 concerns valleys corrupted by a Gaussian blur and Gaussian noise. This type of noise represents an additive noise disturbing gray values in images. Its model is essentially designed and characterized by its probability distribution function (PDF) or histogram normalization with reference to the gray value: = where is the gray value, is the mean, and is the standard deviation. Gaussian noise blurred more the geometrical structure in images. So, as shown in Fig. 19, the general TP for all filters are decreased compared to those with Poisson and Speckle noises. SOAGK, , filters gave better results (TP 0.8). L/L filter showed always the worst result, even at SNR of 9 dB, the TP is only 0.57; when noise becomes stronger, performances decreases. In comparison, with still detects a better quality of ridges than with other non-oriented filters, statistically and visually.

Evaluation with Real Images

After evaluating the filters on synthetic images with different types of noise, the ridge detection on real- world images is presented. These images are from the Ghent University Fungal Images together with their manually annotated ground-truth ridges [START_REF] Lopez-Molina | Unsupervised ridge detection using second order anisotropic Gaussian kernels[END_REF]. This database is extremely challenging. Here, 13 images with their tied ground-truth images are selected randomly for this experiment. The images have very poor contrast and strong noises (see Fig. 20). Regarding the evaluation pixel per pixel, due to the hand-labeled ridge points which create inaccurate ground truth ( ), the overall ridge detection with these filters is around TP 0.2, and they are image-dependent. In the best situation, TP can reach 0.3; otherwise, in worst cases, the TP will drop below 0.5 and be close to 0.06. Oriented filter , , and SOAGK performs well, regarding and , contrary to the L/L and with . Regarding with , its evaluation is better than other non-oriented filters, even though the angle evaluation penalizes the directions perpendicular to the detected ridges (however, the score remains under 0.1 where it was under 0.2 for speckle noise).

Visual Comparison on Real Images

The experiment on real images have been performed on heart images for cartographic analysis and aerial noisy image which are available in Figs. 21 and22, respectively.

For the first image presented in Fig. 21a, the output filtered images demonstrate interesting results with all filtering techniques reviewed in this literature nonmaxima suppression and thresholded images. Here, the thresholded images correspond to the 3500 highest computed points for each filtering techniques. The original image is not specifically noisy but it is blurred.

G 2 SF 4 SF α F N H 4 D W 1 D E
The selected parameters correspond to the detection of ridges of width 3 pixels (see Table 1). Usually, the blood vessels are well extracted, even with the Laplacian, even though ridges are roughly extracted (see Fig. 21b). The discrete filter obtains the worse result with some false positive points and extracted ridges which are doubled. Regarding the Hessian matrix, , , and perform well using with the extraction of main blood vessels without many false positives. Same remark, the Weingarten is reliable contrary to the Hessian matrix with and / , which are corrupted by many false positive pixels and noncontinuous extracted ridges. Now, considering oriented filters in Figs. 21j-21l, performs as well as , , and filters using while obtains the stronger result with continuous extracted ridges without many undesirable extracted pixels. Elsewhere, the SOAGK in Fig. 21l performs well for elongated ridges but too thin blood vessel are not well extracted, while some extracted ridges are tripled, penalizing this ridge detector.

Likewise, the same implementation procedure has been applied on Stanwick aerial image in Fig. 22a for aerial image analysis, where the original image is corrupted by a strong noise. Because the images contain too many ridges to be detected, the thresholded images after filtering correspond to 50% of the highest extracted points after non-maxima suppression (excepted for the discrete filter). The selected parameters correspond to the detection of ridges of width 5 pixels (see Table 1). The extracted information in Figs. 22b-22l, clearly shows that each filtering technique detects different level of information after nonmaximum suppression, as subject to different application analysis. For instance, the LoG (Fig. 22b

) detects Z R σ G 1 D σ G 2 D 4 D 2 SF Z R σ G 1 D 4 SF Fig.

Z R

The Hessian matrix ( ) using Gaussian and in Fig. 22f shows more interesting result with extracting sufficiently ridges mostly connected without need of any further processing even though it contains many false positive points. The results obtained by and / in Figs. 22g and22h show more ridges with sharp noises. The Weingarten in Fig. 22i mathematical properties, driving filter parameters and characterizations. Classical technique s such as DOLP and Laplacian filter are detailed, along with the elegant way using the Hessian matrix . Different techniques exist to compute ridge or valley extraction with the eigenvalues of ; they are detailed through this communication, as for the Weingarten. Three low pass filter are compared, namely, , , which are exponential filters, and the Gaussian one . Elsewhere, method using oriented Gaussian filters are also reported. The evaluation and comparison of filtering techniques have been performed both theoretically and experimentally on synthetic and real images. Each filtering technique has been examined on complex images, where different types of noises have been applied. The acquired comparison and evaluation graphs exhibited which method is reliable as a function of the width feature and the noise type.

Regarding non-oriented filters, the filter performs well when the ridge or the valley are very thin (width of one pixel) and requires the less computational complexity computed the Hessian matrix . On the other hand, associated with the Gaussian and the highest eigenvalue ( ) is a good compromise when the feature widths are growing. Yet, the Weingarten and its eigenvalue give suitable and better continuous detected ridges. Steerable filters of order 2 ( ) and of order 4 ( ) obtain similar results, they are particularly reliable for images corrupted by noise, especially for bended features ( is a little more reliable), contrary to the SOAGK which is well adapted for straight features.

Eventually, this study would serve as ridge/valley optimal parameter configuration and adjustment guide for its interested applied researchers and application tools and domain such as satellite or aerial image analysis (road, river, etc.), medical image analysis (blood vessels, filaments, nerve system, etc.), lines detection, image segmentation, and object detection.

COMPLIANCE WITH ETHICAL STANDARDS

This manuscript is a completely original work of its authors; a preliminary work was published in IMTA VII -Workshop on Image Mining Theory and Applications.

2 .

 2 In the Cartesian space, considering the image surface , thus, defines all pixel coordi-plan of the surface in all points where

1 )

 1 hence, the use of a discrete Laplacian for 1D signals is the convolution of a signal with the vector ±[1 0 -2 0 1]. Regarding the 2-dimensional (2D) signals such as dig-

1 L 2 LFig. 1 . 1 Fig. 2 .

 12112 Fig. 1. Ridge and Valley highlighted on a discrete 1D signal using Laplacian.

Fig. 3 . 9 Fig. 4 .

 394 Fig. 3. Valley detection using a DOLP filter.

Fig. 5 .

 5 Fig. 5. Non-maxima suppression technique and different steps to obtain a thresholded image in 4 directions (usually used for ridge/valley detection with DOLD or Laplacian).

( a )

 a Image 256 u 256 (b) I xx image of (a) I yy image of (a) I xy image of (a) |k 1 | image of H |k 2 | image of H |k 2 | image of W (h) Valley detection, from (e)

  , smoothing the image with a Gaussian and then considering .

1

 1

Fig. 7 .

 7 Fig. 7. Comparison of valley detection on real images. The image in (a) is obtained using scanning electron microscopy of melt ceramic. Here, the valleys are detected with 3 different techniques: Hessian matrix (c) without and (d) with Gaussian , respectively and steerable filter of order 4 ( ) in (e).

G 4 SFFig. 8 .

 48 Fig. 8. Display of , , and computing local contrast as a function of . 15 4 100 80 60 40 20 0 10 8 6 4 2 0 2 4 6 8 10 3 2 1 0 10 8 6 4 2 0 2 4 6 8 10

Fig. 9 .Fig. 10 .Z

 910 Fig. 9. Example of extracted ridges with their tied perpendicular directions.

2

 2

Fig. 11 .

 11 Fig.11. Two-dimensional discrete filters for valley detection in y direction of width 7, parameters are available in Table1. The negative of filters detect ridges.

  the direction. Now, this 2D kernel can be oriented in different directions to capture valleys (or ridges with the opposite filter) in the image, as illustrated in Fig.12.

Fig. 12 .Fig. 13 .Fig. 14 .

 121314 Fig. 12. Representation and visualization of the second derivative of an image computed by convolution with the SOAGK with = 2.88 and .

  (1) , V = 1.81 S*n(1) , V = 1.81 r S*S p , V = 1.81 S*S N , V = 1.81 S*n(3) , V = 1.81

  Fig.[START_REF] Laligant | Noise estimation from digital step-model signal[END_REF]. Complexity schema, depending on the recursive filters order, the number of calculated images (Table3) and the filter rotation.

Fig. 16 .

 16 Fig. 16. Valley detection in green on real image of a dragonfly, with thin, blurred, and very close junctions. The original image is inverted for a better visualization.

1 DFig. 17 .

 117 Fig. 17. Evaluation of the different ridge/valley extraction techniques on synthetic images corrupted by Poisson noise.

R

  and H, D 1 G V and H, D 1 G V and H, D 2 G V and H, D 3 G V and H, D 4 G V and W SNR = 9 dB SNR = 3 dB SNR = 9 dB SNR = 9 dB SNR = 3 dB SNR = 3 dB SNR = 9 dB SNR =

  t

SNR = 9 5 SNR = 9

 959 dB SNR = 3 dB SNR = 9 dB SNR = 3 dB SNR = 3 dB SNR = 9 dB SNR, dB SNR = 3 dB SNR = 9 dB SNR, dB SNR = 3 dB SNR = dBroughly disconnected contours as of blob like structure, and Hessian with discrete filter (Fig.22c) has detected most of the discernible details but the majority of detected ridges are misconnected. The and filters in Figs.22d and 22ehave extracted similar roughly disconnected ridges, which cannot be enhanced even with post morphological processing.

displays less * 1 D * 2 D 4 DFig. 19 .RFig. 20 .R 5 Fig. 21 .

 1241920521 Fig. 19. Evaluation of the different ridge/valley extraction techniques on synthetic images corrupted by Gaussian noise and Gaussian blur.

( a )

 a Original image 256 u 256 (b)LoG, V = 1.81 (c) H with discrete filter and D 1 (d) H with filter Z and s z = 0.63 (e) H with filter R and s r = 0.69 (f) H with G, V = 1.81 and D 1 (g) H with G, V = 1.81 and D 2 (g) H with G, V = 1.81 and D 4 (i) W with G, V = 1.81 and D 1 (j) SF 2 , V = 1.81 (j) SF 4 , V = 1.81 (l) SOAGK, V u = 1.81, V v = 5V u information than with , it is similar to and . Finally, detected ridges extracted by oriented filters in Figs. 22j-22l demonstrate quiet interesting and directly usable information for its application analysis, in particular, the SOAGK for the straight ridges.

  This paper presents an extensive evaluation and comparison of ridge/valley detection with imagebased filtering techniques including the ridge/valley

Fig. 22 .

 22 Fig. 22. Ridge detection on Stanwick aerial image. Thresholded images correspond to 50% of the positive points after non-maxima suppression step, excepted for (c) where the 7500 highest points are reported because it contains too many positive pixels.

Table 1 .

 1 Optimum parameters of the high-pass filters z, r, and g σ as a function of the width of the ridge or valley Kernels based on the derivative of anisotropic Gaussian functions have been successfully applied in edge detection, showing certain advantages compared with the isotropic Gaussian derivatives

	Width 1	3	5	7	9	11 13 15 17
	s

T z 1.696 0.63 0.37 0.27 0.21 0.17 0.15 0.13 0.12 s r 2.05 0.69 0.42 0.31 0.24 0.20 0.17 0.15 0.13 σ 0.58 1.81 2.88 3.91 4.93 5.94 6.95 7.95 8.96 2.4.2. Anisotropic Gaussian filter.

Table 2 .

 2 Recursive orders of the filters and image computations as a function of the chosen technique

	Filter	Z and z	R and r	G σ and g σ

Table 3 .

 3 Image computations as a function of the chosen technique

		Computed basis images	Rotation	Other
	LoG	and		-	-
		,	,

18 .

 18 Evaluation of the different ridge/valley extraction techniques on synthetic images corrupted by speckle noise.R and H, D 1 G V and H, D 1 G V and H, D 2 G V and H, D 3 G V and H, D 4 G V and W

	0.85	900	0.9
	0.80	800
	0.75	700	0.8
	0.70	600
	0.65	500	0.7
	0.60	400
	0.55	300	0.6
	0.50	200
	SNR, dB	SNR, dB	SNR, dB
	0.98	0.08
	0.96	0.07	Z and H, D 1
	0.94	0.06
	0.92	0.05
		0.04
	0.90	0.03
	0.88	0.02	SF 2 SF 4
	0.86	0.01	SOAGK L/L

CONFLICT OF INTEREST

The content of the article does not give grounds for raising the issue of a conflict of interest.


Darshan Venkatrayappa
email: darshan.venkatrayappa@mines-ales.fr

Philippe Montesinos
email: philippe.montesinos@mines-ales.fr

Daniel Diep
email: daniel.diep@mines-ales.fr

Baptiste Magnier
email: baptiste.magnier@mines-ales.fr

RSD-HoG: A New Image Descriptor

Keywords: Image descriptor, Half gaussian kernel, Feature matching, Rotation signal descriptor, HoG

In this paper we propose a novel local image descriptor called RSD-HoG. For each pixel in a given support region around a key-point, we extract the rotation signal descriptor(RSD) by spinning a filter made of oriented anisotropic half-gaussian derivative convolution kernel. The obtained signal has extremums at different orientations of the filter. These characteristics are combined with a HoG technique, to obtain a novel descriptor RSD-HoG. The obtained descriptor has rich, discriminative set of local information related to the curvature of the image surface. With these rich set of features, our descriptor finds applications in various branches of computer vision. For evaluation, we have used the standard Oxford data set which has rotation, brightness, illumination, compression and viewpoint changes. Extensive experiments on these images demonstrates that our approach performs better than many state of the art descriptors such as SIFT, GLOH, DAISY and PCA-SIFT.

Introduction

In computer vision, problems related to object matching, tracking, panorama generation, image classification, structure and motion estimation are effectively addressed by the popular approach of image representation by a set of local image descriptors. The main purpose of local image descriptor is to capture the geometry of a support region around a key-point. In addition to this, the image descriptor should be invariant to certain image transformations such as rotation, brightness, blurring and scale changes. Scanning through the com-puter vision literature, one can come across the term image matching pipeline. The image matching pipeline has four stages. In the first stage, key-points or regions are selected using the popular detectors such as LoG[2], DoG [START_REF] Lowe | Distinctive Image Features from Scale-Invariant Keypoints[END_REF]orHar-ris Affine [START_REF] Mikolajczyk | Scale & Affine Invariant Interest Point Detectors[END_REF]. This is followed by the extraction of features or feature description from the support region around the key-point. Next, various post processing applied. The final block involves matching the descriptor using different distance measures such as euclidean distance [START_REF] Lowe | Distinctive Image Features from Scale-Invariant Keypoints[END_REF], hamming distance [START_REF] Calonder | BRIEF: Computing a Local Binary Descriptor Very Fast[END_REF], Earth Mover's Distance(EMD) [START_REF] Lazebnik | A Sparse Texture Representation Using Local Affine Regions[END_REF].

We are mainly interested in the second stage of the above described pipeline. The popular approach to extract features from a support region is to use the Histogram of Gradients(HoG) [START_REF] Dalal | Histograms of oriented gradients for human detection[END_REF]. SIFT [START_REF] Lowe | Distinctive Image Features from Scale-Invariant Keypoints[END_REF] effectively makes use of HoG to generate a descriptor. In PCA-SIFT [START_REF] Ke | PCA-SIFT: A more distinctive representation for local image descriptors[END_REF] the dimension of the SIFT descriptor is reduced by applying the PCA technique. GLOH [START_REF] Mikolajczyk | A performance evaluation of local descriptors[END_REF] and DAISY [START_REF] Tola | DAISY: An Efficient Dense Descriptor Applied to Wide Baseline Stereo[END_REF] use radial and daisy binning strategies respectively to improve on the existing SIFT descriptor. The standard version of SURF [START_REF] Bay | Speeded-Up Robust Features (SURF)[END_REF] achieves speed-up over SIFT by using Haar wavelets. MROGH [START_REF] Fan | Rotationally Invariant Descriptors Using Intensity Order Pooling[END_REF] and MRRID [START_REF] Fan | Rotationally Invariant Descriptors Using Intensity Order Pooling[END_REF] improve on the rotational invariance aspect by pooling local features based on their intensity orders in multiple support regions. In RIFF [START_REF] Takacs | Rotation-invariant fast features for large-scale recognition and real-time tracking[END_REF], radial and tangential components are extracted from the support region to form the descriptor. CSLBP [START_REF] Heikkilä | Description of interest regions with local binary patterns[END_REF] and MRRID [START_REF] Fan | Rotationally Invariant Descriptors Using Intensity Order Pooling[END_REF] use local binary patterns to encode the image structure in the support region. A detailed evaluation of these descriptors can be found in [START_REF] Mikolajczyk | A Performance Evaluation of Local Descriptors[END_REF].

There is a line of research oriented towards object recognition and image matching on camera-enabled mobile devices (e.g. phones and tablets). Limited computational power and storage space on these devices has enabled the emergence of binary descriptors. The main idea behind binary descriptor is that each bit is independent and Hamming distance can be used as similarity measure. Some of the most popular binary descriptors are BRIEF [START_REF] Calonder | BRIEF: Computing a Local Binary Descriptor Very Fast[END_REF], ORB [START_REF] Rublee | ORB: An efficient alternative to SIFT or SURF[END_REF], BRISK [START_REF] Leutenegger | BRISK: binary robust invariant scalable keypoints[END_REF] and FREAK [START_REF] Alahi | FREAK: fast retina keypoint[END_REF]. A detailed evaluation of these descriptors can be found in [START_REF] Figat | Performance evaluation of binary descriptors of local features[END_REF].

Filter responses has been used in abundance for image description. Schmid and Mohr [START_REF] Schmid | Local grayvalue invariants for image retrieval[END_REF] use differential invariant responses to compute new local image descriptors. Differential invariant responses are obtained from a combination of Gaussian derivatives of different orders, which are invariant to 2-dimensional rigid transformations. Larsen et al. [START_REF] Larsen | Jet-based local image descriptors[END_REF] follow a new approach for the construction of an image descriptor based on local k-jet, which uses filter bank responses for feature description. Palomares et al. [START_REF] Palomares | A new affine invariant method for image matching[END_REF] have come up with a local image descriptor issued from a filtering stage made of oriented anisotropic half-gaussian smoothing convolution kernels. In their work the authors use euclidean invariance is achieved by FFT correlation between the two signatures/descriptors. Moderate deformation invariance is achieved using Dynamic Time Warping (DTW).

The contributions of this paper include a new approach for descriptor generation. We embed the response of an anisotropic half Gaussian kernel in a HoG framework. For each pixel in the support region a rotation signal descriptor(RSD) is extracted by spinning a filter made of oriented anisotropic half-Gaussian derivative convolution kernel. We select two orientations at which a global maxima and a global minima occurs. Thus, for each pixel in the support region we have two orientations and a magnitude. These orientations are binned by weighing them with the magnitude.

Image-Processing

Isotropic Gaussian filters Fig. 1(a), Fig. 1(b) are widely used for smoothing images. But, it is well known that this type of filters blur the image features. Gaussian anisotropic filtering overcome this by preserving the features along a particular direction. Compared to isotropic filters, anisotropic filters are elongated along one of the directions. Fig. 1(c) and Fig. 1(d) shows smoothing and derivative anisotropic filters respectively. We use a filter made of anisotropic half Gaussian derivative kernel as in [START_REF] Magnier | Evolution image with PDEs toward a new anisotropic smoothing based on half kernels[END_REF] [START_REF] Montesinos | A new perceptual edge detector in color images[END_REF]. An anisotropic half-Gaussian derivative filter 1(f) is designed with a Gaussian half-smoothing filter in Y direction and a derivative Gaussian filter in X direction. The filtering in the Y direction acts as a causal integration, as the filtering in X direction is a full derivative. Then, rotating the filter around pixels generates a description around the pixels. Applying such filters on an edge point will produce a minima and a maxima in two directions. On a straight line edge these two direction are exactly opposite (180 • ), but generally the difference of angles reflects the curvature at edge points.

ϕ(x, y, θ) = I θ * C • H(-y) • x • e - ⎛ ⎝ x 2 2λ 2 + y 2 2μ 2 ⎞ ⎠ (1)
Where I θ corresponds to a rotated image of orientation θ, C is a normalization coefficient, (x, y) are pixel coordinates, and (μ, λ) the standard deviation of the anisotropic Gaussian kernel. Only the causal part of this filter along the Y axis is used. This is obtained by cutting the kernel in the middle, in an operation that corresponds to the Heaviside function H. As in [START_REF] Montesinos | A new perceptual edge detector in color images[END_REF] we have chosen to rotate the image instead of the filter, there by reducing the algorithmic complexity by making use of a recursive Gaussian filter [START_REF] Deriche | Recursively implementing the gaussian and its derivatives[END_REF].

Rotational Signal Descriptor(RSD)

RSD is obtained by rotating the above described filter around a key-point. Fig. 2(c) shows a sample RSD obtained by applying the Gaussian derivative half filter at the the pixel location (x p , y p ) in steps of 2 • . We compute the gradient ∇I and the two angles for the key-point/pixel at location (x p , y p ) by considering the global extrema of the function Q(x p , y p , θ). The two angles θ 1 and θ 2 define a curve crossing the pixel(an incoming and outgoing direction). Two of these global extrema are combined to maximize ∇I , i.e :

⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ ∇I = max θ∈[0,360 • [ Q(x p , y p , θ) -min θ∈[0,360 • [ Q(x p , y p , θ) θ 1 = arg max θ∈[0,360 • [ (Q(x p , y p , θ)) θ 2 = arg min θ∈[0,360 • [ (Q(x p , y p , θ)) (2)
We can spin an anisotropic half Gaussian derivative kernel at a key-point and the resulting response(RSD) can be considered as a descriptor. But, RSD alone gives a weak description, as it fails to capture the geometry around the key-point. This is the main motivation for combining the HoG technique with RSD.

Pre-Processing

Prior to feature extraction, the image is smoothed with a Gaussian filter. Then, key-points(regions) need to be localized. We have used Harris affine feature detector. It can also be noted that any of the existing feature detector can be used for this step. As in Fig. 3, we first extract image patches surrounding the key-points. Then, depending on the key-point attributes such as scale and orientations, we rotate the patch by its orientation. This is followed by scale normalizing the oriented patch to a fixed size of 41x41. More details about the process for patch extraction and normalization can be found in [START_REF] Mikolajczyk | Scale & Affine Invariant Interest Point Detectors[END_REF].

RSD-HoG Construction

The framework of the RSD-HOG extraction is illustrated in Fig. 3. On each pixel of the normalized patch, we apply the rotating semi Gaussian filter to obtain the RSD. From this RSD, we extract two angles and a magnitude for each pixel as explained above. Then we bin the two angles separately as in Eq.3 and Eq.4.

The image patch is divided in to 16 blocks. Since the image patch is of size 41x41, most of the blocks are of size 10x10 ( blocks on the extreme right and extreme bottom are of size 11x11 ). Each of these block contributes 8 bins to the final descriptor. We fuse the two intermediate descriptors to form the final descriptor as in Eq.5. The intermediate descriptors in Eq.3 and Eq.4 alone can be used as descriptors. But, fusing these two descriptors as in Eq.5 results in a more robust description. The performance of two intermediate descriptor and the final descriptor for the boat dataset(Rotation changes) can be seen in the first row of Fig. 4.

RSD -HoG -T heta1 = {θ 1 bin1 , θ 1 bin2 , θ 1 bin3 , θ 1 bin4 .....θ 1 bin128 } (3) RSD -HoG -T heta2 = {θ 2 bin1 , θ 2 bin2 , θ 2 bin3 , θ 2 bin4 .....θ 2 bin128 } (4) RSD -HoG = {θ 1 bin1 , θ 1 bin2 , ...θ 1 bin128 , θ 2 bin1 , θ 2 bin2 , ...θ 2 bin128 } (5)
3 Experiments, Discussions and Results

Dataset and Evaluation

We evaluate and compare the performance of our descriptor against the state of the art descriptors on the standard dataset using the standard protocol provided by Oxford group. The binaries and dataset are obtained from website linked to [START_REF] Mikolajczyk | A Performance Evaluation of Local Descriptors[END_REF] (http://www.robots.ox.ac.uk/ ∼ vgg/research/affine/) . The dataset used in our experiments has different geometric and photometric transformations such as change of scale, image rotation, viewpoint change, image blur, JPEG compression and illumination change. For each type of image transformation there is a set of six images with established ground truth homographies. Recall versus precision curves as proposed by [START_REF] Mikolajczyk | A Performance Evaluation of Local Descriptors[END_REF] are used to evaluate our descriptor. This is based on the number of correspondences, correct matches and false matches obtained for an image pair. Number of correspondences is defined as the number of pairs of features whose euclidean distance is below a certain threshold. We vary the threshold in steps to obtain the recall vs precision curves. A correct match is recorded when the pair of features with the smallest euclidean distance coincides with the pair of points/regions matched by the ground-truth. As in Eq.6, recall is defined as the total number of correctly matched regions over the of corresponding regions between two images of the same scene. From Eq.7, 1-precision is represented by the number of false matches relative to the total number of matches.

Parameter Selection

Our descriptor has 4 different parameters that are tabulated in the table 1. The rotation step of the filter is fixed to 5

Descriptor Performance

The performance of RSD-HoG is compared against SIFT-OXFORD, SIFT-PATCH, GLOH, DAISY and PCA-SIFT. For SIFT-OXFORD, PCA-SIFT and GLOH the descriptors are extracted from the binaries provided by Oxford group (http://www.robots.ox.ac.uk/ ∼ vgg/research/affine/) [START_REF] Mikolajczyk | A Performance Evaluation of Local Descriptors[END_REF]. DAISY descriptor is extracted from the code provided by SIFT-PATCH is the SIFT descriptor applied on a patch using the VLFEAT library [START_REF] Vedaldi | VLFeat: An Open and Portable Library of Computer Vision Algorithms[END_REF]. Due to lack of space we restrict ourself to image pairs(1-3),(1-4) and (1-5).

1. Rotation changes. First, we compare the performance of RSD -HoG (Eq.5) with the intermediate descriptors RSD -HOG -T heta1 (Eq.3)and RSD -HOG -T heta2 (Eq.4). From the graphs in the first row of Fig. 4 it is clear that RSD -HoG performs better than the two intermediate descriptors. This is same for all the images with different transformations. Due to lack of space we restrict ourselves to graphs for rotational changes. When RSD -HOG is compared with other descriptors it can be seen that our descriptor performs better than other descriptors. This can be seen in the 2nd row of Fig. 4. 2. Viewpoint and Blur changes.

Graphs in the third row of Fig. 4 represent Recall vs Precision curves for the Viewpoint changes. It can be seen that for the first 2 cases (1-2),(1-3)and (1-4) RSD -HoG performs similar to other descriptors. For the final case the performance of RSD -HoG deteriorates. This is a challenging sequence and all the other descriptors perform badly. Graphs in the fourth row of Fig. 4 represent Recall vs Precision curves for the Blur changes. From the graphs it is clear that RSD -HoG outperforms other descriptors.

Brightness and Compression changes.

Based on graphs in the fifth row of Fig. 4, we conclude that the performance of RSD -HoG is superior to the performance of other descriptors when it comes to handling Brightness variations. Similarly, graphs in the last row of Fig. 4 illustrate the dominance of our descriptor for compression variations.

Conclusion

This paper proposes a new image descriptor called RSD-HoG. It also proposes a new approach to construct the descriptor by interlacing the bins of the two intermediate descriptors. On the standard dataset provided by the Oxford group, RSD-HoG outperforms other state of the art descriptors. Currently, high complexity and the dimension of the descriptor are a major drawback. In the future we would like to reduce both the complexity and dimension of our descriptor.

Here, we have used a fixed set of parameters for the anisotropic half Gaussian kernel. In the future, we would like to experiment with other variations of the anisotropic half Gaussian kernel. We would also like to focus on the real time implementation of our descriptor using parallel programming techniques. Another direction of our future work would be to test our approach on tasks related to object detection, writer classification and shape retrieval.

Fig. 1 .

 1 Fig. 1. (a) Isotropic smoothing filter. (b) Isotropic derivative filter. (c) Anisotropic smoothing filter. (d) Anisotropic derivative filter. (e) Half anisotropic smoothing filter. (f) Half anisotropic derivative filter.

Fig. 2 .

 2 Fig. 2. (a) A thin rotating Gaussian derivative half-filter. (b) Half Gaussian kernel applied to a keypoint xp, yp. (c) Extrema of a function Q(xp, yp, θ). μ = 10, λ = 1 and Δθ = 2 o . Note that the initial orientation of the filter is vertical, upwardly directed and steerable clockwise.

Fig. 3 .

 3 Fig.3. Affine region is normalized to a square of size 41x41. On each pixel of the patch, RSD is generated and the two angles θ1 and θ2 are extracted. These angles are binned separately.

  matches Number of correct matches + Number of false matches[START_REF] Lazebnik | A Sparse Texture Representation Using Local Affine Regions[END_REF]

Fig. 4 .

 4 Fig. 4. Recall vs Precision curves for Oxford dataset using nearest neighbour matching method

Table 1 .

 1 • . Increasing the rotation step results in loss of information. We have fixed the number of bins to 8 per block, resulting in a 8 * 16 = 128 bins for 16 blocks. Increasing the number of bins results in the same performance but, increases the dimensionality of the descriptor. All the parameters are chosen empirically. Parameters

	filter Height (μ) filter Width (λ) Rotation step (Δθ) No of BINS
	6	1	5 •	8