Jean-Hugues Pruvot
email: jhpruvot@greyc.ensicaen.fr

Luc Brun

Hierarchy Construction Schemes within the Scale Set Framework

Segmentation algorithms based on an energy minimisation framework often depend on a scale parameter which balances a fit to data and a regularising term. Irregular pyramids are defined as a stack of graphs successively reduced. Within this framework, the scale is often defined implicitly as the height in the pyramid. However, each level of an irregular pyramid can not usually be readily associated to the global optimum of an energy or a global criterion on the base level graph. This last drawback is addressed by the scale set framework designed by Guigues. The methods designed by this author allow to build a hierarchy and to design cuts within this hierarchy which globally minimise an energy. This paper studies the influence of the construction scheme of the initial hierarchy on the resulting optimal cuts. We propose one sequential and one parallel method with two variations within both. Our sequential methods provide partitions near an energy lower bound defined in this paper. Parallel methods require less execution times than the sequential method of Guigues even on sequential machines.

Introduction

Despite much efforts and significant progresses in recent years, image segmentation remains a notoriously challenging computer vision problem. It's usually a preliminary step towards image interpretation and plays a major role in many applications.

The use of an energy minimisation scheme within the region based segmentation framework allows to define criteria which should be globally optimised over a partition. Several types of methods such as the Level set [START_REF] Lecellier | Region-based active contours with noise and shape priors[END_REF], the Bayesian [START_REF] Geman | Stochastic relaxation, gibbs distribution, and the bayesian restoration of images[END_REF], the minimum description length [START_REF] Leclerc | Constructing simple stable descriptions for image partitioning[END_REF] and the minimal cut [START_REF] Boykov | An experimental comparison of min-cut/max-flow algorithms for energy minimization in vision[END_REF] frameworks are based on this approach. Within these frameworks the energy of a partition P is usually defined as E λ (P) = D(P) + λC(P) where D and C denote respectively the fit to data and the regularising term. The energy E λ (P) corresponds to the Lagrangian of the constraint problem: minimise D(P) subject to C(P) ≤ ǫ. Where ǫ is a function of λ. Under large assumptions, minimising E λ (P) is also equivalent to the dual problem: minimise C(P) subject to D(P) ≤ ǫ ′ , where ǫ ′ is also a function of λ. Therefore λ may be interpreted as the amount of freedom allowed to minimise D (D(P) ≤ ǫ ′ ) while keeping C as low as possible. Since ǫ ′ is a growing function of λ, as λ is growing, the constraint on D is more and more relaxed while the importance of the term C is getting more and more important. This parameter λ may thus be interpreted as a scale parameter which represents the relative weighting between the two energy terms.

In many approaches the parameter λ is fixed experimentally and a minimisation algorithm determines for a value of λ a locally optimal partition from the set P of all the possible partitions on image I. A sequence of λ may also be defined a priori in order to compute the optimal partition on each sampled value of λ [START_REF] Guigues | Scale-sets image analysis[END_REF].

The scale set framework proposed by Guigues [START_REF] Guigues | Scale-sets image analysis[END_REF] is based on a different approach. Instead of performing the minimisation scheme on the whole set P of possible partitions of an image I, Guigues proposes to restrict the search on a hierarchy H. The advantages of this approach are twofold: firstly as shown by Guigues the globally optimal partition on H may be found efficiently while the search on the whole set P of partitions only provides local minima. Secondly, Guigues shown that if the energy satisfies some basic properties, the whole set of solutions on H when λ describes R+ corresponds to a sequence of increasing cuts within the hierarchy H hereby providing a contiguous representation of the solutions for the parameter λ. A method to build the hierarchy H has been proposed by Guigues. Since the research space used by Guigues is restricted to the initial hierarchy H the construction scheme of this hierarchy is of crucial importance for the optimal partitions within H built in the second step.

This paper explores different heuristics to build the initial hierarchy. These heuristics represent different compromises between the energy of the final partitions and the execution times. We first present in Section 2 the scale set framework. The different heuristics are then presented in Section 3. These heuristics are evaluated and compared to the method of Guigues in Section 4.

The Scale Set framework

Given an image I and two partitions P and Q on I, we will say that P is finer than Q (or Q is coarser then P) iff Q may be deduced from P by merging operations. This relationship is denoted by P Q. Let us now consider a theoretic segmentation algorithm P λ parametrised by λ. We will say that P is an unbiased multi-scale segmentation algorithm iff for any couple (λ 1 , λ 2 ) such that λ 1 ≤ λ 2 , and any image I, P λ 1 (I) P λ 2 (I). If P λ is an unbiased multi-scale segmentation algorithm, P λ (I) increases according to λ and the set H = λ∈R+ P λ (I) defines a hierarchy as an union of nested partitions. Note that the set P of partitions on I being finite, H must be also finite.

Unbiased multi-scale segmentation algorithms follow a well known causal principal: increasing the scale of observation should not create new information. In other words any phenomenon observed at one scale should be caused by objects defined at finer scales. In our framework, increasing the scale should not create new contours.

The family of energies considered by Guigues corresponds to the set of Affine Separable Energies (ASE) which can be written for any partition P of I in n regions {R 1 , . . . , R n } as:

E(P) = D(P) + λC(P) = n i=1 D(R i ) + λ n i=1 C(R i ) = n i=1 D(R i ) + λC(R i )
Let us consider a hierarchy H and the sequence (C * λ (H)) λ∈R+ of optimal cuts within H. The approach of Guigues is based on the following result: If E λ (P) is an ASE and if C λ (P) is decreasing within P:

∀(P, Q) ∈ P P⊳Q ⇒ C(P) > C(Q) then the sequence (C * λ (H))
λ∈R+ is an unbiased multi-scale segmentation. The union of all (C * λ (H)) λ∈R+ defines thus a new hierarchy within H. The tree corresponding to the hierarchical structure of λ∈R+ C * λ (H) may be deduced from H by merging with their fathers all the nodes which do not belong to any optimal cuts. Note that an equivalent result may be obtained if no condition is imposed to C but if D is increasing according to λ.

The restriction by Guigues of the research space to a hierarchy may thus be justified by the fact that the set of partitions produced by any unbiased multi-scale segmentation algorithm describes a hierarchy. Conversely, given a hierarchy H, if the energy E λ is an ASE with a decreasing term C the sequence of optimal cuts of H according to E λ : (C * λ (H)) λ∈R+ is an unbiased multi-scale segmentation algorithm. Given a partition P ∈ P, the decrease of C may be equivalently expressed as a sub-additivity relationship:

∀(R, R ′ ) ∈ P | R is adjacent to R ′ C(R ∪ R ′ ) < C(R) + C(R ′ ) (1)
Note that the sub-additivity of the regularising term C in common is many applications. For example, if C is proportional to some quantity summed up along contours, C is sub-additive due to the removal of the common boundaries between the two merged regions. Moreover, the term C may be interpreted within the Minimum Description Length framework [START_REF] Leclerc | Constructing simple stable descriptions for image partitioning[END_REF] as the amount of information required to encode a partition. Therefore, one can expect C to decrease when the partition gets coarser.

Given a hierarchy H, the sequence of optimal cuts C * λ (H) within H has to be computed. Let us consider one region R at the second level of the hierarchy (computed from the base) and its set of sons S 1 , . . . , S n . Let us additionally consider the tree H(R) rooted at R within H (Fig. 1(a)). Since R is a level 2 node, the hierarchy H(R) allows only two cuts: one encoding the partition P 1 made of the sons of R whose energy is equal to E λ (P 1 ) = n i=1 D(S i )+λ n i=1 C(S i ) and one encoding the partition P 2 reduced to the single region R. The energy of P 2 is equal to E λ (P 2 ) = D(R) + λC(R). Due to the sub additivity of C we have n i=1 C(R i ) > C(R). Therefore, using the linear expression of E λ (P 1 ) and

E λ (P 2 ) in λ, if n i=1 D(S i ) < D(R) the line E λ (P 1 ) = n i=1 D(S i ) + λ n i=1 C(R i ) is below the line E λ (P 2 ) = D(R) + λC(R) until a value λ + (
R) of λ for which the two lines cross(Fig. 1

(b)). If n i=1 D(S i ) ≥ D(R), E λ (P 2
) is always greater or equal to E λ (P 1 ) in which case we set λ + (R) to 0. Therefore, in both cases the partition P 1 is associated to a lower energy than P 2 for λ = 0 until λ = λ + (R). Above this value the partition P 2 is associated to the lowest energy. In terms of optimal cuts, P 1 corresponds to the optimal cut of H(R) until λ + (R) and P 2 is the optimal cut above this value(Fig. 1(c)). The value λ + (R) is called the scale of appearance of the region R.

Guigues shown that the above process may be generalised to the whole tree. Each node of H is then valuated by a scale of appearance. Some of the nodes of H may get a greater scale of appearance than their father. Such nodes do not belong to any optimal cut and are removed from H during a cleaning step which merges them with their fathers. Each node R of the resulting hierarchy belongs to an optimal cut from λ = λ + (R) until the scale of appearance of its father λ + (F (R)), where F (R) denotes the father of R in H. The value λ + (R) may be set for each node of the tree using a bottom-up process. The optimal cut C * λ (H) for a given value of λ may then be determined using a top-down process which selects in each branch of the tree the first node with a scale of appearance lower than λ. The set of selected nodes constitutes a cut of H which is optimal by construction according to E λ . The function E λ (C * λ (H)) corresponds to a concave piecewise linear function whose each linear interval corresponds to the energy of an optimal cut within H (Fig. 1(d)).

S 1 . . . S n R (a) H(R) E T E λ (P 1 ) n i=1 D(S i ) D(R) E λ (P 2 ) = D(R) + λC(R) λ + (R) λ E λ (b) E λ (P 1 ), E λ (P 2 ) E T λ + (R) λ E λ E λ (C * λ ) (c) E λ (C * λ (H(R))) E T ¨2 2 2 2 2 2 2 A E λ λ (d) E λ (C * λ (H))
Given a hierarchy H and the function E λ (C * λ (H)) encoding the energy of the sequence of optimal cuts, the optimality of H may be measured as the area under the curve E λ (C * λ (H)) for a given range of scales or as the area of the surface A (Fig. 1(d)) between E λ (C * λ (H)) and the energy of the coarsest cut E λ (P max ). Where P max denote the partition composed of a single region encoding the whole image. We propose in Section 4 an alternative measure of the quality of a hierarchy which allows to reduce the influence of the initial image.

Guigues proposed to build a hierarchy H by using an initial partition P 0 and a strategy called the scale climbing. This strategy merges at each step the two adjacent regions R and R ′ such that:

λ + (R ∪ R ′ ) = D(R ∪ R ′ ) -D(R) -D(R ′ ) C(R) + C(R ′ ) -C(R ∪ R ′ ) = min (R 1 ,R 2 )∈P 2 ,R 1 ∼R 2 D(R 1 ∪ R 2 ) -D(R 1 ) -D(R 2 ) C(R 1 ) + C(R 2 ) -C(R 1 ∪ R 2 )
(2) where P denotes the current partition and R 1 ∼ R 2 indicates that R 1 and R 2 are adjacent in P.

This process merges thus at each step the two regions whose union would appear at the lowest scale. Such a construction scheme is coherent with the further processes applied on the hierarchy. However, there is no evidence that the resulting hierarchy may be optimal according to any of the previously mentioned criteria. We indeed show in the next section that other construction schemes of a hierarchy may lead to lower energies.

Construction of the initial hierarchy

Many energies have been designed in order to encode different types of homogeneity criteria (piecewise constant [START_REF] Leclerc | Constructing simple stable descriptions for image partitioning[END_REF][START_REF] Mumford | Optimal approximation by piecewise smooth functions and associated variational problems[END_REF], linear or Polynomial [START_REF] Leclerc | Constructing simple stable descriptions for image partitioning[END_REF] variations,. . . ). This paper being devoted to the construction schemes of the hierarchy, we restrict our topic to the piecewise constant model described by Leclerc [3] and Mumford and Shah [START_REF] Mumford | Optimal approximation by piecewise smooth functions and associated variational problems[END_REF]. The energy of this model may be written as:

E λ (P) = D(P) + λC(P) = n i=1 SE(R i ) + λ|δ(R i )| (3)
where

P = {R 1 , . . . , R n } represents the partition of the image, SE(R i ) = p∈R c p -µ R 2
is the squared error of region R i and |δ(R i )| is the total length of its boundaries. Within the Minimum Description Length framework, S E(R i ) may be understood as the amount of information required to encode the deviation of the data against the model, while |δ(R i )| is proportional to the amount of information required to encode the shape of the model. Within the statistical framework, the squared error may also be understood as the log of the probability that the region satisfies the model (i.e. is constant) using a Gaussian assumption while |δ(R i )| is a regularising term.

Our approach follows the scale climbing strategy proposed by Guigues (equation 2). Given a set W of regions within a partition P we thus consider the scale of appearance of the region R defined as the union of the regions in W. The heuristics below use this basic approach but differ on the sets W which are considered and on the ordering of the merge operations.

Sequential Merging

Given a current partition P, let us consider for each region R of P, its set V(R) defined as {R} union its set of neighbours and the set P * (V(R)) of all possible subsets of V(R) including R. Each subset W ∈ P * (V(R)) encodes a possible merging of the region R with at least one of its neighbour. Let us denote by R W = R ′ ∈W R ′ the region formed by the union of the regions in W. Note that the region R W is connected since R belongs to W and all the regions of W are adjacent to R. Let us additionally consider the two partitions of R W : P R W = {R W } and P W = W. The energies associated to these partitions are respectively equal to E λ (P R W ) = D(R W ) + λC(R W ) and:

E λ (P W ) = D(W) + λC(W) = R ′ ∈W D(R ′ ) + λ R ′ ∈W C(R ′ )
where D(W) and C(W) denote respectively the fit to data and the regularising terms of the partition P W .

Since C is sub additive (equation 1) we have C(W) > C(R W ). The energy E λ (P W ) is thus lower than E λ (P R W ) until a value λ + (R W ) called the scale of appearance of R W (Section 2). Using the scale climbing principle, our sequential merging algorithm computes for each region R of the partition the minimal scale of appearance of a region R W :

λ + min (R) = arg min W∈P * (V(R)) D(R W ) -D(W) C(W) -C(R W )
the set W ∈ P * (V(R)) which realises the min is denoted W min (R).

Given the quantities λ + min (R) and W min (R), our sequential algorithm iterates the following steps:

1. Let P denotes the current partition initialised with an initial partition P 0 , 2. For each region R of P compute λ + min (R) and W min (R) 3. Compute R min = arg min R∈P λ + min (R) and merge all the regions of W min (R min ). 4. If more than one region remains go to step 2, 5. Output the final hierarchy H encoding the sequence of merge operations.

This algorithm performs thus one merge operation at each step of the algorithm. Note that all the regions of W min (R min ) are adjacent to R min . Therefore, within the irregular pyramid framework, the merge operation may be encoded by a contraction kernel of depth one composed of a single tree whose root is equal to R min . The computation of λ + min (R) for each region R of the partition requires to traverse P * (V(R)) whose cardinal is equal to 2 |V(R)|-1 . Therefore, if the partition is encoded by a graph G = (V, E), the complexity of each step of our algorithm is bounded by O(|V|2 k ) where |V| denotes the number of vertices (i.e. the number of regions) and k represents the maximal vertices's degree of G. The cardinal of V is decreased by |W min (R min )| -1 at each iteration. Since |W min (R min )| is at least equal to 2, the cardinal of V decreases by at least 1. The computation of λ + min (R) for each region R of the partition may induce important execution times when the degree of the vertices of the graph is important. However, experiments presented in Section 4 show that the cardinal of the subsets W ∈ P * (R) may be bounded without altering significantly the energy of the optimal cuts. Let us finally note that this algorithm includes the scale climbing approach proposed by Guigues. Indeed, the merge operations studied by Guigues (Section 2) correspond to the subsets W ∈ P * (V(R)) with |W| = 2 which are considered by our algorithm.

Parallel Merge algorithm

Our parallel merge algorithm is based on the notion of maximal matching. A set of edges M of a graph G = (V, E) is called a maximal matching if each vertex of G is incident to at most one edge of M and if M is maximal according to this property. Moreover, we would like to design a maximal matching M such that the scale of appearance of the regions produced by the contraction of M is as low as possible. Let us denote by ι(e), the two vertices incident to e. Using the same approach as in Section 3.1 we associate to each edge e of the graph the scale of appearance λ + (ι(e)) (equation 2) of the region R ι(e) defined as the union of the regions encoded by the two vertices incident to e. Following, the same approach as Haxhimusa [START_REF] Haxhimusa | Constructing stochastic pyramids by midesmaximal independent directed edge set[END_REF] we define our maximal matching as a Maximal Independent Set on the set of edges of the graph. The iterative process which builds the maximal independent set selects at each step edges whose scale of appearance is locally minimal. This process may be formulated thanks to two boolean variables p and q attached to each edge such that: This iterative process stops when no change occurs between two iterations. If n denotes the final iteration, the set of edges such that p n e is true defines a maximal matching [START_REF] Haxhimusa | Constructing stochastic pyramids by midesmaximal independent directed edge set[END_REF] M which encodes the set of edges to be contracted. Moreover, the set of selected edges corresponds to local minima according to the scale of appearance λ + (e). Roughly speaking if λ + (e) is understood as a merge score, one edge between two vertices will be marked (p k e = true) at iteration k, if among all the remaining possible merge operations involving these two vertices, the one involving them is the one with the best merge score. Note that the construction of a maximal matching is only the first step of the method of Haxhimusa which completes this maximal matching in order to get a decimation ratio of order 2. The restriction of our method to a maximal matching allows to restrict the merge operations to edges which become locally optimal at a given iteration. We thus favour the energy criterion against the reduction factor. As shown by Bield [START_REF] Biedl | Tight bounds on maximal and maximum matching[END_REF], the reduction factor in terms of edges induced by the use of a maximal matching is a least equal to 2 k-1 2k-1 where k is the maximal vertex's degree of the graph. The edge's decimation ratio may thus be very low for graphs with important vertices's degrees. Nevertheless, experiments performed on 100 natural images of the Berkeley database 1 shown that the mean vertex's decimation ratio between levels on this database is equal to 1.73 which is comparable to the 2.0 decimation ratio obtained by Haxhimusa.

The local minima selected in equation 4 are computed on decreasing sets along the iterations in order to complete the maximal matching. We can thus consider that the detected minima are less and less significants as the iterations progress. We thus propose an alternative solution which consists in contracting at each step only the edges selected at the first iteration (p 1 e = true). These edges correspond to minima computed on the whole neighbourhood of each edge. This method may be understood as a combination of the method proposed by Haxhimusa [START_REF] Haxhimusa | Constructing stochastic pyramids by midesmaximal independent directed edge set[END_REF] and the stochastic decimation process of Jolion [START_REF] Jolion | Data driven decimation of graphs[END_REF] which consists in merging immediately vertices corresponding to local minima.

Experiments

The different heuristics presented in this paper have been evaluated on the Berkeley database. The evaluated heuristics include our parallel merge heuristic based on a maximal matching (MM) and the variation of this method(MM 1 ) which merges at each step the edges selected during the first iteration (Section 3.2). We also evaluated our sequential method (SM) and two variations of this method: the first variation (S M 2 ), considers for each region R of the partition the subsets of cardinal 2 of V(R). This method corresponds to the heuristic proposed by Guigues. We also evaluated an intermediate method (S M 5 ) which restricts the cardinal of the subsets of V(R) including R to an upper threshold fixed to five in these experiments. All the experiments have used an initial partition obtained by a Watershed algorithm [START_REF] Brun | Hierarchical watersheds within the combinatorial pyramid framework[END_REF].

Fig. 2 shows 5 optimal cuts obtained for increasing values of λ on the Mushroom and Fisherman images of the Berkeley database 1 . The heuristics used to build the hierarchies are displayed on the first column of Fig. 2. The original images are displayed in Fig. 4(a). 2) using an initial partition with a varying number of regions. (b) mean energies of optimal cuts obtained by our heuristics on the Berkeley database. (c) bounds of the optimal cut's energies.

Fig. 3(a) shows the influence of the number of initial regions on the execution time. These curves have been obtained on the Mushroom image with different initial partitions obtained by varying the smoothing parameter of the gradient within our Watershed algorithm.

Fig. 3(b) allows to compare the performance of each heuristic on the whole Berkeley database. However, a direct comparison of the energies obtained by the different heuristics on different images would be meaningless since the shape of the function E λ (C * λ (H)) depends both of the intrinsic performances of the heuristic used to build H and of the image I on which H has been built. We have thus to normalise the energies E λ (C * λ (H)) produced by the different heuristics before any comparison. Given a hierarchy H, since C * λ (H) is an unbiased multi-scale segmentation (Section 2), the hierarchy H obtained by each of our methods may be associated to a value λ H max above which the optimal partition P max is reduced to a single region encoding the whole image. The energy of P max is defined as: E λ (P max ) = D I + λC I where D I = SE(I) denotes the global image's squared error and C I = |δ(I)| the perimeter of the image. Since the energy of the optimal cuts E λ (C * λ (H)) of a hierarchy H is a piecewise linear concave function of λ, the function E λ (C * λ (H)) is below the energy E λ (P max ) associated to the coarser partition(Fig. 3(c)). Moreover, if P 0 denotes the initial partition, the two points (0, E 0 (P 0 )) and (λ max , E λ max (P max )) belong to the curve. Therefore, E λ (C * λ (H)) being concave, it should be above the line connecting these two points. Finally, the line connecting (0, 0) to (λ max , E λ max (P max )) being below the line joining (0, E 0 (P 0 )) and (λ max , E λ max (P max )) we have for any hierarchy H and any scale λ (Fig. 3(c)):

λ λ max E λ max (P max ) ≤ E λ (C * λ (H)) ≤ E λ (P max )
We obtain from this last inequality and after some calculus the following equation: Therefore, using the normalised energy,

∀λ ∈ R+ x λ ≤ 1+ x λ -1 1 + x λ E I ≤ E λ (C * λ (H)) E λ (P max ) ≤ 1 with x λ = λ λ max and E I = λ max C I D I ( 5
) (a) Original Images (b) D(R) = SE(R) (c) D(R) = S E(R)(1 + f ( Int(R) Ext(R) )
E λ (C * λ (H))
E λ (P max ) and the normalised scale x λ = λ λ max , any curve

E λ (C * λ (H))
E λ (P max ) lies in the upper left part of the unit cube [0, 1] 2 . Note that this result is valid for any hierarchy H and thus any heuristic.

Using our piecewise constant model (equation 3), the energy E λ (P max ) is roughly equal to the squared error of the image for small values of λ and may be interpreted as the global variation of the image. The normalised energy allows thus to reduce the influence of the global variation of the images on the energy and to compare energies computed with a same heuristic but on different images. Note however, that the use of the normalised scale x λ = λ λ max discards the absolute value of λ max . We thus do not take into account the range of scales for which the optimal cut is not reduced to the trivial partition P max . However, the absolute value of λ max varies according to each image and each heuristics. The normalised scale allows thus to remove the influence of the image. Moreover, our experiments shown thus that for each image, our different heuristics obtain close λ max values.

Fig. 3(b) represents for each value of x λ and each heuristic, the mean value of the normalised energy

E λ (C * λ (H))
E λ (P max ) computed on the whole set of images of the Berckley database.

As shown in Fig 3(b) the energy of the optimal cuts obtained by the heuristic MM 1 (--) is lower than the one obtained by the maximal matching heuristic (-• -). This result is confirmed by Fig. 2 (lines MM and MM 1 ) where the heuristic MM removes more details of the mushroom at a given scale. This result is connected to the greater decimation ratio of the MM heuristic. The MM heuristic merges at each step regions with important scale of appearance without considering regions which may appear at further steps. The algorithms MM and MM 1 induce equivalent execution times on a sequential machine. The execution times of the method MM 1 (--) are overlayed by the ones of the method MM (-• -) in Fig. 3(a) due to the vertical scale of this figure.

The subjective quality of the partitions obtained by the heuristics MM 1 and S M 2 (Fig. 2) seems roughly similar. We can notice that the heuristic MM 1 seems to produce slightly coarser partitions at each scale. However, considering Fig. 3(b), the optimal energy obtained by the heuristic S M 2 (--) are lower than the one obtained by MM 1 (--). Note that the heuristic MM 1 produces lower execution times than S M 2 even on a sequential machine(Fig. 3(a)).

As shown by Fig. 3(b) the optimal energies produced by the heuristic S M (-+ -) are always below the one produced by the heuristic S M 2 (--). Note that, the curve (-+ -) is close to the diagonal of the square [0, 1] 2 . This last point indicates that on most of the images of the Berkeley database the hierarchies produced by the S M heuristic provide optimal cuts whose normalised energy is closed from the lower bound of the optimal cut's energies (equation 5). This result is confirmed by Fig. 2 where the heuristic S M preserves more details of the image at each scale. However, the heuristic S M is the one which requires the more important execution times on a sequential machine (Fig. 3(a)).

The heuristic S M 5 may be understood as a compromise between S M 2 and S M. As shown by Fig. 3(b) the optimal energies obtained by the heuristic S M 5 ( ) are close to the one obtain by S M(-+ -) and below the one obtained by S M 2 (--). Moreover, as shown by Fig. 3(a), the execution times required by S M 5 are between the one required by the heuristics S M 2 and S M. Finally, the partitions obtained by the S M 5 heuristic in Fig. 2 are closed from the one obtained by the heuristic S M. Fig. 4 shows results obtained using an other fit to data criterion based on the intuitive notion of contrast. The basic idea of this criterion [START_REF] Felzenszwalb | Image segmentation using local variation[END_REF] states that a region should have a higher contrast with its neighbours (called external contrast) than within its eventual subparts (called internal contrast). Let us denote by G e the mean gradient computed along the contour associated to an edge e. The internal and external contrasts of a region R are then respectively defined as Int(R) = max e∈CC(R) G e and Ext(R) = min e∈E|v∈ι(e) G e . Where CC(R) denotes the set of edges which have been contracted to define R and e ∈ E|v ∈ ι(e) denotes the set of edges incident to v. Our new energy combines the contrast and the squared error criteria as follows:

E λ (P) = n i=1 S E(R i ) 1 + f Int(R i ) Ext(R i ) + λ|δ(R i )| (6)
where f () denotes a sigmoid function.

A contrasted region will thus have a low ratio between its internal and external contrast. Conversely, a poorly contrasted region may have a fit to data term close to twice its squared error. As shown by Fig. 3(b) and (c) this energy favours highly contrasted regions. For example, the cloud merged with the sky in Fig. 3(b) remains in Fig. 3(c). Moreover, experiments not reported here, shown us that the same type of discussion about the advantages and drawbacks of the different heuristics may be conducted on this new energy with the same conclusions.

Conclusion

The Scale Set framework is based on two steps: the determination of a hierarchy according to an energy criterion and the determination of optimal cuts within this hierarchy.

We have presented in this article parallel and sequential heuristics to build such hierarchies. The normalised energy of the optimal cuts, associated with these hierarchy are bounded bellow by the diagonal of the unit square [0, 1] 2 . Our experimental results suggest that our sequential heuristic S M provides hierarchies whose normalised energies are closed from this lower bound. This methods may however require important execution times. We thus propose an alternative heuristic providing lower execution time at the price of generally slightly higher optimal cut's energies. Our parallel methods provide greater energies than the one produced by Guigues's heuristic. However, these methods require less execution times even on sequential machine.

Hierarchies encoding a sequence of optimal cuts are usually composed of a lower number of levels and regions than the initial hierarchies built by our merge heuristics. In the future, we would like to use these hierarchies of optimal cuts in order to match two hierarchies encoding the content of two images sharing a significant part of a same scene.

Fig. 1 .

 1 Fig. 1. (a) a node R of the hierarchy whose sons {S 1 , . . . , S n } correspond to initial regions. (b) the energies of the partitions associated to R and {S 1 , . . . , S n } plotted as functions of λ. (c) the energy of the optimal cuts within H(R) (a). (d) an example of concave piecewise linear function encoding the energy of the optimal cuts within a global hierarchy H.

p 1 e

 1 = λ + (e) = min e ′ ∈Γ(e) {λ + (e ′ )} q1 e = e ′ ∈Γ(e) e ∨ q k e ∧ λ + (e) = min e ′ ∈Γ(e) | q k e ′ {λ + (e ′ )} q k+1 e = e ′ ∈Γ(e) p k+1 e ′ (4) where Γ(e) denotes the neighbourhood of the edge e and is defined as Γ(e) = {e} ∪ {e ′ ∈ E|ι(e) ∩ ι(e ′ ) ∅}.

8 Fig. 2 .

 82 Fig. 2. Partitions of the mushroom and the fisherman images at different scales. Each line of the array corresponds to an heuristic whose acronym is indicated on the first column.

Fig. 3 .

 3 Fig. 3. (a) execution times of the different heuristics on the Mushroom image (Fig. 2) using an initial partition with a varying number of regions. (b) mean energies of optimal cuts obtained by our heuristics on the Berkeley database. (c) bounds of the optimal cut's energies.

Fig. 4 .

 4 Fig. 4. (a) Original images. (b) and (c), partitions of the tower image built with a same heuristic(SM) at a same normalised scale (x λ = .8) but with energies defined using two different fit to data terms. (b) is defined using the squared error D(R) = S E(R) while (c) is defined using the formula defined by equation 6.

available at http://www.eecs.berkeley.edu/Research/Projects/CS/vision/bsds/

Color plates are available at the following url: http://www.greyc.ensicaen.fr/∼jhpruvot/Cut/


Xavier Roynard
email: xavier.roynard@mines-paristech.fr

¸ois Franc

Goulette

Paris-Lille-3D: A Point Cloud Dataset for Urban Scene Segmentation and Classification

This article presents a dataset called Paris-Lille-3D 1 . This dataset is composed of several point clouds of outdoor scenes in Paris and Lille, France, with a total of more than 140 million hand labeled and classified points with more than 50 classes (e.g., the ground, cars and benches). This dataset is large enough and of high enough quality to further research on techniques regarding the automatic classification of urban point clouds. The fields to which that research may be applied are vast, as it provides the ability to increase productivity in regards to the management of urban infrastructures. Moreover, this type of data has the potential to be crucial in the field of autonomous vehicles.

Introduction

LiDAR mobile mapping systems were developed in the early 2000s ([3], [START_REF] Abuhadrous | Digitizing and 3d modeling of urban environments and roads using vehicle-borne laser scanner system[END_REF]). Since then, these systems also called Mobile Laser Scanning (MLS) have improved in regards to robustness. They are now marketed by various companies (such as Riegl and Optech) and are increasingly used by surveyor companies. MLS, unlike fixed LiDAR significantly improves productivity and makes it possible to digitize entire cities in 3D. One example is the dataset TorontoCity [START_REF] Wang | TorontoCity: Seeing the World with a Million Eyes[END_REF] with LiDAR data on the entire city of Toronto (with high-precision maps as ground truth but without ground truth regarding point wise classification).

Key aspect of MLS technology is the ability to process large volumes of product data, including segmentation and classification steps. Segmentation is the ability to group points into homogeneous regions. Classification is the ability to give a semantic class to a group of points or to each point (e.g., car, bench and tree).

The two main fields of LiDAR application that require 3D data semantics are infrastructure's management by cities and the creation of precise maps for autonomous vehicles. Cities are increasingly using MLS scans to inventory their street furniture and to update their Geographic Information Systems (GIS). MLS's production capacity allows cities to update their GISs at higher frequencies and with greater accuracies (measurement standardizations) than other systems.

Another major application of 3D outdoor scene classification is the ability to create maps for the navigation of autonomous vehicles. Segmented and classified point clouds can be used to locate vehicles ( [START_REF] Maddern | Leveraging experience for large-scale lidar localisation in changing cities[END_REF]) and to provide critical information regarding urban scenes (e.g., traffic lights, poles and traffic signs) to assit such vehicles with real-time navigation. Currently, work is underway on the ability of the infrastucture to send up-to-date point clouds to vehicles ( [START_REF] Chu | Index coding of point cloud-based road map data for autonomous driving[END_REF]).

The Paris-Lille-3D dataset can further research on this topic and allow deep learning methods to be developed that involve large volumes of data.

Dataset

The Paris-Lille-3D dataset were acquired with the MLS prototype of the center for robotics of MINES ParisTech. The system is a van equipped with a precise localisation system (GPS Novatel FlexPak 6 + IMU Ixsea LANDINS) and a Velodyne HDL-32e LiDAR mounted at the rear of the truck. We have taken care to keep the origin of the Lidar for each measured point and its timestamp. The dataset comprises three point clouds files: Lille1.ply, Lille2.ply, Paris.ply. Each serve as a training base for tests on the two others. There are more than 140 million manually classified points and 50 classes. Indeed, the separation of classes was furthered to note vehicles parked, stopped on the road, and in motion. Figure 1 shows the different classes that were used in the data. The added difficulty is that some classes contain very few examples which makes the task of learning very challenging. The article [START_REF] Roynard | Paris-Lille-3D: a large and high-quality ground truth urban point cloud dataset for automatic segmentation and classification[END_REF] can be referred to for details on the different files and annotations that were used. For the dataset, great care was taken to ensure quality segmentations and classifications were used. Most other urban datasets (such as [START_REF] Serna | Detection, segmentation and classification of 3d urban objects using mathematical morphology and supervised learning[END_REF] and [START_REF] Vallet | Terramobilita/iqmulus urban point cloud analysis benchmark[END_REF]) used semi-automatic methods to create ground truth, and this led to inaccuracies in regards to semantization that could block learning, especially for deep methods.

Figure 2 shows part of one of the point clouds in the dataset. It is raw data that was obtained using the mobile mapping system. The colors correspond to the intensity of the LiDAR return signal. Figure 3 shows a point cloud with labels corresponding to a manual segmentation of this same part of the point cloud. The numerous objects present in the urban scene are shown. Finally, Figure 4 shows the same point cloud with manual classification that can be used as ground truth.

The dataset closest to the Paris-Lille-3D dataset is the Semantic3D dataset [START_REF] Hackel | SEMANTIC3D.NET: a New Large-Scale Point Cloud Classification Benchmark[END_REF]. It was formed using a set of acquisitions of fixed scanners and was also hand labeled. A benchmark with test data is available online 2 . However, the number of classes available in the dataset is relatively limited as it has only eight classes, unlike the Paris-Lille-3D dataset. Moreover, the data from the Semantic3D dataset 2 http://www.semantic3d.net/ is quite different than the data that can be gathered using mobile mapping systems. MLS systems are better suited to digitizing large areas, such as cities, from a terrestrial point of view, and thus, they should be popularized. More details about shortcomings of other existing benchmarks are presented in [START_REF] Roynard | Paris-Lille-3D: a large and high-quality ground truth urban point cloud dataset for automatic segmentation and classification[END_REF].

We propose two evaluation metrics : one for the object detection task and one for the point wise classification. To evaluate detection of objects, we use the same metric as used in iQmulus/TerraMobilita contest [START_REF] Vallet | Terramobilita/iqmulus urban point cloud analysis benchmark[END_REF]. For the point wise classification, we take classical precision, recall, F1 score and IoU to evaluate performance [START_REF] Roynard | Paris-Lille-3D: a large and high-quality ground truth urban point cloud dataset for automatic segmentation and classification[END_REF].

The quality of the manual classification and the quantity of data in the Paris-Lille-3D dataset will enable deepnetwork learning. The publication [START_REF] Roynard | Classification of Point Cloud Scenes with Multiscale Voxel Deep Network[END_REF] presents a new deepnetwork trained on the Paris-Lille-3D dataset. In addition to detecting numerous objects present in an urban scene, the quality of a good classification method lies in its ability to demonstrate that what works in one part of a city works not only in another part of the same city but also in a different city. This is why two different cities were used in the dataset: Paris and Lille in France. A classification benchmark for the dataset will soon be released based on nonpublished data in Paris, Lille and other cities.

Conclusion

We have presented in this article the Paris-Lille-3D dataset and the interest of having this type of dataset for the scientific community to improve the techniques of automatic classification of urban scenes. Important characteristics of our benchmark are high quality of segmentation and classification (the annotation of each point was done completely manually), the wide variety of object classes and that the dataset is quite large, making it useful for deep learning techniques.

Figure 1 .

 1 Figure 1. The classes used in the dataset.

Figure 2 .

 2 Figure 2. A section of the Paris-Lille-3D point cloud (reflectance ranges from 0 [blue] to 255 [red]).

Figure 3 .

 3 Figure 3. The manually labeled Paris-Lille-3D point cloud (one color was used for each object instance).

Figure 4 .

 4 Figure 4. The classes of the Paris-Lille-3D point cloud (one color was used for each class).

http://caor-mines-paristech.fr/fr/ paris-lille-3d-dataset/