Adrian Ion
email: ion@ins.uni-bonn.de

Nicole M Artner
email: artner@prip.tuwien.ac.at

Gabriel Peyré
email: peyre@ceremade.dauphine.fr

Walter G Kropatsch

Laurent D Cohen
email: cohen@ceremade.dauphine.fr

Matching 2D and 3D Articulated Shapes Using the Eccentricity Transform

Keywords: eccentricity transform, shape matching, articulation, geodesic distance

This paper presents a novel method for 2D and 3D shape matching that is insensitive to articulation. It uses the eccentricity transform, which is based on the computation of geodesic distances. Geodesic distances computed over a 2D or 3D shape are articulation insensitive. The eccentricity transform considers the length of the longest geodesics. Histograms of the eccentricity transform characterize the compactness of a shape, in a way insensitive to rotation, scaling, and articulation. To characterize the structure of a shape, a histogram of the connected components of the level sets of the transform is used. These two histograms make up a highly compact descriptor and the resulting method for shape matching is straightforward. Experimental results on established 2D and 3D benchmarks show results similar to more complex state of the art methods, especially when considering articulation. The connection between the geometrical modification of a shape and the corresponding impact on its histogram representation is explained. The influence of the number of bins in the two histograms and the respective importance of each histogram is studied in detail.

Introduction

The recent increase in available 3D models and acquisition systems has created the need for efficient retrieval of stored models, making 3D shape matching gain attention also outside the computer vision community. Together with its 2D counterpart, 3D shape matching is useful for identification and retrieval in classical vision tasks, but can also be found in Computer Aided Design/Computer Aided Manufacturing (CAD/CAM), virtual reality (VR), medicine, molecular biology, security, and entertainment [START_REF] Bustos | Featurebased similarity search in 3D object databases[END_REF].

Shape matching requires to set up a signature that characterizes the properties of interest for the recognition [START_REF] Veltkamp | Properties and performance of shape similarity measures[END_REF]. Depending on the task, the invariance of this signature to local deformations such as articulation is important for the identification of 2D and 3D shapes. Matching can then be carried out over this (usually lower dimensional) space of signatures.

Related work

Most shape descriptors are computed over a transformed domain that amplifies the important features of the shape while throwing away ambiguities such as translation, rotation or local deformations.

For 2D shapes, the Fourier transform of the boundary curve [START_REF] Zahn | Fourier descriptors for plane closed curves[END_REF] is an example of such a transformed-domain descriptor adapted to smooth shapes. Shape transformations computed with geodesic distances [START_REF] Bronstein | Analysis of two-dimensional non-rigid shapes[END_REF][START_REF] Hamza | Geodesic matching of triangulated surfaces[END_REF] lead to signatures invariant to isometric deformations such as bending or articulation. To capture salient features of 2D shapes, local quantities such as curvature [START_REF] Mokhtarian | A theory of multiscale, curvaturebased shape representation for planar curves[END_REF] or shape contexts [START_REF] Belongie | Shape matching and object recognition using shape contexts[END_REF] can be computed. They can be extended to bending invariant signatures using geodesic distances [START_REF] Ling | Shape classification using the inner-distance[END_REF][START_REF] Nasreddine | Shape geodesics for boundarybased object recognition and retrieval[END_REF]. Another possibility is to represent a shape as a collection of modestly overlapping disk components, which contain both local geometric and structural information [START_REF] Xu | Shape matching using morphological structural shape components[END_REF]. More global features include the Laplace spectra [START_REF] Reuter | Laplace-spectra as fingerprints for shape matching[END_REF] and the skeleton [START_REF] Siddiqi | Shock graphs and shape matching[END_REF]. Contour flexibility [START_REF] Xu | 2d shape matching by contour flexibility, Pattern Analysis and Machine Intelligence[END_REF] is a novel shape descriptor of planar contours, which obtains both local and global features from a contour. It represents the deformable potential at each point along a contour. The rolling penetrate descriptor [START_REF] Chen | Rolling penetrate descriptor for shape-based image retrieval and object recognition[END_REF] combines the advantages of contour-based and region-based methods, and provides a unified scheme to handle various shapes, geometrical transforms, noise, distortion and occlusion. Some transformations involve the computation of a function defined on the shape, for instance the solution to a linear partial differential equation [START_REF] Gorelick | Shape representation and classification using the poisson equation[END_REF] or geometric quantities [START_REF] Osada | Shape distributions[END_REF].

Among approaches matching 3D shapes, existing methods can be divided into [START_REF] Bustos | Featurebased similarity search in 3D object databases[END_REF]: statistical descriptors, like for example geometric 3D moments employed by [START_REF] Elad | Content based retrieval of vrml objects: an iterative and interactive approach[END_REF][START_REF] Paquet | Description of shape information for 2-D and 3-D objects[END_REF], 3D moment invariants [START_REF] Sadjadi | Three-dimensional moment invariants, Pattern Analysis and Machine Intelligence[END_REF], and the shape distribution [START_REF] Osada | Shape distributions[END_REF][START_REF] Ip | Automated learning of model classifications[END_REF]. A novel shape representation based on [START_REF] Sadjadi | Three-dimensional moment invariants, Pattern Analysis and Machine Intelligence[END_REF] is 3D gray level moment invariants [START_REF] Guo | 3D gray level moment invariants: A novel shape representation[END_REF], which are independent of translation, scaling and rotation. Extension-based descriptors are calculated from features sampled along certain directions from a position within the shape [START_REF] Vranic | Description of 3D-shape using a complex function on the sphere[END_REF][START_REF] Vranic | 3D model retrieval with spherical harmonics and moments[END_REF]. Volumebased descriptors use the volumetric representation of a 3D shape to extract features (examples are Shape histograms [START_REF] Ankerst | 3D shape histograms for similarity search and classification in spatial databases[END_REF], Model Voxelization [START_REF] Vranic | 3D shape descriptor based on 3D fourier transform[END_REF], and point set methods [START_REF] Tangelder | Polyhedral model retrieval using weighted point sets[END_REF]). The 3D shape impact descriptor [START_REF] Mademlis | 3D object retrieval using the 3D shape impact descriptor[END_REF] considers a 3D object as a distributed 3D mass and it is indirectly computed from the resulting fields. The field is described using both Newton's and general relativity laws. In [START_REF] Ruggeri | Spectral-driven isometry-invariant matching of 3D shapes[END_REF], 3D shapes are sampled using a technique based on critical points of the eigenfunctions of the Laplace-Beltrami operator. A point-based statistical descriptor is used that incorporates an approximation of the geodesic shape distribution and other geometric information describing the surface at that point. Matching is carried out using Bipartite graph matching. Descriptors using the surface geometry compute curvature measures and/or the distribution of surface normal vectors [START_REF] Paquet | Nefertiti: A tool for 3-D shape databases management[END_REF][START_REF] Zaharia | Three-dimensional shape-based retrieval within the mpeg-7 framework[END_REF]. Image-based descriptors reduce the problem of 3D shape matching to an image similarity problem by comparing 2D projections of the 3D shapes [START_REF] Ansary | A bayesian framework for 3D models retrieval based on characteristic views[END_REF][START_REF] Cyr | A similarity-based aspect-graph approach to 3D object recognition[END_REF][START_REF] Chen | On visual similarity based 3D model retrieval[END_REF]. Reeb graphs have been used to match the topology of two shapes [START_REF] Hilaga | Topology matching for fully automatic similarity estimation of 3D shapes[END_REF][START_REF] Shinagawa | Surface coding based on morse theory[END_REF]. Skeletons are intuitive shape descriptions and can be obtained from a 3D shape by applying a thinning algorithm on the voxelization of a solid object like in [START_REF] Sundar | Skeleton based shape matching and retrieval[END_REF]. Descriptors using spin images work with a set of 2D histograms of the shape geometry and a search for point-to-point correspondences is done to match 3D objects [START_REF] Johnson | Using spin images for efficient object recognition in cluttered 3D scenes[END_REF]. In [START_REF] Reuter | Hierarchical shape segmentation and registration via topological features of laplace-beltrami eigenfunctions[END_REF], 3D shapes are automatically decomposed into parts using topological features of the Laplace-Beltrami Eigenfunctions. The parts of near isometric shapes are registered to each other.

A generic class of approaches for 2D and 3D shape retrieval is to define a metric between pairs of points on the shape, and compare either directly these metric spaces or compare features extracted from these spaces. An important goal of this metric space design is to make the shape retrieval more or less invariant to bending and articulations. The Gromov-Hausdorff framework directly compares metric spaces [START_REF] Mémoli | A theoretical and computational framework for isometry invariant recognition of point cloud data[END_REF], and can be used for shape retrieval in conjunction with geodesic metric spaces [START_REF] Bronstein | Analysis of two-dimensional non-rigid shapes[END_REF] or diffusion spaces that are more robust to topological noise [START_REF] Mémoli | Spectral gromov-wasserstein distances for shape matching[END_REF][START_REF] Bronstein | A gromov-hausdorff framework with diffusion geometry for topologically-robust non-rigid shape matching[END_REF]. Associated to the diffusion metric, it is possible to define a metric using a dimensionality reduction within a few eigenvectors of the Laplacian [START_REF] Rustamov | Laplace-beltrami eigenfunctions for deformation invariant shape representation[END_REF]. To speed-up retrieval applications, one can consider low dimensional features extracted from these metric spaces. This can be for instance the set of geodesic distances between critical points of the Laplace eigenfunction [START_REF] Ruggeri | Spectral-driven isometry-invariant matching of 3D shapes[END_REF], local distributions of geodesic curves [START_REF] Ling | Shape classification using the inner-distance[END_REF], statistical moments of a diffusion distance [START_REF] Gorelick | Shape representation and classification using the poisson equation[END_REF]. A popular class of approaches considers histograms, such as the distribution of Euclidean distances [START_REF] Osada | Shape distributions[END_REF], of the mean geodesic distance [START_REF] Hamza | Geodesic matching of triangulated surfaces[END_REF] or the maximum geodesic distance [START_REF] Ion | Shape matching using the geodesic eccentricity transform -a study[END_REF][START_REF] Ion | 3D shape matching by geodesic eccentricity[END_REF], or bags of features [START_REF] Bronstein | Shape google: geometric words and expressions for invariant shape retrieval[END_REF]. This article elaborates on this idea of building compact geodesic descriptors using two different kinds of histograms to represent faithfully the distribution of geodesic distances.

Considering the underlying transform on which our descriptor is built, probably the most similar works in shape matching are [START_REF] Hamza | Geodesic matching of triangulated surfaces[END_REF][START_REF] Hilaga | Topology matching for fully automatic similarity estimation of 3D shapes[END_REF], where instead of the length of the geodesic to the point furthest away (this work), the mean or the sum over all points are considered. Distance based transforms are also used in [START_REF] Gorelick | Shape representation and classification using the poisson equation[END_REF][START_REF] Siddiqi | Shock graphs and shape matching[END_REF][START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF], where the length of random or shortest paths to an (existing) boundary are computed. Describing a shape as a normalized histogram of the values of a function at all points of the shape, as in the case of the normalized ECC histogram in this paper, is conceptually similar to the shape distributions in [START_REF] Osada | Shape distributions[END_REF], where different functions and histogram comparison methods are considered.

In [START_REF] Ling | Shape classification using the inner-distance[END_REF] a model of articulated objects is presented. It is defined as a union of (rigid) parts O i and joints (named 'junctions' by the authors). An articulation is defined as a transformation that is rigid when limited to any part O i , but can be non-rigid on the junctions. An articulated instance of an object is an articulated object itself (actually the same object) that can be articulated back to the original one. The term articulated shape refers to the shape of an articulated object in a certain pose. In the context of shape matching the concept of articulated shape means that shapes that belong to articulations of the same object, belong to the same class. Assuming that the size of the junctions is very small compared to the size of the parts O i , it is shown that the variation of the geodesic distance1 during articulation is small and that geodesic distances are articulation insensitive.

Similar to the method in [START_REF] Ling | Shape classification using the inner-distance[END_REF], our method does not explicitly involve any part models. In [START_REF] Ling | Shape classification using the inner-distance[END_REF] the part based articulation model was used to support the analysis of the properties of the geodesic distance under articulation. We have briefly recalled it in Section 1.1 to aid the definition of an articulated shape. Finding the correspondences between all the parts of two shapes is an NP -complete problem in graph theory (known also as the 'matching' of two graphs) and requires the correct decomposition of the unknown object into parts. A one-to-one correspondence (bijection) for all parts is not always possible as some parts might be missing (e.g. due to segmentation errors).

Contribution

The eccentricity transform of a shape associates to each of its points the distance to the point furthest away. It is based on the computation of geodesic distances and thus robust with respect to articulation. It is robust against minor segmentation errors, and Salt and Pepper like noise [START_REF] Kropatsch | The eccentricity transform (of a digital shape)[END_REF], and stable in the presence of holes (i.e. it is defined in the same way for shapes with and without holes, and does not require pre-selection of for e.g. a single closed boundary to be processed).

In a common framework, we propose histograms built on the eccentricity transform as descriptors for 2D and 3D shape matching. The descriptor consists of two histograms: the ECC histogram h which is the normalized histogram of the eccentricity transform of the shape, and the ECC structure histogram s which is the histogram of the number of connected components of the level sets of the eccentricity transform. The descriptor is invariant to changes in orientation, scale, and articulation. It requires only a simple representation and can be efficiently matched. We present an in-depth study of the properties of the approach (relation between shape and descriptor, parameters), supported by experimental results, and an analysis of the results and possibilities for improvement. The descriptor is computed on four different domains using geodesic distances: inside 2D shapes, inside the volume, inside the border voxels, or on the surface meshes of 3D shapes, and compared to state of the art methods. These numerical results support that the proposed approach performs similarly and in some case better than the state of the art on databases of articulated 2D and 3D shapes. Initial results using only the ECC histogram h have been presented for 2D in [START_REF] Ion | Shape matching using the geodesic eccentricity transform -a study[END_REF], and for 3D (volumetric representation only) in [START_REF] Ion | 3D shape matching by geodesic eccentricity[END_REF].

To the best of our knowledge, this is the first approach applying the eccentricity transform to the problem of shape matching.

Overview of the Paper

The paper is organized as follows: Section 2 recalls the eccentricity transform and discusses used variants and computation. Section 3 explains the proposed matching method and discusses pros and cons of the descriptor (Section 3.4). In Section 4 experiments, and the effect of the parameters are presented, followed by future work (Section 5) and conclusion (Section 6). The Appendix recalls the algorithm used to compute the eccentricity transform.

Eccentricity Transform

The following definitions and properties follow [START_REF] Kropatsch | The eccentricity transform (of a digital shape)[END_REF][START_REF] Ion | Euclidean eccentricity transform by discrete arc paving[END_REF], and are extended to n-dimensional domains.

Let the shape S be a closed set in R n . A path π in S is the continuous mapping from the interval [0, 1] to S. Let Π(p 1 , p 2 ) be the set of all paths between two points p 1 , p 2 ∈ S within the set S.

The geodesic distance d(p 1 , p 2 ) between two points p 1 , p 2 ∈ S is defined as the length λ(π) of the shortest path π ∈ Π(p 1 , p 2 )

d(p 1 , p 2 ) = min{λ(π) | π ∈ Π(p 1 , p 2 )}, (1)
where the length λ(π) is The eccentricity transform of S is defined as, ∀p ∈ S ECC(S, p) = max{d(p, q) | q ∈ S}.

λ(π(t)) = 1 0 || π(t)||dt, π(t)
(

) 2
To each point p it assigns the length of the geodesic path(s) to the points farthest away from it. The definition above accommodates n-dimensional objects embedded in R n as well as n-dimensional objects embedded in higher dimensional spaces (e.g. the 2D manifold given by the surface of a closed 3D object). The distance between any two points whose connecting segment is contained in S, is computed using the L 2 -norm, i.e. distances are not computed on a graph, but are a discretization of the continuous geodesic distance. For a definition of the ECC of a graph see [START_REF] Kropatsch | The eccentricity transform (of a digital shape)[END_REF].

The ECC is quasi-invariant to articulated motion and robust against salt and pepper noise, which creates small (typically 1 pixel) holes in the shape [START_REF] Kropatsch | The eccentricity transform (of a digital shape)[END_REF]. An analysis of the variation of geodesic distance under articulation can be found in [START_REF] Ling | Shape classification using the inner-distance[END_REF].

An eccentric point is a point q that reaches a maximum in Equation 2, and for most shapes, all eccentric points lie on the border of S [START_REF] Kropatsch | The eccentricity transform (of a digital shape)[END_REF]. An eccentric path of a point p is a geodesic to one of its eccentric points. The (geodesic) center is the set of points that have the smallest eccentricity (global minimum). The diameter of a shape S is the maximum ECC, which is the length of the longest geodesic path in S.

The classes of 2n-connected discrete shapes S defined by points on a square grid Z n , n ∈ {2, 3}, as well as connected triangular meshes representing the surface of the 6-connected 3D shapes are considered in the paper. Table 1 shows the types of manifolds used in this article, for which ECC is computed. For ECCobj2D, ECCobj, and ECCborder, paths need to be contained in the area of R n defined by the union of the support squares/cubes for the pixels/voxels of S. For ECCmesh, paths need to be contained in the 2D manifold defined by the union of the triangles of the mesh (including the interior of the triangles). The used (approximated) metric is in all cases the Euclidean based geodesic distance d ε , i.e. the distance between the endpoints of any line segment included in the shape is computed using the Euclidean distance. When the resolution of the shapes increases ECCborder [START_REF] Manay | Integral invariants for shape matching[END_REF]), random missing pixels [START_REF] Ling | Shape classification using the inner-distance[END_REF], and adding a large hole [START_REF] Nasreddine | Shape geodesics for boundarybased object recognition and retrieval[END_REF]. Bottom: variation of the values of DT, GGF, and ECC due to the 8 modifications. The mean and max values give the mean/max difference in the values of pixels in the original and altered shapes, normalized by the difference between the highest and lowest value on the original shape. The "% same" values give for cases [START_REF] Veltkamp | Properties and performance of shape similarity measures[END_REF][START_REF] Zahn | Fourier descriptors for plane closed curves[END_REF][START_REF] Bronstein | Analysis of two-dimensional non-rigid shapes[END_REF][START_REF] Hamza | Geodesic matching of triangulated surfaces[END_REF][START_REF] Mokhtarian | A theory of multiscale, curvaturebased shape representation for planar curves[END_REF] the ratio of pixels in the altered shape with values unchanged as compared to the whole hand and to the values in the previous shape [START_REF] Bustos | Featurebased similarity search in 3D object databases[END_REF][START_REF] Veltkamp | Properties and performance of shape similarity measures[END_REF][START_REF] Zahn | Fourier descriptors for plane closed curves[END_REF][START_REF] Bronstein | Analysis of two-dimensional non-rigid shapes[END_REF][START_REF] Hamza | Geodesic matching of triangulated surfaces[END_REF], respectively.

and ECCmesh converge to the same value.

Similarity to other transforms

The eccentricity transform is part of a greater class of distance-based transforms along with:

• the distance transform (DT) [START_REF] Soille | Morphological Image Analysis[END_REF] defined as:

DT (S, p) = min{d(p, q) | q ∈ Y ⊆ S}
where the set Y is called marker set and is usually taken as Y = ∂S;

• the global geodesic function (GGF) [START_REF] Hilaga | Topology matching for fully automatic similarity estimation of 3D shapes[END_REF] defined as:

GGF (S, p) = q∈S d(p, q).
Shortest geodesics (DT) have the advantage that they locally characterize the shape and are invariant to deformations in all other parts of the shape. On the other side, the shorter a path, the higher the perturbation that is created by a hole (obstacle).

ECC uses longest geodesics and has thus the highest stability with respect to small holes. With respect to deformations, ECC is invariant to nonisometric deformations that do not affect the parts where eccentric points lie, i.e. as long as a certain part does not contain any eccentric points it can even be removed, without changing the values for the remaining points of the shape [START_REF] Ion | The eccentricity transform of n-dimensional shapes with and without boundary[END_REF] (e.g. for a shape like a "T" with the vertical line connected to the center of an at least twice longer horizontal line, removing the vertical line will have no effect on the eccentricity values on the horizontal line).

Figure 1 shows the changes in the mean and maximum values of the three transforms when the shape of a hand from the Kimia 99 database [START_REF] Sebastian | Recognition of shapes by editing their shock graphs[END_REF] is changed. The DT is more local and thus insensitive to the missing fingers, but comparable in the case of boundary noise, and the highest variation with respect to missing points inside. Due to the averaging that the GGF does, it has a maximum variation which is in all cases smaller that of the ECC, but almost all pixels change even for a minor change in the shapes. The variation of the ECC is larger then the one of the GGF, but with a considerable amount of pixels keeping the same value before and after a finger is removed.

Both the DT and the GGF have previously been used as a basis for shape matching methods. The DT is closely related to the skeleton/medial axis/medial surface [START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF][START_REF] Zhang | Retrieving articulated 3-D models using medial surfaces and their graph spectra[END_REF] of a 2D/3D shape, on which shock graphs [START_REF] Siddiqi | Shock graphs and shape matching[END_REF] are built. Matching of shapes can be carried out by matching their shock graphs or medial surfaces. A continuous definition of the GGF has been used in [START_REF] Hilaga | Topology matching for fully automatic similarity estimation of 3D shapes[END_REF] to match the topology of 3D shapes based on the similarity of their Reeb graphs. In [START_REF] Hamza | Geodesic matching of triangulated surfaces[END_REF] a discrete version of the squared GGF is used to build so called "geodesic shape distributions". Kernel density estimation is used to get a smooth function from the computed samples, and the Jensen-Shannon divergence is used to match two distributions. Our normalized eccentricity histograms are in spirit similar to the geodesic shape distributions. Differences lie in the used function (ECC vs squared GGF) and in the normalization and matching of histograms (L 2 -norm). In addition we use a second histogram to characterize the structure of a shape.

Computation

Furthest point computation has been approached in [START_REF] Suri | The all-geodesic-furthest neighbor problem for simple polygons[END_REF], where an algorithm is presented which finds for each vertex of a simple polygon the vertex that is furthest away (eccentric). Later, [START_REF] Maisonneuve | An efficient algorithm to compute the hexagonal and dodecagonal propagation function[END_REF][START_REF] Schmitt | Propagation function: Towards constant time algorithms[END_REF] proposed an efficient algorithm for simply connected shapes on the hexagonal and dodecagonal grid. The concept of the eccentricity of a vertex can be found in classical graph theory books [START_REF] Harary | Graph Theory[END_REF][START_REF] Diestel | Graph Theory[END_REF], and the concept of the eccentricity transform2 in recent discrete geometry [START_REF] Klette | Digital Geometry[END_REF] and mathematical morphology books [START_REF] Soille | Morphological Image Analysis[END_REF]. Computation is not discussed and no references to holes in a shape are made.

The straightforward approach to compute ECC(S) is: for each point of S, compute the distance to all other points and take the maximum. A faster computation and efficient approximation algorithms are presented in [START_REF] Ion | Euclidean eccentricity transform by discrete arc paving[END_REF]. For this paper the fastest one, algorithm ECC06, is used.

ECC06 relies on the computation of the shape bounded single source distance transform3 D S (p) (Figure 2(b)), which is computed for estimated eccentric point candidates in an iterative manner (see the Appendix for more details). D S (p) associates to each point q ∈ S the geodesic distance to p. D S can be computed using Fast Marching [START_REF] Sethian | Level Sets Methods and Fast Marching Methods, 2nd Edition[END_REF] (FM), without the need to explicitly build a neighborhood graph. The runtime complexity of FM is O(N log(N)) steps, for N = |S| grid points, where | | denotes set cardinality.

The complexity for computing ECC(S) using ECC06 and Fast Marching is O(KN log(N)), where 2 K |∂S|, K ∈ N depends on the shape and is the number of eccentric point candidates that are evaluated.

Figure 2 shows a comparison of the geodesic and Euclidean distances. Figures 3 and4 show the eccentricity transform of a 2D, respectively 3D, shape. For the 3D shape, the eccentricity transform is presented for the whole shape (ECCobj), for the border voxels (ECCborder), and the surface mesh (ECCmesh). Figure 5 shows the difference between ECCobj and ECCborder, both using distances computed on a voxel description of S.

Eccentricity Histogram Matching

To match two shapes we first create a shape descriptor for each of them and then match these descriptors to obtain a similarity measure. The pro- posed shape descriptor is made out of two components: the first one characterizes the compactness (geometry) of the shape (Section 3.1) and the second one its structure (Section 3.2). The descriptor is highly compact, which is an advantage for real time retrieving and low memory devices, it is invariant under many natural deformations, it can handle shapes without as well as with holes4 (Figure 8 (g) and (h)), and gives good results, comparable to the presented state of the art methods (experiments follow in Section 4).

ECC histogram

The first component the shape descriptor is the histogram h of the eccentricity transform ECC of the shape S. We use k h bins for the histogram. The eccentricity histogram is the vector h ∈ R k h defined by: ∀i

= 1, . . . , k h h(S, i) = 1 |S| # p ∈ S | i -1 k h ECC(S, p) -m M -m < i k h ,
where |S| is the number of pixels/voxels/vertices in S, and m and M are the smallest, respectively, largest eccentricity values over S. A discussion about choosing the number of bins k h follows in Section 4.4. The obtained histogram only contains bins for the values which exist in the eccentricity transform, i.e. from minimum to maximum eccentricity, and the sum over all bins is 1. Figures 6 and7 show examples of eccentricity histograms for 2D and 3D shapes with different geometric features. Note that the histogram h is invariant under Euclidean transformations, scaling and isometric bending of S (Figure 8 shows examples).

ECC structure histogram

The second component of our shape descriptor is a histogram of the number of connected components of the discrete level-sets of the eccentricity transform. A discrete level-set of ECC(S) is the set of points of the shape having their eccentricity value in a certain domain (a, b]. More formally:

L(a, b) = {p ∈ S | a < ECC(S, p) b}.
The eccentricity structure histogram is the vector s ∈ R ks defined by: ∀i

= 1, . . . , k s s(S, i) = # C {L(e i-1 , e i )} ,
where # C denotes the number of connected components. We use 8 connectivity in 2D, 26 connectivity for the 3D voxel representations, and vertex adjacency for the triangle meshes. For a function which is based on Euclidean distances, to obtain at least 8 respectively 26 connected isolines/isosurfaces, the thresholds e i have to satisfy that e ie i-1 g [START_REF] Andres | The discrete analytical hyperspheres[END_REF]. Where g = 1 is the distance between two closest grid points. For the case of the triangle meshes, we use g = √ 2, which is the length of the longest edge in our meshes. If e ie i-1 < g then the pixels covering the points of the same (connected) continuous isoline will not be connected in the discrete level-sets. The values for e i are computed as:

e i = m + (M -m) i k s ,
where as before m, M are the minimum respectively maximum eccentricity.

We take e 0 = m -, > 0 to ensure that the geodesic center is included in the first bin. If (Mm)/g < k s , which means that e ie i-1 < g, we compute s for k s = (Mm)/g and then resize the obtained histogram to the required k s bins using bilinear interpolation.

Comparison of histograms

To match the descriptors of two shapes S and S, it is necessary to compute the distance between the corresponding histograms. Let h, s and h, s be the two histograms of S respectively S computed as in Sections 3.1 and 3.2. The distance between two histograms h, h ∈ R k is measured using the L 2 -norm:

δ(h, h) def. = h -h . ( 3
)
One could use more elaborate metrics such as the χ 2 statistic [START_REF] Belongie | Shape matching and object recognition using shape contexts[END_REF] or those defined in [START_REF] Osada | Shape distributions[END_REF]. In numerical experiments we found that all these metrics give results similar to δ, which is the easiest and fastest to compute (discussion follows in Section 5).

The dissimilarity Δ(S, S) of two shapes S and S is computed as a weighted distance of their eccentricity histograms h, h and structure histograms s, s:

Δ(S, S) def. = w • δ(h, h) + (1 -w) • δ(s, s) min( s , s ) , ( 4
)
where w is a mixing parameter. The histograms h ∈ R k h lie in the k h -1 dimensional simplex with vertices h r = (h 1 , . . . , h k h ), h i = 1 for i = r and 0 otherwise. Thus we have 0

δ(h, h) √ 2 for any h, h ∈ R k h .
The values of the bins s(S, i) 1 and can be arbitrarily large depending on the number of "parts" of the object. The expression min( s , s ) √ k s , and δ(s, s) 0. For similar shapes (same or similar looking classes) both terms in Equation 4have values around 0 (in the same range), and extend the descriptive power of the other one. Matching very different shapes (e.g are horses more similar to scissors then to cars?) is not considered.

In the remainder of the paper, we will use eccentricity histogram to denote h, structure histogram to denote s, and eccentricity based histograms to refer to both.

Characteristics of the ECC based histograms: from S to h, s and back

The eccentricity histogram h characterizes the compactness of the shape (e.g. a flat histogram characterizes a very elongated shape, a histogram with monotonically increasing values characterizes a rather compact shape). The structure histogram s characterizes the structure of the shape by representing the evolution of the number of parts of the shape when going from the geodesic center towards larger eccentricity values (e.g. for a spider the number will range from 1 to 8, and for an ant from 1 and 6). Intuitively, one could say that h looks at the widths of parts and s at their number. Figure 8 shows the ECC based histograms for basic shapes, with and without holes and articulation.

Effects of basic changes of the shape to the eccentricity histogram 5 h. The

q 3 q 3 ECC S c q 1 h i d(e 1 , q 1 ) d(e 1 , q 3 ) (c) S c e 1 e 2 c q 1 q 2 q 3 q 2 q 3 ECC S c q 1 q 2 h i d(e 1 , q 1 ) d(e 1 , q 3 ) (d) S d e 1 e 2 c c q 1 q 2 ECC S c, c q 1 q 2 h i d(e 1 , q 2 )
Figure 9: Behavior of eccentricity histograms for basic changes in the shape. Column S: where possible, straight lines where used for illustration, but only the length of the curves is relevant, not whether they are straight or not. For the shown shapes (composed of curves) s(S) = h(S) as each point of S contributes with "1" to both of them: one eccentricity value and one connected component. (d) thick line in ECC(S): between q 1 and q 2 two points have the same eccentricity.

histogram of the ECC of a simple open curve6 S a with length l = d(e 1 , e 2 ) (Figure 9(a)) is flat with a possibly smaller value in the first bin. The continuous formula is: Consider S b obtained by adding a simple open curve of length d(q 1 , q 3 ) < l/2 connected at the point q 1 to S a (Figure 9(b)). Let q 3 ∈ S b s.t. d(q 1 , q 3 ) = d(q 1 , q 3 ) and d(q 3 , e 1 ) = d(q 3 , e 1 ). For the points with eccentricity between d(e 1 , q 1 ) and d(e 1 , q 3 ), the eccentricity histogram of S b has increased by 50% (there is one additional point having each of the values in the domain). A shape without cycles (e.g. S a , S b , S c ) has only one center point (ECC minimum) and the histogram value for the center is always one. All other histogram values can be changed by adding branches as above.

h(S a , i) = 1 l if i = min(ECC(S a ))
The histogram value for the center can be changed by introducing cycles. Consider S d obtained by adding a simple open curve q 1 c q 2 of length λ(q 1 , q 2 ) to S a (Figure 9(d)). The length d(e 1 , e 2 ) is kept the same and q 1 q 2 has the same length if going over c or c . Also d(e 1 , c) = d(e 1 , c ) = d(e 1 , e 2 )/2. Two center points exist (c and c ), and for the eccentricity values [d(c, e 1 ), d(q 2 , e 1 )) there is one additional point.

For a given histogram, the steps used to create S b and S d , can be iterated to grow the continuous shape (for geodesics computed along thin lines). For discrete shapes, the number of points is finite 7 , which limits the number of curves that can be put close to each other and do not intersect. If the maximum shape size (number of pixels/voxels) and the number of bins k h is fixed, not all (real valued) histograms can result as ECC histograms (it can also be seen as a discretization problem: the lower the resolution/maximum size, the higher the dependence between neighboring histogram bins).

Equivalence classes of ECC based histograms.

A histogram has a smaller dimension (in our case 1D) than the shape, and a whole class of shapes is projected into the same histogram. Two shapes S and S with the same eccentricity based histograms satisfy Δ(S, S) = 0, and are thus considered to be the same according to our matching algorithm. Consider S c in Figure 9(c) obtained from S a , similar to S b , but with two curves s.t. d(q 1 , q 2 ) = d(q 1 , q 2 ), d(q 2 , q 3 ) = d(q 2 , q 3 ), and d(q 1 , q 3 ) is equal in both S b and S c . Hence, the two shapes S b and S c have the same eccentricity based histograms h, s and cannot be differentiated using only the histograms.

When the following operations are applied to a shape they create a shape in the same equivalence class, i.e. having the same or similar histograms:

• scaling;

• rotation;

• isometric deformation that does not change the structure of the shape (e.g. moving a finger, without touching another finger);

• moving certain points: simultaneous thinning and thickening of two parts with the same eccentricity values (e.g. making one finger thinner while making another one thicker);

• moving certain parts: taking a part and attaching it in another point with the same eccentricity while keeping the eccentricity values the same (e.g. moving a finger from the left hand to the right hand, obtaining a human with 4 respectively 6 fingers on his hands, also S b and S c in Figure 9);

• certain changes in topology without disconnecting the shape, when the contact points in the shape without the hole8 have the same eccentricity value and no eccentric (longest) paths go over the created/destroyed connection (e.g. if considering the whole human as a shape, touching the tip of a finger with the neighboring finger, with the fingers as straight as possible). This change affects at least one bin of s, but could be only minimally (one point contact as opposed to keeping two fingers touched along the whole length). Disconnected shapes cannot be handled by our descriptor.

One can say that the eccentricity based histograms are influenced by the geometry and the structure of shapes, but they do not uniquely characterize them. This can be both a positive or a negative aspect, as one would like filter out for example noise, but probably not "attaching a second head".

Matching Experiments in 2D and 3D

This section shows results on popular benchmarks and comparison with state of the art methods. When comparing the results, keep in mind that the proposed method is simple and matching the computed descriptors is fast. An approximation of the ECC can be computed for many shapes with as few as 50 distance propagations (e.g. the average number for the ECCmesh on the McGill database is 54), and determining δ between two computed descriptors (L 2 -norm) has practically no CPU time consumption 9 . Two fixedlength vectors as descriptors can be a very efficient indexing method. The approaches compared with are more complicated and include the decomposition of shapes, aligning or finding correspondences between features, etc.

2D Shape Matching

For the experiments with 2D shapes we have used four shape databases: Kimia 25 [START_REF] Sharvit | Symmetry-based indexing of image databases[END_REF] (Figure 13), Kimia 99 [START_REF] Sebastian | Recognition of shapes by editing their shock graphs[END_REF] (Figure 10), MPEG7 CE-Shape-1 [START_REF] Latecki | Shape descriptors for non-rigid shapes with a single closed contour[END_REF] (Figure 11), and the articulated shape database of [START_REF] Ling | Shape classification using the inner-distance[END_REF] (Ling articulated) (Figure 12). As a baseline we have added the results obtained by using the Shape Index, which is computed as the length of the boundary of the shape over its area.

A shape database is composed of q shapes {S i } q i=1 and each shape S i has a label L(i) ∈ {1, . . . , l max }. Each label value 1 l l max defines a class of shapes Q(l) = {S i | L(i) = l}. The first columns of the blocks in Figure 13 show the shapes from the Kimia 25 database, ordered by classes (such as fish, planes, rabbits, etc.). Any shape matching algorithm α assigns to each shape S i a vector of best matches Φ i , where Φ i (1) is the shape the most similar to S i , Φ i (2) is the second hit, and so on. Depending on the benchmark, Φ i contains all shapes including the query shape S i (MPEG7), or leaves S i out, i.e. the shape S i is not matched to itself and Φ i has q -1 elements (all other benchmarks presented).

For the Kimia 25 database l max = 6 and q = 25, and for the Kimia 99 database, l max = 9 and q = 99. The efficiency of various matching algorithms on Kimia databases is measured by the number of correct matches for each ranking position r:

Match r (Φ) def. = q i=1 1 L(Φ i (r))=L(i) q.
(5)

Tables 2 and3 give the value of Match k for various shape matching algorithms. We also show results when using only one of the two: h or s.

In the case of the MPEG7 database, which contains l max = 70 classes with 20 images each (q = 70 × 20 = 1400), the efficiency of matching algorithms is computed using the standard Bullseye test:

Bullseye(Φ) def. = 1 20q 40 r=1 q i=1 1 L(Φ i (r))=L(i) = 1 20q 40 r=1 Match r (Φ). (6)
This test counts the number of correct hits (same class) in the first 40 hits.

For each image there can be at most 20 correct hits and a maximum of 20×1400 hits can be obtained during the benchmark and thus Bullseye(Φ) 1. Table 4 gives the value of Bullseye for various shape matching algorithms. The Ling articulated shape database was created specifically to show results in the presence of articulation. The database contains q = 40 shapes of articulated objects from l max = 8 classes. The efficiency of an algorithm is measured by the number of correct matches Match r (Φ) for each ranking position r (Equation 5). Table 5 gives the value of Match k for various shape matching algorithms.

Case study -Kimia 25. Figure 13 shows the retrieval results for Kimia 25 when using only h, only s or both terms. The first columns show the 25 shapes S i (1-4 fish, 5-8 greebles, 9-13 hands, 14-17 airplanes, 18-21 rabbits, 22-25 tools). The following columns r show Φ i (r), the rank-r shape associated to S i . 6) for various algorithms on the MPEG7 database (q = 1400 shapes from l max = 70 classes). Algorithm α Bullseye(Φ) random 2.86% Shape Index 25.46% ECCobj2D 54.56% Shape Context [START_REF] Belongie | Shape matching and object recognition using shape contexts[END_REF] 64.59% ID-Shape Context [START_REF] Ling | Shape classification using the inner-distance[END_REF] 68.83% The Kimia 25 database has shapes from 6 classes: 5 classes with 4 images each, and one (hands) with 5 images (1 simulating a segmentation error). If considering only h (Figure 13 top) the class with the best results are rabbits, followed by tools, hands, fishes, airplanes and greebles. Two questions immediately rise when looking at these results:

1. Why are the greebles considered to be more similar to the hands than to other greebles? 2. Why does a rabbit appear in so many cases when the matching has failed?

For the first question, consider the ECC histograms h of the greebles and the not occluded hands (Figure 14). The histograms are similar even though the shapes are of different classes, e.g. the histogram of the first greeble (Figure 14 top-left) looks more similar to the hands, than the second and third greeble. This is due to the abstraction of a 2D shape to a 1D histogram which, in the case of h, disregards certain structural properties of distances/paths (studied in detail in Section 3.4). The structure histogram s can better discriminate between these shapes, and adding it to the descriptor compensates for this effect.

For the second question, consider the shapes in Figure 15 (a rabbit -S 19 , and two tools -S 25 and S 22 ). When matching S 25 , the rabbit has a smaller distance of h than S 22 , even though one might say that the histograms h of S 25 and S 22 reveal more similar distance characteristics than the histogram h of S 19 (see Figure 15). Both S 25 and S 22 have more long distances than medium, and short, while S 19 has a peak in the medium. This effect is due to typical histogram matching methods, which are inherently low level and fail to capture the high level context of the task. Note that similar examples could be constructed also using s. Discussion follows in Section 5. Global geometrical and structural statistics of the shapes are well captured by our low dimensional descriptors. However, more detailed structural properties like part decomposition, ordering and adjacency of parts, and the geometrical features corresponding to these parts are not considered by our signature extraction.

3D Articulated Shape Matching

A widely used 3D object retrieval database is the Princeton Shape Benchmark [START_REF] Shilane | The princeton shape benchmark[END_REF]. It contains 1,814 3D object models organized by class and is effective for comparing the performance of a variety of methods. However, the majority of the models corresponds to rigid, man-made objects. Only a limited number of shapes in the database have articulated parts. As one of the main advantages of using the eccentricity transform is its robustness with respect to articulation, we have turned to the McGill Shape Benchmark [START_REF] Zhang | Retrieving articulated 3-D models using medial surfaces and their graph spectra[END_REF]. It contains several models from the Princeton repository and others added by the authors. The main advantage of this benchmark is that from its 455 3D shapes, 255 have significant part articulation. We show the results on the q = 255 shapes grouped into the l max = 10 classes of articulated shapes (Figure 16). Shapes are not matched to themselves and so Φ i contains q -1 shapes.

We compare results with:

• medial surfaces (MS) [START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF];

• spherical harmonic descriptor (HS) [START_REF] Kazhdan | Rotation invariant spherical harmonic representation of 3D shape descriptors[END_REF];

• shape distributions (SD) [START_REF] Osada | Shape distributions[END_REF].

MS computes the medial surfaces of a voxelized shape and decomposes them in parts. The similarity of two shapes is obtained from the matching of the directed acyclic graphs describing their parts. HS transforms rotation dependent spherical shape descriptors into rotation invariant ones. It describes spherical functions on the shape in terms of the amount of energy they contain at different frequencies. The L 2 -norm is used to match two descriptors. SD computes the similarity of two shapes by the comparison of two probability distributions sampled from a shape function measuring geometric properties of the 3D model. Best results are obtained with a function called D2 which represents the distribution of Euclidean distances between pairs of randomly selected points on the surface of a 3D model.

Three ECC based descriptors are evaluated (Figure 4):

ECCobj

ECCborder ECCmesh

ECCobj

ECCborder ECCmesh 1. ECCobj -eccentricity of the whole shape (all object voxels); 2. ECCborder -eccentricity of the border/boundary voxels; 3. ECCmesh -eccentricity of the triangular mesh of the surface of the shape.

ECCborder is the eccentricity transform ECC(∂ 6 S), where ∂ 6 S is the 6 connected voxel boundary of S. ECCmesh is computed on the 2D manifold defined on the boundary of the 3D shapes. ECCborder uses distance computation in the 3D volume, ECCmesh in the 2D surface. If the resolution of the shapes is increased, ECCborder and ECCmesh converge to the same value. For a similar resolution, ECCmesh needs less memory, as cells not part of the boundary do not have to be stored (e.g. interior of the shape), and it can be more accurate when approximating the eccentricity of the surface, as the computation is done on the surface itself, not on an approximating (thicker) volume.

The size of the voxel models used for ECCobj and ECCborder is below 128 × 128 × 128 (the size of the binary 3D images is 128 × 128 × 128).

Experimental results. In the following, results of the three variants on shapes of the 10 articulated classes of the McGill Shape Benchmark are given. The notation from Section 4.1 is used. The following measures are considered.

Recall(Φ i , t) = 1 |Q(L(i))| -1 t r=1 1 L(Φ i (r))=L(i)
The recall computes the ratio of models in the database in the same category as the query, with indexing rank t, to the total number of shapes in the same category (never including the query itself). The average results and standard deviation for several rank thresholds (t = 10, 20, . . . ), over all classes, are given in Figure 17.

AvgRank(Φ i ) = 1 |Q(L(i))| -1 q-1 r=1 r • 1 L(Φ i (r))=L(i)
For all queries in a class, the average of the ranks of all other shapes in that class are computed. Figure 18 shows the average and the standard deviation of the ranks for each class (lower average is better).

Table 6 shows the average score for all pairs of classes. Each shape in the database is matched against all other shapes and each cell shows the average of the score (Equation 4) between all combinations of shapes of the two classes defined by the row and column.

Precision(Φ i , t) = 1 t t r=1 1 L(Φ i (r))=L(i)
Precision refers to the ratio of the relevant shapes retrieved, to the total number retrieved. Figure 19 shows the precision-recall curves for each of the 10 classes. Precision-recall curves are produced by varying the parameter t.

Better results are characterized by curves closer to the top, i.e. recall = 1 for all values of precision. Precision and recall are common in information retrieval for evaluating retrieval performance. They are usually used where static document sets can be assumed. However, they are also used in dynamic environments such as web page retrieval [START_REF] Fawcett | An introduction to ROC analysis[END_REF].

As can be seen in Figures 17,18, and 19, and Table 6, ECCobj does in most cases a better job than ECCborder and ECCmesh. The recall of the three methods is very similar, with slightly better results from ECCobj (noticeable on the first three bars). With respect to the average ranks, EC-Cobj does better with the ants, octopus, spiders, teddy, is equal to one of ECCborder and ECCmesh with the crabs, humans, pliers, snakes, and is slightly worse than one of the other methods with the humans and spectacles. None of the three variants produces an average class rank higher than 21% of 255 (average rank 42, 53, and 52, for the octopus, for ECCobj, ECCborder, and ECCmesh, respectively). ECCobj has the smallest average class distance (highest similarity) correct for all 10 classes. For both EC-Cborder and ECCmesh the smallest average class distance is correct for 9 classes, while the correct class is the second smallest one for the remaining class, the octopus (see Table 6). A discussion considering the differences between the three ECC variants follows at the end of this section.

Figure 19 shows comparative precision-recall results of ECCobj, ECCborder, ECCmesh, MS, HS, and SD. ECCobj, ECCborder, and ECCmesh are close, except for the teddy bears, where ECCobj is superior to the other two. The best results (higher precision vs. recall) are reached by the ECC variants for the snakes, by MS for the crabs, and HS and SD for teddy. For these best results the MS and the ECCobj have almost identical precision-recall, followed by HS and SD (Figure 20). The worst results are achieved by EC-Cobj, ECCborder, and ECCmesh for the octopus, MS for the pliers, HS for snakes, and SD for the hands. The results of MS for the pliers are superior to ECCobj for the octopus, which are in turn superior to the HS, superior to SD (Figure 20). In comparison to all other three methods (MS, HS, SD), the eccentricity based methods score better on the pliers, spectacles and snakes, and ECCobj is on par with the best on ants, hands, and spiders. MS searches for a concrete mapping between the parts of articulated shapes, which explains the overall superior results. The ECC variants are articulation insensitive and able to characterize the compactness (geometry) of the shapes quite well. However, describing the structure of the shape using a histogram does not capture the connections between the parts (e.g. one could detach a finger from a hand and attached it to the other one, without changing the descriptor of the human). The descriptor computed by HS is invariant to the independent rotation of parts, even if the shape gets disconnected. SD uses a function which is based on Euclidean distance and not invariant to articulation.

The differences in the results of ECCobj vs. ECCborder and ECCmesh, can be linked to the compactness of the shapes and the width of their joints. During articulation, the variation of the geodesic distances is larger when computed on the 'skin' (boundary) compared to computed inside the shape (smaller). In the case of 2D shapes the eccentricity of the boundary is a constant. In 3D it manages to capture some of the properties of the shape, but it is more unstable, e.g. the eccentricity transform of a simply connected volume has in most cases a single stable center (minimum), while the eccentricity transform of its boundary will have a disconnected center or at least one with a more complex structure. A more concrete example of such instability with only minor deformation of the shape is the ECC of the surface of:

• a sphere: all points make up the center;

• an ellipsoid with two equal radiuses and a slightly longer third one: the center is one circle.

As a last experiment we show results of ECCmesh on the Robustness Benchmark from the SHREC'10 dataset [START_REF] Bronstein | Shrec 2010: robust large-scale shape retrieval benchmark[END_REF] 10 . The dataset consists of 1184 shapes given as triangle meshes. The query set consists of 13 shape classes taken from a subset of the dataset with simulated transformations applied to them. For each shape, transformations are split into 10 classes (isometry, topology, small and big holes, global and local scaling, noise, shot noise, partial occlusion, sampling, and a combination of all transformations). In each class, the transformation appears in five different strength levels. The total number of transformations per shape is 55 plus the null shape (neutral pose), and the total query set size is 728 shapes. The reported results are mean average precision (mAP).

The given triangle meshes are possibly disconnected and non-manifold (e.g. more than 2 triangles can share an edge). As the ECC is computed on a single connected component and our geodesic distance computation requires manifolds, for each shape we select the largest subset of triangles (based on area), which are connected and approximate a manifold.

Table 7 gives the results of ECCmesh on the SHREC'10 robustness benchmark. The transformation "partial" reduces the shape to a disconnected subset of the original triangles, consisting of many components. Our method cannot handle such disconnected shapes and computes the descriptor of one of the connected mesh parts instead. Thus we also give the average score over only the transforms which keep the query shape connected (without "partial" or "all"). When computing the shape descriptors we did not refine the meshes or sample from them uniformly, but took the ECC values for the given vertices of the mesh. The resolution of the meshes is not uniform and Table 7: Retrieval results on the SHREC'10 Robustness Benchmark (total q = 1184 shapes, with a query set of 728 shapes from l max = 13 classes). The reported performance is mean average precision. Rows indicate the type of transform that the query shape was altered with, columns indicate the maximum strength of the applied transform. "avg": results over all queries in the dataset. "avg (connected)": results only with connected shapes (without transformations "partial" and "all"). explains why the results for "scale" are not 100%, as expected from the normalization. The sensitivity of the histograms to the non-uniform sampling of the vertices is reflected in the results for the transformation "sampling" 11 . For an evaluation of different methods on the SHREC'10 benchmark see [START_REF] Bronstein | Shrec 2010: robust large-scale shape retrieval benchmark[END_REF].

Discussion: 2D and 3D

The computed shape similarities are robust with respect to scaling, rotation, and part articulation. The matching results are good, especially when considering the straightforward approach, and databases of shapes undergoing articulation. Our histogram representation is invariant to articulation, thus achieving state of the art performance on these databases. In contrast, the most efficient shape matching algorithms [START_REF] Ling | Shape classification using the inner-distance[END_REF][START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF] are more complicated and require extraction of salient features and local signatures that need to be aligned or registered. E.g. the method in [START_REF] Ling | Shape classification using the inner-distance[END_REF] describes shapes by a collection of local descriptors sampled at different locations on one selected boundary of the shape. Each descriptor is a 2-dimensional histogram, and is based on the computation of geodesic distances and the angle made by the tangent to the geodesic with the tangent to the shape boundary at the point of the descriptor. Similarities are computed using dynamic programming as the cost of the best matching of the descriptors along the boundary. In [START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF], first a medial surface is obtained by doing topology preserving thinning on the average outward flux computed on the distance transform of the shape. The obtained medial surface is segmented into components and using saliency a directed acyclic graph (DAG) is built to describe the shape. To match a query DAG with a database, first an indexing strategy is used to quickly select a few candidates, followed by a bipartite graph matching to compute node correspondences and the matching score. In addition to the matching score both methods [START_REF] Ling | Shape classification using the inner-distance[END_REF][START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF] can provide correspondences between the used parts/points. Also the higher quality results of [START_REF] Ling | Shape classification using the inner-distance[END_REF] can be explained by the usage of local descriptors together with a score computation which considers their order -this gives more detailed geometric information with a higher tolerance to missing parts.

The limitations of our approach include: (1) Eccentricity structure histograms s are not able to fully capture the part structure of the shape and thus histograms of different shapes can be very similar. (2) Histogram 'matching' (whether using the L 2 -norm or more sophisticated methods) is inherently low level and does not consider the higher level context in which it is applied.

One can identify the limitations discussed before (see Figure 19 and Table 6): for classes with simple topology (e.g. snakes and spectacles), the results are very good. For classes where part decomposition and detailed geometry play an important role (e.g. octopus, hands, humans), the discrimination capabilities are reduced.

The Parameters k h , k s and w

The approach has three parameters: the number of bins k h , k s , of the histograms h, s, and the mixing factor w. In all presented experiments 2D and 3D, k h = 200 and k s = 100 where used. The values where chosen based on a few initial trials on a small set of shapes (Kimia25). Table 8 shows the results on Kimia99 for different values of the parameter k h when using only h as a shape descriptor, i.e. when w = 1. For s, a difference of at least 1 (2D, and voxel representations in 3D) respectively √ 2 is required between the centers of two neighboring bins. Thus any value of k s larger than half the average of the geodesic diameters of the shapes is enough. The used value for w was 0.94 and was selected based on a few trials on the Kimia25, Kimia99, and Ling databases (note that h and s use a different normalization scheme, which explains the apparently large bias toward h). Table 9 shows results for different values of w.

As the shapes are discrete, the number of distance values of the ECC is finite. Let h c be the ordered set of eccentricity values computed for a shape S, i.e. each distinct value that exists in the ECC of the discrete shape S. We have min(h c ) equal to the ECC value of the center (minimum ECC) and greater or equal to half the diameter of the shape (max(h c ) = max(ECC)). The largest distance between two neighboring (grid) points is equal to one (shapes are required to be 4 respectively 6 connected). For the ECC histogram of a shape not to contain any empty bins, the number of bins k h has to satisfy: k max(ECC(S))min(ECC(S)). Depending on the shape, k h could be much higher and still have no empty bins in h, e.g. for S a disk with radius u in Z 2 and the Euclidean distance, there are more distinct values than u (consider the discrete approximation of the Euclidean circle). An absolute upper bound is k h = |S|. If this number is exceeded, there will be empty bins in h.

As k h decreases, the description capability of the histogram also decreases. In the extreme case, a single bin would just contain |S|, and for the normalized histogram it would contain the value 1. Two bins can give the equivalent of a simple compactness measure (similar to the circularity ratio or shape index ), which relates the area of the shape to the area of the circle with the same diameter). Three bins could be considered as a relative measure for short/long/medium distances and can characterize more than the simple compactness measure.

A higher number of bins increases the dimension of the space in which distances are computed, and gives more flexibility in the relations, e.g. in 2D there can be maximum 3 points s.t. they are pairwise at the same Euclidean distance (equilateral triangle), and this number increases to 4 in 3D (regular tetrahedron). Assuming that the number of classes and their pairwise (average) distances are known, a lower bound for the number of bins is equal to the smallest dimension in which the classes can be embedded s.t. the distances computed in the histogram space are equal or close to the given distances.

Potential Extensions

Topological changes in the shape. Geodesic distances are sensitive to changes in the topology of the shape. For example, it is enough to touch the index fingers of two hands in one point to drastically change the geodesic distance between the points of the two palms: before we had to go over the arms and torso, now we only have to travel over two fingers. This problem has been approached in [START_REF] Rustamov | Laplace-beltrami eigenfunctions for deformation invariant shape representation[END_REF][START_REF] Bronstein | A gromov-hausdorff framework with diffusion geometry for topologically-robust non-rigid shape matching[END_REF] by using diffusion distance as an alternative to geodesic distance. Diffusion distances consider all possible paths by which two points can be connected and are expected to be more stable to changes in topology (in the previous example, the high number of paths going through the thicker arms will still be dominant over the path(s) going over the small size contact area between the two fingers). The definition of the eccentricity transform (Equation 2) is defined over a metric space and diffusion distances could also be used.

Describing and matching part structure. One of the problems identified in Section 3.4 and during the experiments (Sections 4.1 and 4.2) is that the histograms s do not capture the exact structure/topology of the shape. Classical methods to describe the topology of a shape (e.g. Reeb graphs [START_REF] Reeb | Sur les points singuliers d'une forme de pfaff complément intégrable ou d'une fonction numérique[END_REF], and homology generators [START_REF] Munkres | Elements of Algebraic Topology[END_REF]) do not capture geometrical aspects. An approach to deal with this problem is presented in [START_REF] Aouada | Geometric modeling of rigid and non-rigid 3D shapes using the global geodesic function[END_REF]. Two descriptors are used to describe a shape: a geometric one, based on the Global Geodesic Function (GGF) [START_REF] Hilaga | Topology matching for fully automatic similarity estimation of 3D shapes[END_REF], which is defined for a point as the sum of the geodesic distances to all points of the shape multiplied by a factor, and a topological one, the Reeb graph of the shape using the GGF as the Morse function.

Initial steps in combining the eccentricity transform with Reeb graphs have been presented in [START_REF] Ion | Eccentricity based topological feature extraction[END_REF].

A better histogram matching. The problem of having a matching function that is aware of the context in which it is applied can be approached in two ways: (1) use expert knowledge about the context to create an algorithm that considers the proper features, or (2) learn the important features by giving a set of representative examples (e.g. [START_REF] Bronstein | Shape google: geometric words and expressions for invariant shape retrieval[END_REF]). In [START_REF] Yang | Distance metric learning: A comprehensive survey[END_REF][START_REF] Yang | An efficient algorithm for local distance metric learning[END_REF], a survey of current distance metric learning methods is given. The purpose of distance metric learning is to learn a distance metric for a space, from a given collection of pairs of similar/dissimilar points. The learned distance is supposed to preserve the distance relation among the training data. Example training data would be: S 1 is more similar to S 2 than to S 3 . The result is a distance function that would replace the L 2 -norm in Equation 3with a new measure, which is adapted to the task of computing the distance of eccentricity histograms as given by the training examples.

Partial shape similarity. The presented approach matches "whole" shapes and will work in the presence of only minor occlusion. Even though the eccentricity values remain the same when certain parts of a shape are removed, removing parts will change the histograms, proportional to the size of the removed parts. Partial shape similarity could be approached by first decomposing a shape into parts (like in the approaches in [START_REF] Siddiqi | Shock graphs and shape matching[END_REF][START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF]) and then matching the parts to each other. This way one can obtain part correspondences and higher robustness with respect to occlusion.

Higher dimensional data. 4D data has started to be available in the medical image processing community (e.g. 3D scans of a beating heart, over time). The presented method is general and should work in any metric space. This includes 4D, but also gray scale images (e.g. gray values can determine the distance propagation speed in the respective cells). A study in this direction is planned.

Conclusion

We have presented a method for matching 2D and 3D shapes. The method is based on the eccentricity transform, which uses maximal geodesic distances and is insensitive to articulation. Descriptors are composed of two terms: a normalized histogram of the eccentricity transform and a histogram of the connected components of the level-sets of the eccentricity transform. They characterize the compactness and structure of the shape, are compact and easy to match. The method is straightforward but still efficient, achieving state of the art results on databases of shapes undergoing articulation. Articulation is indeed the main invariance provided by our histogram representation. Experimental results on popular 2D and 3D shape matching benchmarks are given, with computation on binary 2D images, binary 3D voxel shapes, and 3D triangular meshes. The parameters of the method, and the relation between changes to a shape and its corresponding eccentricity Algorithm 1 ECC06(S) -Estimate eccentricity transform by progressive refinement. Input: Discrete shape S.

1: for all q ∈ S, ECC(q) ← 0 /*initialize distance matrix */ 2: p ← random point of S /*find a starting point*/

Figure 2 :Figure 3 :Figure 4 :Figure 5 :

 2345 Figure 2: Euclidean (a) and geodesic (b) distance function, for starting point p 0 . Gray values are distances modulo a constant.

Figure 6 :

 6 Figure 6: Top: ECCobj2D for some 2D shapes. Middle and bottom: corresponding histograms h respectively s.

Figure 7 :

 7 Figure 7: Top: example 3D shapes. Middle and bottom: corresponding ECCobj histograms h respectively s.

Figure 8 :

 8 Figure 8: Basic shapes and their histograms h, s. (vertical axes independently scaled)

2 l

 2 if i > min(ECC(S a )) , where min(ECC(S a )) = d(e 1 , c) = d(e 2 , c).

Figure 9

 9 Figure9).6 The term curve is used to denote a one dimensional and continuous manifold, and includes both straight and non-straight lines.7 Depends on the discretization and maximum shape size.

Figure 10 :

 10 Figure 10: The 99 shapes from the Kimia 99 shape database. Shapes in the same row belong to the same class.

Figure 11 :

 11 Figure 11: One shape from each of the 70 classes in the MPEG7 shape database.

Figure 12 :

 12 Figure 12: The 40 shapes from the Ling articulated database. Shapes in the same column belong to the same class.

Figure 13 :Figure 14 :

 1314 Figure 13: Retrieval results on Kimia 25 using top: h, middle: s, bottom: both h and s.Each row shows the query shape, followed by the first three matches.

Figure 15 :

 15 Figure 15: Three shapes from the Kimia 25 database and their eccentricity histograms h.

Figure 16 :

 16 Figure 16: The object classes from the McGill 3D shape database having significant part articulation. The number of instances in each class are given.

Figure 17 :

 17 Figure 17: Recall for several rank thresholds

Figure 18 :

 18 Figure 18: Average ranks for each class. The first three letters of each class name are printed.

Figure 19 :

 19 Figure 19: Precision-recall for the articulated shapes from the McGill dataset (q = 255 shapes from l max = 10 classes). Left two columns: ECCobj, ECCborder, ECCmesh. Right two columns (image taken from[START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF], with kind permission of Springer Science and Business Media): results of three other methods on the same database: medial surfaces (MS)[START_REF] Siddiqi | Retrieving articulated 3-D models using medial surfaces[END_REF], spherical harmonic descriptor (HS)[START_REF] Kazhdan | Rotation invariant spherical harmonic representation of 3D shape descriptors[END_REF], and shape distributions (SD)[START_REF] Osada | Shape distributions[END_REF]. Precision: horizontal axis, recall: vertical axis. (Best visualized in color.)

Figure 20 :

 20 Figure 20: Best and worst precision-recall curves for the results of ECCobj, MS, HS, and SD (Picture synthesized from Figure19). (Best visualized in color.)

3 : 4 : 6 :

 346 /*Phase 1: find a diameter */ 5: while p not computed do ECC ← max{ECC, D S (p)} /*accumulate & mark p as computed */ 7: p ← arg max{ECC(p) | p ∈ S} /*highest current ECC (farthest away)*/ 8: end while 9: 10: /*Phase 2: find center points and local maxima*/ 11: pECC ← 0 /*make sure we enter the loop*/ 12: while pECC = ECC do 13: pECC ← ECC14: C ← arg min{ECC(p) | p ∈ S} /*points with min. ECC */ 15: for all c ∈ C, c not computed do 16: D ← D S (c) /*compute distances from the center */ 17: ECC ← max{ECC, D} /*accumulate & mark c as computed */ 18: 19: M ← {q ∈ S | D(q) local max. in S & q not computed} 20: for all m ∈ M, m not computed do 21: ECC ← max{ECC, D S (m)} /*accumulate & mark m as computed *Distances ECC.

Table 1 :

 1 Types of manifolds used for matching.

	Name	input computing on	S (d ε is used)
	ECCobj2D	2D	2D: whole shape	4-connected binary 2D shape
	ECCobj	3D	3D: whole shape	6-connected 3D voxel shape
				6 connected voxel surface in 3D, made
	ECCborder	3D	3D: border voxels	out of voxels of the shape that are 26
				connected to a background voxel
	ECCmesh	3D	2D: triangular mesh	connected triangular mesh of the sur-face of the 3D shape

Table 2 :

 2 The value of Match r (Φ) for various algorithms on the Kimia 25 database (q = 25 shapes from l max = 6 classes).

	Algorithm α	r=1 2	3
	Shape Index	14	10 10
	ECCobj2D, s only	18	19 17
	ECCobj2D, h only	20	16 14
	ECCobj2D	22	20 17
	Sharvit et al. [62]	23	21 20
	Gdalyahu et al. [64]	25	21 19
	Shape Context [7]	25	24 22
	ID-Shape Context [8]	25	24 25

Table 3 :

 3 The value of Match r (Φ) for various algorithms on the Kimia 99 database (q = 99 shapes from l max = 9 classes).

	Algorithm α	r=1 2	3	4	5	6	7	8	9 10
	Shape Index	43	51 58 52 52 49 51 45 44
	ECCobj2D, s only	84	68 65 67 56 57 51 41 31
	ECCobj2D, h only	87	74 66 64 49 52 45 33 33
	ECCobj2D	94	85 81 73 81 73 64 56 35
	Shape Context [7]	97	91 88 85 84 77 75 56 37
	Gen. Model [65]	99	97 99 98 96 96 94 75 48
	Shock Edit [52]	99	99 99 98 98 97 96 93 82
	ID-Shape Context [8]	99	99 99 98 98 97 97 94 79

Table 4 :

 4 The value of Bullseye(Φ) (Equation

Table 5 :

 5 The value of Match r (Φ) for various algorithms on the Ling articulated shape database (q = 40 shapes from l max = 8 classes). See[START_REF] Ling | Shape classification using the inner-distance[END_REF] for a description of these algorithms.

	Algorithm α	r=1 2	3	4
	Shape Context+DP [8]	20	10 11 5
	Shape Index	24	17 17 18
	ECCobj2D	40	33 29 22
	ID-Shape Context+DP [8]	40	34 35 27

Table 6 :

 6 Average matching results multiplied by 100 (smaller means more similar). For each row, the first and second smallest value are printed in bold.

		ants	crabs	hands	humans	octopus	pliers	snakes	spectacles	spiders	teddy
		ants 2.99	7.74	7.93	8.15	8.31	7.90 15.73 19.49 7.25 10.57
		crabs 7.74	2.60	9.09	9.85	5.88 10.37 20.85 20.17 4.75 11.24
		hands 7.93	9.09	4.39 5.34 10.15 5.36 10.83 13.02 9.76	5.72
	ECCobj	humans 8.15 octopus 8.31 pliers 7.90 10.37 5.36 9.85 5.34 5.88 10.15 10.11 5.10 10.46 19.94 17.55 5.60 11.97 3.62 10.11 4.08 8.88 10.60 9.51 5.21 4.08 10.46 1.45 6.58 9.44 10.44 4.53 snakes 15.73 20.85 10.83 8.88 19.94 6.58 1.40 6.29 20.78 7.80
		spectacles 19.49 20.17 13.02 10.60 17.55 9.44	6.29 3.00 19.60 9.29
		spiders 7.25	4.75	9.76	9.51	5.60 10.44 20.78 19.60 2.98 11.91
		teddy 10.57 11.24 5.72	5.21 11.97 4.53	7.80	9.29 11.91 3.47
		ants 2.40	5.62	6.51	6.39	4.83	6.26 14.33 16.84 4.93	8.32
		crabs 5.62	2.58	7.83	9.27	5.78 10.14 19.65 22.26 4.71 11.20
	ECCborder	hands 6.51 humans 6.39 octopus 4.83 pliers 6.26 10.14 5.01 7.83 3.97 4.91 9.27 4.91 3.01 5.78 7.97 8.73 3.62 snakes 14.33 19.65 10.09 7.73 17.74 6.45 7.97 5.01 10.09 12.33 8.30 8.73 3.62 7.73 10.37 8.68 4.49 8.48 17.74 18.12 4.44 10.57 5.60 4.53 8.48 1.57 6.45 8.84 8.85 4.08 1.49 5.63 18.89 6.35 spectacles 16.84 22.26 12.33 10.37 18.12 8.84 5.63 2.60 19.70 8.92
		spiders 4.93	4.71	8.30	8.68	4.44	8.85 18.89 19.70 2.80 10.72
		teddy 8.32 11.20 5.60	4.53 10.57 4.08	6.35	8.92 10.72 3.74
		ants 2.45	5.90	5.86	6.31	4.52	6.01 13.76 15.62 4.77	7.31
		crabs 5.90	2.86	8.02	9.62	5.52 10.42 20.03 22.16 4.49 10.77
	ECCmesh	hands 5.86 humans 6.31 octopus 4.52 pliers 6.01 10.42 4.97 8.02 4.20 9.62 5.19 5.52 7.39 snakes 13.76 20.03 10.38 7.41 17.63 6.43 5.19 7.39 4.97 10.38 11.98 7.74 3.01 8.67 3.61 7.41 9.78 8.74 8.67 4.19 8.39 17.63 18.20 4.16 3.61 8.39 1.71 6.43 8.30 8.76 1.56 4.73 18.76 6.87 5.50 4.53 9.61 3.87
		spectacles 15.62 22.16 11.98 9.78 18.20 8.30	4.73 2.58 19.69 8.78
		spiders 4.77	4.49	7.74	8.74	4.16	8.76 18.76 19.69 2.65	9.75
		teddy 7.31 10.77 5.50	4.53	9.61	3.87	6.87	8.78	9.75	3.70

Table 8 :

 8 The value of Match r (Φ) for ECCobj2D, h only (w = 1 in Equation4), on the Kimia 99 database using descriptors with different number of bins. See Table3for comparison with other methods.

	number of bins k r=1 2 3 4 5 6 7 8 9 10
	10	85 69 66 59 52 51 40 46 42 30
	25	87 74 68 63 48 53 45 38 35 33
	50	87 74 67 68 45 51 43 38 37 31
	100	87 74 66 65 48 53 44 38 34 32
	200	87 74 66 65 49 51 45 39 33 33
	500	87 74 67 64 48 53 45 38 33 33
	1000	87 74 67 64 48 54 44 37 34 33
	2000	87 74 67 64 48 53 44 38 34 33

Table 9 :

 9 The value of Match r (Φ), r = 1, 2, 3 for ECCobj2D on the Kimia 25, Kimia 99, and Ling databases using k h = 200, k s = 100 and different values of w. See Section 4.1 for comparison with other methods.

	database	w	r=1 2 3
	Kimia 25 0.93 22 19 16
	Kimia 99 0.93 95 84 79
	Ling	0.93 40 33 27
	Kimia 25 0.94 22 20 17
	Kimia 99 0.94 94 85 81
	Ling	0.94 40 33 29
	Kimia 25 0.95 22 20 18
	Kimia 99 0.95 93 86 81
	Ling	0.95 40 33 29

Called 'inner-distance' in[START_REF] Ling | Shape classification using the inner-distance[END_REF].

Known in the mathematical morphology community as the propagation function.

Also called geodesic distance function with marker set p.

As opposed to methods selecting a single boundary, usually the outer/longest one, which cannot differentiate for example a disk from a 2D torus.

For the example class of shapes composed of 1D curves, h(S) and s(S) are equal (see

The term curve is used to denote a one dimensional and continuous manifold, and includes both straight and non-straight lines.

Depends on the discretization and maximum shape size.

Points that will be connected to make the hole appear.

For example, on an Intel Xeon with 3 GHz, computing ECCmesh for the 255 shapes of the McGill database (see the Appendix) takes around 2.5 hours in total i.e. an average of 35 seconds per shape (code is partly matlab, partly c++). Computing Δ for all pairs of descriptors (255*254 matches) takes less than 1 second in total.

See http://tosca.cs.technion.ac.il/book/shrec_robustness.html for images of the shapes in the dataset.

Refining the meshes to obtain a uniform sampling for the locations of the vertices could in principle be added as a pre-processing step.

histogram are discussed in detail. Adding extended structural information to the descriptor and searching for correspondences between parts would help to overcome the cases where shapes of different classes have similar histograms. Matching can be improved by adding context in the histogram similarity computation e.g. through metric distance learning.

Acknowledgments

The presented work was partially supported by the Austrian Science Fund under grants S9103-N13 and P18716-N13. Adrian Ion was supported, in part, by the European Commission, under project MCEXT-025481.

Appendix

For completeness, the algorithm ECC06 [START_REF] Kropatsch | The eccentricity transform (of a digital shape)[END_REF] used to compute the eccentricity transform for the shapes in our experiments is included (see [START_REF] Ion | Euclidean eccentricity transform by discrete arc paving[END_REF] for an analysis of the speed/error performance). ECC06 (see Algorithm 1) tries to identify points of the geodesic center (minimum ECC) and use those to find eccentric point candidates. Computing D S (c) for a center point c ∈ C(S) is expected to create local maxima where eccentric points lie. In a first phase, the algorithm identifies at least two diameter ends by repeatedly 'jumping' (computing D S (p)) for the point that had the highest value in the previous estimation. In the second phase, the center points c i are estimated as the points with the minimum eccentricity and all local maxima m of D S (c) are marked as eccentric point candidates. For all m, D S (m) is computed and accumulated. When no new local maxima are found (i.e. with D S (m) not previously computed), the algorithm stops.