Emna Chebbi
email: emna.chebbi@univ-littoral.fr

Patrick Sondi

Eric Ramat

Lucas Rivoirard

Martine Wahl

Simulation of a Clustering Scheme for Vehicular Ad Hoc Networks Using a DEVS-based Virtual Laboratory Environment

Keywords: Vehicular Networks, Routing Protocol, Intelligent Transport Systems, DEVS, Event-B Vehicular Networks, Routing Protocol, Intelligent Transport Systems, DEVS, Event-B Vehicular Networks, Routing Protocol, Intelligent Transport Systems, DEVS, Event-B

come

Introduction

Wireless communication technologies are one of the key factors for the development of the Internet of the future (connected objects, autonomous vehicles and smart cities). Ad hoc network paradigm, particularly, provides high flexibility by offering the possibility of short-range communications between all connected nodes without any preexisting communication infrastructure. The intrinsically decentralized coordination of an ad hoc network makes its safety, robustness and performance highly dependent from the protocols used. The goal is to ensure the discovery of the nodes in the network, their self-organization and routing, and also the quality of service and security of their exchanges. The connected objects can adopt an opportunistic functioning by communicating via an infrastructure when it is operational or switching to ad hoc mode when the infrastructure fails. The Internet protocols initially de-

Introduction

Wireless communication technologies are one of the key factors for the development of the Internet of the future (connected objects, autonomous vehicles and smart cities). Ad hoc network paradigm, particularly, provides high flexibility by offering the possibility of short-range communications between all connected nodes without any preexisting communication infrastructure. The intrinsically decentralized coordination of an ad hoc network makes its safety, robustness and performance highly dependent from the protocols used. The goal is to ensure the discovery of the nodes in the network, their self-organization and routing, and also the quality of service and security of their exchanges. The connected objects can adopt an opportunistic functioning by communicating via an infrastructure when it is operational or switching to ad hoc mode when the infrastructure fails. The Internet protocols initially de-

Introduction

Wireless communication technologies are one of the key factors for the development of the Internet of the future (connected objects, autonomous vehicles and smart cities). Ad hoc network paradigm, particularly, provides high flexibility by offering the possibility of short-range communications between all connected nodes without any preexisting communication infrastructure. The intrinsically decentralized coordination of an ad hoc network makes its safety, robustness and performance highly dependent from the protocols used. The goal is to ensure the discovery of the nodes in the network, their self-organization and routing, and also the quality of service and security of their exchanges. The connected objects can adopt an opportunistic functioning by communicating via an infrastructure when it is operational or switching to ad hoc mode when the infrastructure fails. The Internet protocols initially de-signed for best-effort applications are confronted with the requirements of traditionally more sensitive applications, such as logistics or tactical units, transport, industrial platforms, e-health, etc. Facing the variety of the requirements that all these applications may present, it is necessary to develop fast and effective evaluation and validation tools, from the earliest steps of the protocol design. The objective of our work is to extend a virtual laboratory based on multi-modeling in order to simulate the communication systems dedicated to transport, and especially, the protocols dedicated to ad hoc communications between the vehicles. The targeted approach consists in three phases from the design to the evaluation of the vehicular ad hoc network protocols:

• Step 1: formal modeling of the components that support formal specification (eg. the routing protocol, etc).

This modeling makes it possible to verify and prove the properties of the related components (for the protocol this may concern convergence, quality of service or security properties) by resorting to formal methods and tools (eg. B and Event-B); • Step 2: integrating the formal models into a larger simulation, by the means of multi-modeling: transforming the formal models (from Event-B) into DEVS models (Discrete Event System Specification), and connecting them to the models of the other components of the transport system (eg. the vehicular ad hoc network); • Step 3: integrating real-world data into DEVS simulation and managing the interactions with other specialized simulators for the different components (eg. MATLAB for propagation models, OPNET or NS3 for communications, SUMO for mobility models, etc). The goal is to achieve realistic evaluations of the entire system.

In this paper, the work focuses on the presentation of the models of a vehicular ad hoc network realized with a DEVS-based virtual laboratory environment (VLE). The rest of this paper is organized as follows. After a brief review of the literature in section 2, the proposed approach is described in section3. The models currently developed are explained in section 4. The results of their validation are presented in section 5.

Related work

Several work propose various approaches to design and evaluate the protocols for ad hoc communications in vehicular networks. This review only illustrates the main trends in these approaches through representative examples.

Ad hoc routing protocols

The numerous ad hoc routing protocols 1 identified in the literature2 can be grouped into a limited number of approaches. The first one concerns reactive protocols (DSR, AODV, etc.). In this approach, the vehicles compute the routes on-demand according to the traffic to send 3 . These protocols do not continually maintain network topology. As a result, they generate less routing traffic. However, since the route computation starts only when the traffic has to be sent, this often leads to higher delays, which is not suitable for real-time applications and high mobility environment such as vehicular ad hoc networks. Another family of protocols identified as proactive (OLSR, etc.) permanently maintains a structure in the network topology. These protocols compute periodically the routes ready to use when traffic arrives for one of the known destinations [START_REF] Shenbagapriya | A survey on proactive routing protocols in manets[END_REF] . Though they are convenient for highly dynamic environments, these protocols are more prone to routing traffic overhead when the network is dense. Thus, they generally implement a clustering scheme in order to reduce the effects of broadcast transmissions, including that related to routing traffic, and avoid congestion. The last approach developed concerns the so called geographical protocols (GRP, etc.). These protocols use information on the geographic location of the vehicles in order to organize the routing of the data [START_REF] Maghsoudlou | A survey on geographic routing protocols for mobile ad hoc networks[END_REF] . However, they assume the existence of a centralized location service which provides to each vehicle the positions of all the other vehicles. A such assumption do not cope with pure ad hoc networks where no infrastructure should be needed. They are more likely to serve in infrastructure-based or mixed vehicular networks.

Routing protocols evaluation

Most of the evaluations performed on ad hoc routing protocols focus on measuring the performance obtained based on a set of metrics in specific scenarios. The issues mainly addressed in these studies are the following:

• In particular conditions of density, mobility and given a communication technology, what performance profile can the protocol guarantee to the various vehicular network applications? • For a set of possible configurations of the vehicular network, does the protocol allow all the time to satisfy a set of requirements according to a set of metrics (delay, bandwidth, loss rate, etc.)? • By varying both network configurations and the application requirements, how both the protocol behavior and the application performance vary?

The results obtained through this type of evaluations generally lead to a qualitative and quantitative appreciation on the performance of the protocol. However, they are not performed to obtain the proofs on the protocol intrinsic properties. One of the studied properties concerning ad hoc network protocols is the convergence of the protocol, otherwise its ability to complete its tasks and deliver a stable network structure in a given time. Another property is the robustness which reflects its ability to rebuild its structure and maintain its functionalities despite the changes in the network. A third property studied is its scalability (in terms of vehicle density or traffic load). We are particularly interested by the properties related to quality of service 6 or security [START_REF] Hu | SEAD: Secure efficient distance vector routing for mobile wireless ad hoc networks[END_REF] . These latter are the most crucial for guaranteeing an increased interest for ad hoc communications in future and effective vehicular ad hoc networks.

Contribution of formal approaches

The adoption of wireless technologies in transportation systems has been concerned with formal methods in their earlier stages, since most of these systems impose stringent safety requirements. The European Train Control System (ETCS), to mention few, contains a set of safety applications modeled and verified using B formalism. The testings for verifying that the GSM-R (Global System for Mobile communications in Rail) could meet ETCS requirements were conducted over four years before a final adoption. Many formal approaches have been also proposed for the analysis of the routing protocols for vehicular ad hoc networks. Other authors [START_REF] Gyesik | How to formally model features of network security protocols[END_REF] propose a methodology for verifying the properties related to security in network protocols [START_REF] Diaz | A formal methodology for integral security design and verification of network protocols[END_REF] . Singh et al. [START_REF] Singh | Formal Modeling of Distance Vector Routing Protocol using Event-B[END_REF] present a formal model of AODV with Event-B. Another study [START_REF] Méry | Analysis of DSR protocol in e Event-B. Stabilization[END_REF] concerns a similar model of the DSR protocol in order to prove its properties related to security. Finally, Kamali et al. [START_REF] Kamali | Formal development of wireless sensor-actor networks[END_REF] describe a set of Event-B refinements of a formal model of the OLSR protocol. Two problems persist with formal approaches in network protocol design for intelligent transport systems, and suggest an approach combining formal methods with discrete event simulation:

• The first concerns the joint evaluation of components which are prone to formal modeling with the other components of the ITS which are not adapted to such approach. The tools such as Event-B can only evaluate formal models, while the simulation tools based on multi-modeling allow connecting heterogeneous models and devices in a single simulation-based evaluation process; • The second problem is related to the scalability of the formal tools. Though they provide some extensions that can animate formal models, they do not support large scale animation of a great number of interacting objects, which would be mandatory in order to verify the behavior of systems such as ITS and VANETS.

Therefore, two approaches can apply: introducing DEVS mechanisms into formal tools or coupling formal models using DEVS based tools. Yacoub et al. [START_REF] Yacoub | Using DEv-PROMELA for Modelling and Verification of Software[END_REF] propose to integrate the DEVS mechanisms into the formal modeling tool PROMELA. They partially solve the problems related to scalability, but not those involving interactions with nonformal models since all the components should be modeled with PROMELA. We investigate the other solution.

Projecting formal models into DEVS multi-modeling

The proposed approach consists in developing a methodology in order to integrate formal models in DEVS-based multi-modeling. As shown in figure 1, this methodology operates in two phases:

• The first phase consists in modeling some components using formal tools such as Event-B in order to obtain a set of proofs using the automatic prover, and a set of proof obligations that necessitate an interactive proving process involving an expert (denoted here by RPO for Residual Proof Obligations); • In the second phase, the formal models, the proven properties and the RPOs are transferred into a DEVS multimodeling, which integrates the models of the components that do not fit to formal modeling. Then a simulator generated from the multi-modeling allows evaluating the entire system through a discrete event simulation.

This approach will allow obtaining proven properties through formal tools. It will also address the issue of interacting with the components that are not prone to formal modeling, and that of large-scale scenarios. Moreover, it will allow detecting design errors in formal models when the simulation results are in contradiction with some theorems. Finally, it may help increasing the number of proven properties if the simulation results bring new data that allow solving the RPOs. However, the implementation of this approach raises some issues that need to be addressed :

• Building automatically a DEVS representation of the models, the related proven properties and the RPOs obtained using the formal tool (Event-B). This issue raises itself several implementation problems. • Designing the DEVS multi-modeling in such a way that the simulation results allow verifying the properties that were proven in the formal tool, and also producing data that could be used in an interactive proving process.

DEVS models of VANET components

A first step in the development of the proposed approach is to implement a DEVS-based multi-modeling of all the components used to simulate an ITS (e.g. a VANET). The DEVS modeling that we propose, is based on two variants of the DEVS formalism: P-DEVS [START_REF] Chow | Parallel DEVS: A parallel, hierarchical, modular modeling formalism[END_REF] (Parallel-DEVS) and DS-DEVS [START_REF] Barros | Modeling formalisms for dynamic structure systems[END_REF] (DEVS Dynamic Structure). The first variant manages simultaneous external events and internal transitions by introducing the conflict function. The notion of transient state is also implemented by resorting to zero lifetime events. DS-DEVS and its improvements like DS-DE [START_REF] Barros | Abstract simulators for the DSDE formalism[END_REF] introduce the possibility to modify the graph model during the simulation. For example, it is possible to create and destroy atomic or coupled models, or to create and destroy connections between the models. In our case, all these possibilities, which are not available in the classical DEVS formalism, are fundamental. Indeed, we chose, for the moment, to represent a vehicle as an atomic model whose connections represent the communication channels of the vehicles in the ad hoc network. The second important aspect is the management of the vehicles movements in a 3D continuous space (the road traffic lanes). Several space management options exist: discretization of the space which raises the problem of the discretization step, distributing the space definition within each model (vehicle) or centralizing the definition and the management of the space into a specialized model. We chose the third option. As shown in Figure 2 , the model "space" collaborates with the model "controller" which has a special type: it is an Fig. 2. Ad hoc network: a DEVS model executive from the point of view of DS-DEVS. An executive is an atomic model, unique within a coupled model, that can modify the structure of the coupled model. All these operations are performed by the abstract algorithm of the associated coordinator so that we can guarantee the causality. The couple "space"-"controller" is responsible for: location management of the vehicles, the detection and dynamic creation of potential connections between the vehicles and the appearance and disappearance of the vehicles in the studied section of traffic lane according to their respective trajectory. The "space" model is notified by the "controller" when a vehicle enters or leaves the section. The "space" model calculates the connections based on the changes in the speed and direction sent by each vehicle. For each new connection, the "controller" is notified and it updates the connection graph. Depending on its connections, the "vehicle" model transmits and receives messages that allow it to execute the routing protocol and the clustering method. As discussed in 2.1, proactive routing protocols offer a bigger range of possibilities because they maintain a local topology of the network. So, we chose a variant of the OLSR protocol that implements the clustering algorithm Chain-Branch-Leaf (CBL) [START_REF] Rivoirard | Chain-Branch-Leaf: A clustering scheme for vehicular networks using only V2V communications[END_REF] . As described in Figure 3, CBL builds a chain formed by particularly stable vehicles called "branch" to which attach vehicles located in their coverage area, called "leaf". Through OLSR routing messages, the vehicles exchange information that allows each one deciding in a completely distributed way if it is a branch or a leaf. Each vehicle which is a leaf designates the branch to which it is attached. As for the MPR in native OLSR, when a broadcast message is sent, it is retransmitted only by the vehicles called "branch" A, B, C. In addition, CBL realizes an additional optimization which makes it possible to indicate that a broadcast message must be flooded only upstream (blue arrow) or downstream (green arrow) in the chain. The "vehicle" component including its routing protocol (OLSR implementing the clustering scheme CBL) is a component that can be modeled using formal tools. All the models were developed using VLE, a DEVS-based multi-modeling tool that allows creating coupled models. Except the "vehicle" model which will be modified according to the projection of its related model realized with a formal tool, the other models should remain unchanged. Therefore, it is necessary to validate their design and verify that they behave correctly according to the corresponding vehicular network when compared to other specialized and well-established simulation tools. We performed this validation by simulation using the following configuration:

Validation of the models

• Network size: the network includes a total of 358 vehicles on the A27 highway in France. The trajectory data have been generated based on real-world traffic data of the A27. • Mobility: each vehicle is associated with one of the available trajectories. The trajectory determines the entry date of the vehicle on the road, its movement realized by successive segments with heterogeneous constant speeds (segments of different lengths and different speeds respectively) and its exit date from the road section. • The OLSR protocol operates according to referenced settings [START_REF] Rivoirard | Chain-Branch-Leaf: A clustering scheme for vehicular networks using only V2V communications[END_REF] . For the moment, only the sending and receiving of HELLO messages required by the CBL scheme are modeled. The communication range is fixed to 500 m.

Figure 4 shows the number of vehicles entering to the simulation, the number of vehicles leaving the simulation and the number of vehicles simultaneously present in the road section. We also see the total number of active connections at each observation step. The highest value is reached when the maximum number of simultaneous vehicles is reached (155), then it decreases and stabilizes around 1300 when the number of vehicles stabilizes itself about 130. The total number of branch nodes is 15% up to 35% of the total number of nodes (Figure 5). Therefore about 70% of the vehicles are leaf nodes and do not retransmit broadcast traffic, which confirms CBL performance results [START_REF] Rivoirard | Chain-Branch-Leaf: A clustering scheme for vehicular networks using only V2V communications[END_REF] . The number of branch nodes per chain indicates that there are 1 to 2 connected subsets in each traffic direction, which confirms that the vehicular network is entirely connected, at least in the same traffic direction. Figure 6 shows that each vehicle has an average of 80 neighbors in the entire VANET, thus 40 in the same traffic direction. Those which are branch nodes are selected by 20% to 40% of their neighbors which attach as leaf nodes. Each of these leaf nodes remains attached to a branch node 90% of the time (Figure 7), and a node remains without any branch less than 10% of the time. This demonstrates that isolated nodes remains only few time out of the vehicular network. All these results

show that these DEVS-based models of a VANET implemented using VLE allow obtaining the same performance for OLSR with the CBL scheme than those obtained with well-established tools such as MATLAB and OPNET [START_REF] Rivoirard | Chain-Branch-Leaf: A clustering scheme for vehicular networks using only V2V communications[END_REF] . In addition to the testings on each component, these additional results contribute to validate the designed models.

Conclusion

This paper presented a new approach in development which consists in projecting a set of formal models and proven properties on these models through formal tools, in a DEVS-based multi-modeling. Our goal is to put in interaction all these formal models with DEVS models of other components not prone to formal modeling in order to perform the evaluation of the global transport system in a single simulation process. This will allow verifying, by simulation, that proven properties on formal models that might not be sufficiently refined, are not contradicted in certain scenarios. It would also allow producing the results of the simulation as data that could be used in an interactive formal proving loop instead of a human expert. This paper presents the DEVS-based models already realized of an ad hoc network using the OLSR routing protocol and the CBL clustering scheme. These models will be used for the development and the proof of concept of the approach we are developing. These models were implemented using VLE, and validated through a scenario of an ad hoc vehicular network built based on actual real-world traffic data on the A27 highway in France. Future work will concern the design of a formal model of OLSR implementing the CBL scheme in Event-B.

4Fig. 1 .

 1 Fig.1. Main steps of a methodology for enhancing the proven properties of ITS components using Event-B and DEVS simulation[START_REF] Chebbi | Enhancing Formal Proofs of Network Protocols for Transport Systems using Discrete Event Simulation[END_REF]

Fig. 3 .

 3 Fig. 3. VANET topology and traffic structuring with the CBL (Chain-Branch-Leaf) scheme

Fig. 4 .Fig. 5 .

 45 Fig. 4. Evolution of the number of vehicles in the simulation

Fig. 6 .Fig. 7 .

 67 Fig.6. The number of one-hop neighbors and leaf nodes per branch in the simulation

Chebbi, Sondi, Ramat, Rivoirard and Wahl / Procedia Computer Science 00 (2016) 000-000

Chebbi, Sondi, Ramat, Rivoirard and Wahl / Procedia Computer Science 00 (2016) 000-000

Chebbi, Sondi, Ramat, Rivoirard and Wahl / Procedia Computer Science 00 (2016) 000-000

Acknowledgments

The authors acknowledge the CPER ELSAT2020, the European Union, the FEDER and Hauts de France Region.


Isabelle Sivignon
email: sivignon@lis.inpg.fr

Rodolphe Breton

Florent Dupont

Eric Andrès

Discrete Analytical Curve Reconstruction without patches

Keywords:

Invertible Euclidean reconstruction methods without patches for 2D and 3D discrete curves are proposed. From a discrete 4-connected curve in 2D, or 6-connected curve in 3D, the proposed algorithms compute a polygonal line which digitization with the standard model is equal to all the pixels or voxels of the curve. The framework of this method is the discrete analytical geometry and parameter spaces are used in order to simplify the algorithms. Moreover, the reconstructed polyline is more compact than classical methods such as the Marching Cubes.

Introduction

The reconstruction of discrete objects is mainly performed in practice with the "Marching Cubes" method [START_REF] Lorensen | Marching cubes : a high resolution 3D surface construction algorithm[END_REF] (and all its follow ups). For a couple of years another approach, based on discrete analytical geometry, is investigated in the discrete geometry community. The aim is to decompose the boundary of a discrete object into discrete analytical polygons and then these polygons into Euclidean polygons. The method has to be invertible, i.e. the discretization of the reconstructed boundary has to be equal to the original discrete object. We do not want any information to be added nor lost. The aim of this new approach is to provide a more compact reconstruction. Several other attempts have already been made in this direction that are not satisfying and usually not invertible (see [START_REF] Coeurjolly | Algorithmique et géométrie discrète pour la caractérisation des courbes et des surfaces[END_REF] for details). Indeed, most the algorithms consist in decomposing a discrete curve (surface) into discrete segments (planes) [START_REF] Debled-Rennesson | Étude et reconnaissance des droites et plans discrets[END_REF] [4] without looking for an exact Euclidean equivalent of the discrete object.

Our method is based on Vittone's recognition algorithm [START_REF] Vittone | Recognition of digital naive planes and polyhedrization[END_REF] for the decomposition of the discrete boundary into discrete line pieces in 2D. The analytical framework is provided by the standard discrete analytical model that defines 2D and 3D discrete polygons [START_REF] Andrès | Discrete linear objects in dimension n: the standard model[END_REF].

In [START_REF] Breton | Towards an invertible euclidean reconstruction of a discrete object[END_REF], we proposed a working algorithm for the 2D case and some hints on how to tackle the 3D surface reconstruction problem. The method we proposed for 2D curves worked basically as follows: a discrete 4-connected boundary is decomposed with Vittone's algorithm [START_REF] Vittone | Recognition of digital naive planes and polyhedrization[END_REF] into discrete line segments in 2D. For each discrete segment a corresponding Euclidean line segment was chosen in the set of all possible solutions (set provided by Vittone's algorithm). The problem we faced then was that these different Euclidean lines could intersect outside the discrete line segments losing the reversibility property. To avoid that, patches (small line segments) where added that force the different Euclidean lines to intersect inside the vertex pixels of the discrete line segments preserving the reversibility property. Nevertheless, extension to dimension 3 is not straightforward. Indeed, the equivalent of patches in 3D is small polygons which are a lot more difficult to define than segments. Although some 3D solutions where proposed in [START_REF] Breton | Towards an invertible euclidean reconstruction of a discrete object[END_REF], they are not very easy to set up.

For that reason we propose in this paper another way of reconstructing a discrete 4-connected curve without patches. The method works basically as follows: we still use Vittone's algorithm. Instead, this time, of performing a recognition of discrete line segments and then replacing them with Euclidean line segments, both are performed at the same time. Indeed, the main idea is to force the first extremity of the segments of the reconstructed polyline to lie inside the discrete curve. Hence, the discrete segment recognition algorithm is constrained by this fixed point. Finally, post-treatments are no more needed to ensure the reversibility of the solution.

A second part of this work shows how to extend this algorithm to the polygonalization of 6-connected 3D discrete lines. We first propose an algorithm to recognize 3D standard line segments and show how the polygonalization process induces many simplifications of this algorithm. Basically, the recognition algorithm relies on three simultaneous 2D recognitions on the projections of the 3D curve. We show that an extra condition has to be settled in order to ensure the correctness of the algorithm. As to the polygonalization process for 3D curves, it works mainly as for the 2D curves. This paper is organized as follows: Section 2 presents some recalls on the tools of discrete geometry we need. In Section 3, the new polygonalization algorithm for 2D standard curves is presented and results are given. Section 4 deals with the 3D curve polygonalization problem: we present a new algorithm to polygonalize in an invertible way any 3D 6-connected curve. Some results are proposed.

Recalls on discrete geometry

In this section, we present some basic notions of discrete geometry that are used in this work. We first introduce the standard digitization model and then present the notion of duality we use to recognize discrete segments.

The standard model

The standard digitization of a Euclidean object consists in all the pixels (resp. voxels) that are cut by the object. The standard lines (resp. planes) can be defined arithmetically [START_REF] Andrès | Defining discrete objects for polygonalization : the standard model[END_REF].

Definition 1 A discrete standard line (resp. plane) of parameters (a, b, µ) (resp. (a, b, c, µ)) is the set of integer points (x, y) (resp. (x, y, z)) verifying:

-ω ≤ ax + by(resp. + cy) + µ < ω where ω = |a|+|b| (resp.+|c| )

2

. A standard line (resp. plane) is a 4-connected line (resp. 6-connected plane).

Remark that the standard digitization requires a so-called orientation convention since the inequalities are not strictly symmetric. This is done by ensuring that a > 0 or if a = 0 that b > 0 or, in the case of planes, if a = b = 0 that c > 0. An illustration of this definition is given in Figure 1. If we denote St(O) the standard digitization of the object O, the following useful properties can be derived from the geometrical definition of this model:

St(O 1 ∩ O 2 ) ⊆ St(O 1 ) ∩ St(O 2 ) and St(O 1 ∪ O 2 ) = St(O 1 ) ∪ St(O 2 ).

Parameter space

In image processing and especially pattern recognition, the dual transformation called Hough Transform [START_REF] Hough | Method and means for recognizing complex patterns[END_REF] is classically used (see [START_REF] Maître | Un panorama de la transformation de Hough[END_REF] for an overview on this method). This transformation is indeed very efficient to recognize parametric shapes in a given image. The general principle of this transformation is to work in a parameter space where each point represents a shape of given parameters. Thus, each selected point in an image is represented by all the shapes that go through it in the parameter space.

In this work, we use a parameter space (0αβ) where a point (α 0 , β 0 ) stands for the Euclidean line of equation α 0 xy + β 0 = 0. Then, each point (x 0 , y 0 ) in the Cartesian space maps to the line αx 0y 0 + β = 0 in the parameter space. Figure 2 illustrates the properties [START_REF] Duda | Use of the Hough transformation to detect lines and curves in pictures[END_REF] which link the Cartesian space with this parameter space. More recently, many geometrical properties on Hough transforms have been studied (see for instance [START_REF] Rosenfeld | A convex polygon is determined by its Hough transform[END_REF] and [START_REF] Milanfar | On the Hough transform of a polygon[END_REF]).

This parameter space P is also very useful in the framework of discrete geometry. Indeed, consider a straight line y = α 0 x+ β 0 (0 ≤ α 0 ≤ 1), the digitization of this line along the Object Boundary Quantization (see [START_REF] Jonas | Digital representation schemes for 3D curves[END_REF]) on an n × n grid is the set of grid point L = {(x, y) ∈ n × n | ⌊α 0 x + β 0 -y⌋ = 0}. Then we can define the domain of a set of discrete points.

Definition 2 Let S be a set of discrete points. The domain of S denoted Dom(S) is the set of parameters (α, β) verifying:

Dom(S) = {(α, β) | ∀(x, y) ∈ S, 0 ≤ αx -y + β < 1}
Then, the domain of a set of pixels is either empty (if the pixels are not collinear) or a convex set since it is defined as the intersection of linear inequalities in the parameter space. Many works have been achieved in order to characterize this domain [START_REF] Dorst | Discrete representation of straight lines[END_REF][START_REF] Mcilroy | A note on discrete representation of lines[END_REF][START_REF] Lindenbaum | On recursive, O(n) partitioning of a digitized curve into digital straight segments[END_REF], and an important result is that if S is a connected discrete segment, then the domain is a (3 or 4)-vertex convex polygon that can only have one of the five shapes illustrated in Figure 3. 3 Reconstruction of a 2D discrete simple curve

Principle

The reconstruction algorithm works on 4-connected curves P = (P 1 , . . . , P n ).

The algorithm starts with a Euclidean point p inside the starting point pixel P 1 of the discrete curve. A recognition direction is chosen once and for all. There is a need to choose a recognition direction especially in case of a closed curve since you can choose two directions from a given starting point. A recognition process determines if a set of pixels is a discrete line segment and if there exists a Euclidean line d that passes through p and through all the pixels of the discrete segment. As long as that holds, a new point P k of the curve, in the recognition direction, is added. Once it does not, a new Euclidean point on the Euclidean line d and inside the last but one pixel is chosen, and the recognition process is started over along the curve until we reach P n .

Starting point and Recognition direction

The method presented in this paper depends on the starting point and on the recognition direction. For each different choice of starting point and recognition direction we get a different reconstruction. We propose a convention which defines a starting point and a recognition direction for any curve. Thus, applying this convention, the reconstruction process always provides the same solution for a given discrete curve. Open and closed simple curves lead to different conventions:

• an open simple curve: the starting point is the curve end-point with lowest abscissa and lowest ordinate in the case of two end-points having the same abscissa. Moreover, the starting point chosen actually induces the recognition direction. • a closed simple curve: the starting point is the curve point with lowest abscissa and lowest ordinate in the case of multiple points having the same abscissa. The clockwise direction is chosen as recognition direction.

This strategy can of course be improved. Different alternative solutions are under investigation. One improvement could be provided by choosing discrete cusps as starting points as in [START_REF] Breton | Towards an invertible euclidean reconstruction of a discrete object[END_REF]. Other cases such as closed non simple curves or open curves with vertices of degree 3 or 4 (a discrete point has at most 4 4-connected neighbours) can be treated pretty much in the same way, as what is presented in section 3.3.2 except that there are more constraints. This is however still under investigation because for multiple neighbouring regions, topological information is helpful.

Recognition process

The recognition process is based on the recognition algorithm proposed by Vittone [START_REF] Vittone | Caractérisation et reconnaissance de droites et de plans en géométrie discrète[END_REF]. This algorithm tells if a set of 8-connected pixels is a discrete line segment and provides the domain of this set of pixels in the parameter space (see Section 2 and [START_REF] Veelaert | Geometric constructions in the digital plane[END_REF] for more details).

The recognition algorithm proposed by Vittone updates the domain in the parameter space as new points of the discrete curve are added. When the parameter space polygon is empty, the given set of pixels in not a discrete straight line segment anymore.

Different modifications to Vittone's algorithm needed however to be done to adapt it to our problem. Firstly, Vittone's algorithm works with 8-connected curves whereas we consider 4-connected curves. A simple shear transform maps the pixels of a 4-connected curve to the pixels of a 8-connected curve: for instance, a pixel (x, y) of a curve lying in the first octant maps to the pixel (x+y, y). Secondly, Vittone's algorithm provides parameter values that do not correspond to what's usual for a discretization. Let's consider the first octant. ) and not ( a b , µ b ) as expected. A translation of the solutions is needed to correct this. Table 1 presents these corrections. These transformations lead to translations of the domain's vertices, and in the following, the so-called "domain" stands for this translated polygon.

Polygonalization algorithm

The polygonalization algorithm is described in Algorithm 1. We detail each step of this algorithm referring to its line numbers.

-→ 0 ≤ ax -(b -a)y + µ < b ↓ µ = µ -b 2 ↓ µ = µ -b 2 After translation 0 ≤ ax-by + µ-b 2 < b b=b-a -→ 0 ≤ ax-(b-a)y + µ-b 2 < b Table 1
Résumé of the different translations from a 8-connected curve to a 4-connected one.

The aim is to reconstruct a Euclidean polyline (p 1 , . . . , p j ) from a 4-connected simple curve given by the ordered set of pixels P = (P 1 , P 2 , . . . , P n ). The first point P 1 is the starting point of the discrete curve. A Euclidean point p 1 is chosen inside the pixel P 1 (line 2). Usually the center of the pixel is chosen unless other constraints exist.

Let us suppose that we have already reconstructed the discrete points (P 1 , . . . , P i ) into the Euclidean polyline (p 1 , . . . , p k ). Then the discrete points (P i , . . . , P n ) remain (line 4). The Euclidean point p k is inside the pixel P i .

A few transformations are needed before the recognition step. Vittone's algorithm works on lines lying in the first octant (for lines axby + µ = 0 where b ≥ a ≥ 0). But the serie of 4-connected discrete points can have an arbitrary direction. The general direction (octant) of the discrete point serie is determined using the Chain code of the set of pixels. Indeed, it is well known that the Chain code of a discrete segment is composed of at most two different directions [START_REF] Freeman | Boundary encoding and processing[END_REF]. From those two directions we derive the octant in which the line lies (line 5). A symmetry is then applied in order to transpose the set of pixels into the first octant. Next the shear transform defined in Section 3.3 is applied and transforms the 4-connected set of pixels (P i , P i+1 , . . . , P n ) into the 8-connected one P ′ i , P ′ i+1 , . . . , P n . The same transformation is also applied to the Euclidean point p k which results in p ′ k (line 6). The Euclidean point

p ′ k = (u k , v k ) corresponds to the Euclidean line d k : u k x + y -v k = 0 in P (line 7). Each point of d k corresponds to a Euclidean line in C that goes through p ′ k .
We can now use Vittone's algorithm. The points P ′ i+j , j > 0 are added one by one in the discrete segment s k until the domain of s k does not intersect the Euclidean line d k anymore (lines 8 to 12). When that happens, if we denote P ′ i+m the last solution pixel, we know that the discrete points (P ′ i , . . . , P ′ i+m ) are recognized as a discrete straight line segment. The intersection between the polygon and the straight line d k in P (denoted I in Algorithm 1) corresponds to all the Euclidean lines in C going through p ′ k and containing (P ′ i , . . . , P ′ i+m ) in their digitization (see Fig. Euclidean line that contains p ′ k in C (Fig. 4). The transformations described in Table 1 are applied back so that we get a Euclidean straight line D that passes through p k and all the pixels P i , . . . , P i+m (line 18). The straight line D intersects the pixel P i+m as an interval. The center of this interval is the point p k+1 in the reconstruction process (line 19). The recognition steps start then all over with the discrete points (P i+m , P i+m+1 , . . . , P n ) and the Euclidean point p k+1 (line 20).

End of the polygonalization process

For a simple curve the recognition process ends simply when we have reached the last discrete point P n of the curve. A Euclidean vertex is set inside the last pixel as last point of the reconstructed polyline. The matter becomes however a little more difficult when we deal with a closed curve. For a closed curve we have P 1 = P n . It is of course reasonable to expect a close Euclidean polyline as result of a closed discrete curve. This means that the last discrete straight line segment that is reconstructed has to have p k+1 = p 1 . This corresponds to an additional condition that has to be setup is Vittone's algorithm. When we reach the last discrete points of the curve (when P n = P 1 is one of the possible points that will be recognized at this step), we have to ensure that both the Euclidean line d k (corresponding to p ′ k ) and the Euclidean line d 1 (corresponding to p ′ 1 ) intersect the polygon in the parameter space (see Figure 5). The recognition algorithm progresses as long as the intersection point of both lines belongs to the parameter polygon. If this goes as far as point P n = P 1 then the reconstruction process is finished. If this condition is not verified it means that the Euclidean line p k p 1 does not contain the discrete points (P i+m-1 , . . . , P n ) and we start over with p k+1 . Figure 6 presents a result obtained with this algorithm.

Algorithm 1 Polygonalization of a 2D 4-connected curve Polygonalization 2Dcurve(ordered set of pixels P )

1: i ← 1, k ← 1 2: Choose a real point p 1 inside P 1 . 3: while (i ≤ n) do 4: s k ← {P i } 5:
Find the current octant with Chain codes.

6:

{P i , . . . P n } become {P ′ i , . . . P ′ n } and p k becomes p ′ k with the shear transform.

7:

I = Dom(s k ) ∩ d k 8:
while (I = ∅ and i ≤ n) do 9:

i ← i+1 10:

s k ← s k ∪ {P ′ i } 11:
Compute the reduction of I according to the constraints of P ′ i .

12:

end while 13:

if (I = ∅) then 14:

s k ← s k -{P ′ i }
and retrieve I before the addition of P ′ i .

15:

i ← i -1 16: end if 17:
Choose a solution π ′ k in I.

18:

Apply the inverse shear transform: π ′ k becomes π k .

19:

Choose a real point p k+1 such that p k+1 belongs to the line of parameters π k and p k+1 ∈ P i .

20:

k ← k + 1 21: end while (a) (b) (c) d 1 d 1 d k d k p ′ 1 p ′ k

3D Discrete curve polygonalization

In the previous section, we have presented a new algorithm to polygonalize any 2D discrete 4-connected curve without post-processing treatments in an invertible way. In this section, we extend this method to 6-connected 3D curves in a three-dimensional cubic grid.

Principle

The principle of the algorithm is exactly the same as in the 2D case. The reconstruction algorithm works on 3D 6-connected curves

V = (V 1 , . . . , V n ).
The algorithm starts with a Euclidean point p inside the starting voxel V 1 of the discrete curve. A recognition direction is chosen, and a 3D digital line recognition process is applied: it determines if a set of voxels is a discrete 3D line segment and if there exists a Euclidean line d that passes through p and through all the voxels of this segment. The voxels are added one by one along the recognition direction as long as this holds. Once this condition fails, a new Euclidean point is chosen inside the last but one voxel and the recognition process starts over along the curve until the voxel V n .

Parameter Spaces and Domains

The parameter space used in the 2D case has been presented in Section 2, and we present here the parameter spaces (extensions of the 2D parameter space) used in this part.

We denote by C the usual Cartesian space and by (x, y, z) a point in this space. In 3D, we can define three parameter spaces P x , P y and P z and the three polymorphic operators C x , C y and C z which link the parameter spaces with the Cartesian space. Their definitions are given in Table 2. The operators C from P to C link one basic geometrical object of P with its representation in C: a point maps to a plane, a plane maps to a point and finally a line maps to a line.

C x : P x → C (α, β, γ) → x + αy + βz + γ = 0 x + αy + βz + γ = 0 → (x, y, z) C y : P y → C (α, β, γ) → αx + y + βz + γ = 0 αx + y + βz + γ = 0 → (x, y, z) C z : P z → C (α, β, γ) → αx + βy + z + γ = 0 αx + βy + γ + z = 0 → (x, y, z)
Table 2 Definitions of the three three-dimensional parameter spaces.

Let us also consider the intersection between those spaces and the planes α = 0 and β = 0. Those intersections define six two-dimensional parameter spaces (Table 3) included in the three three-dimensional parameter spaces.

For instance, the two-dimensional parameter space P xz is equal to P x ∩ (α = 0). Those spaces can be considered either as restrictions of a 3D parameter space or as two-dimensional parameter spaces. Thus, one point of P xz can be considered either as a plane perpendicular to the plane y = 0 in C or as a line in the 2D Cartesian space (0xz). An example which illustrates those definitions is depicted in Figure 7.

The preimage or domain of a set of voxels is defined according to those spaces and the standard digitization scheme St presented in Section 2.

C xz : P xz → C (β, γ) → x + βz + γ = 0 C xy : P xy → C (α, γ) → x + αy + γ = 0 C yz : P yz → C (β, γ) → y + βz + γ = 0 C yx : P yx → C (α, γ) → αx + y + γ = 0 C zy : P zy → C (β, γ) → βy + z + γ = 0 C zx : P zx → C (α, γ) → αx + z + γ = 0
Table 3 Definitions of the six two-dimensional parameter spaces in 3D.

Definition 3 (domain) Consider a set of voxels V . The preimage or domain of V along the x coordinate, denoted Dom x (V ) is:

Dom x (V ) = {(α, β, γ) ∈ P x | V ⊂ St(C x (α, β, γ))}
where St denotes the standard digitization scheme. Dom y (V ) and Dom z (V ) can be defined in the same way according to P y and P z .

The domain of a given set V is the set of Cartesian planes containing V in their standard digitization. This set of planes can be represented in the three three-dimensional parameter spaces up to a projection operation (see [START_REF] Veelaert | Geometric constructions in the digital plane[END_REF]) and in the following, we use the notation Dom(V ) when no particular parameter space needs to be precised or when we deal with the domains of 2D standard lines (see Section 2).

3D Standard line segment recognition

Let us first recall the definition of a standard 3D line [START_REF] Andrès | Discrete linear objects in dimension n: the standard model[END_REF][START_REF] Andrès | Defining discrete objects for polygonalization : the standard model[END_REF].

Definition 4 (3D standard line) Consider a 3D straight line of directional vector (a, b, c), and going through the point (x 0 , y 0 , z 0 ). Then the standard digitization of this line is the set of integer points fulfilling the conditions given by the following double inequalities:

-|a|+|b| 2 ≤ bx -ay + ay 0 -bx 0 < |a|+|b| 2 -|a|+|c| 2 ≤ cx -az + az 0 -cx 0 < |a|+|c| 2 -|b|+|c| 2 ≤ cy -bz + bz 0 -cy 0 < |b|+|c| 2
where the double inequalities are oriented along the standard convention (see Section 2 and [START_REF] Andrès | Discrete linear objects in dimension n: the standard model[END_REF][START_REF] Andrès | Defining discrete objects for polygonalization : the standard model[END_REF]).

In order to design a 3D standard line recognition algorithm, we rewrite this definition using the operators of projection in the Cartesian and discrete spaces and the notion of compatible parameters.

Consider a set of voxels V in Z 3 and an Euclidean object F ⊂ R 3 . We define three projection operators denoted by π for the discrete space and three others denoted by p for the Euclidean space as:

π x (V ) = {(y, z) ∈ Z 2 | ∃(x, y, z) ∈ V } π y (V ) = {(x, z) ∈ Z 2 | ∃(x, y, z) ∈ V } π z (V ) = {(x, y) ∈ Z 2 | ∃(x, y, z) ∈ V } p x (V ) = {(y, z) ∈ R 2 | ∃(x, y, z) ∈ F } p y (V ) = {(x, z) ∈ R 2 | ∃(x, y, z) ∈ F } p z (V ) = {(x, y) ∈ R 2 | ∃(x, y, z) ∈ F }
An illustration of a 3D standard line, a 3D straight line and their projections is given in figure 8. of π c (S). In other words, the parameters of p c (L) have to belong to the 2D domain (see Section 2) of the set of pixels π c (S).

Hence, a 6-connected 3D discrete curve S is a standard 3D line segment if and only if:

(1) the three projections π x (S), π y (S) and π z (S) are 2D standard line segments (2) the parameters of those 2D standard line segments are compatible.

As a consequence, the recognition of a 3D standard line is done in two steps: first check that the projections of this line are 2D standard lines, and next, ensure that there exists a solution 3D straight line.

We use the 2D standard line recognition algorithm presented in the previous section to compute the domains related to the three projections of S. Each point of the 3D curve may induce a modification of two out of the three projections.

Each domain of the projections is a convex polygon (see Section 2) that can be represented in two out of the six two-dimensional parameter spaces presented in part 4.2. For instance, the domain of the projection π z (S) can be represented in P xy and in P yx . Hence, two out of the three projections' domains can be represented in each three-dimensional parameter space. Figure 9 illustrates the embedding of the projections' domains in the parameter spaces for a given set of voxels S.

In order to ensure the compatibility property in the parameter spaces (i.e. during the recognition process), a characterization of the preimage of a 3D Euclidean line is required. This Proposition is illustrated in Figure 10 (a) and (b). We can now give the following result :

Proposition 7 Let S be a 3D discrete segment. Then any line in the parameter space P x that crosses both Dom(π z (S)) and Dom(π y (S)) represents a 3D line L in the Cartesian space such that p z (L) belongs to Dom(π z (S)) and p y (L) belongs to Dom(π y (S)).

The proof of this proposition is straightforward using Proposition 6. Figure 10 (c) and (d) illustrates this result.

Then, according to Proposition 7, any 3D line which crosses the domains of the projections of S in the parameter space maps to a 3D line L C such that two out of its three projections are solutions for the projections of S.

But to ensure the compatibility of the domains, p x (L C ) must also belong to Dom(π x (S)). It is not straightforward to check this condition algorithmically in the parameter spaces we have presented by now. Indeed, such a verification implies to check every couple of points in Dom(π z (S)) and Dom(π y (S)) -this defines a line in the Cartesian space -, then to compute the parameters of the third projection of this line and finally to check whether those parameters belong to Dom(π x (S)) or not. If the answer is "yes" for one couple, then the domains are compatible. Otherwise, S is not a 3D standard line segment. This checking has to be exhaustive to give the exact result.

In the next part, we see how the compatibility checking can actually be simplified in a 3D curve segmentation process.

Polygonalization of a 3D curve

The segmentation process aims at decomposing the 3D discrete curve into 3D discrete segments s k with corresponding straight line segments r k of extremal points r 1 k and r 2 k . As for the 2D case, the main idea is to fix the first extremity of each real segment r k to a given real point v k . In the parameter spaces, a point v is represented by a plane as illustrated in Figure 11. Seeing that, for a given segment s k , all the solution 3D lines have to contain v k , this condition is transposed in the parameter space saying that we only consider the lines included in C -1

x (v k ). Consequently, the domains are no more polygons but simply segments (see Figure 11).

Consider the two domains Dom(π z (S)) and Dom(π y (S)). As said previously, any couple of points in those domains defines a Cartesian line L C . To ensure the compatibility, we have to check that there exist two points m 1 and m 2 in the parameter space, the first one in Dom(π z (S)) and the second one in Dom(π y (S)) such that the third projection of the line L C defined by m 1 and m 2 lies in Dom(π x (S)). The computation of the third projection's parameters from the two points m 1 and m 2 is given by the following property : Then, the projection p x (L C ) maps to the point (-α 2 α 1 , β 1 -β 2 α 1 ) in the parameter space P yz . PROOF. The line L C is defined by the two points m 1 and m 2 , which correspond to the system of equations

     x + α 1 y + β 1 = 0 x + α 2 z + β 2 = 0
. Then we get α 1 y -

α 2 z + β 1 -β 2 = 0.
If α 1 = 0, we can divide by α 1 , and this equation becomes

y -α 2 α 1 z + β 1 -β 2 α 1
= 0, which corresponds to the point (-α 2 α 1 , β 1 -β 2 α 1 ) in P yz . Otherwise, if α 1 = 0 and α 2 = 0, we can divide by α 2 and get the same result in the parameter space P zy . Finally, the case α 1 = 0 and α 2 = 0 is a degenerate case which can be handled algorithmically using infinite domains. 2 Since all the domains are embedded in a plane, all the parameters of the third projection computed from any m 1 and m 2 in P xy and P xz lie on the same line in P yz . Moreover, this plane is not perpendicular to the (0αβ) plane. Then, all the domains, which are segments, can be defined by the α-coordinates of their extremities. For instance, Dom(π z (S)) which is equal to the segment [A z , B z ] in P xy is defined by the interval I z = [α zA , α zB ] where α zA is the α-coordinate of point A z and α zB is the α-coordinate of point B z . In the same way, if Dom(π y (S)) = [A y , B y ], we denote it by the interval I y = [α yA , α yB ] in P xz . Then we have the following property: Proposition 9 Consider a discrete segment S and the domains of its projections in P xy , P xz and P yz denoted I z , I y and I x . Let f be the function defined as: f (α, α ′ ) = -α ′ α . Then the domains are compatible if and only if

I x ∩ [min(f (α 1 , α 2 )), max(f (α 1 , α 2 ))] = ∅ where α 1 = α zA or α 1 = α zB and α 2 = α yA or α 2 = α yB .
The proof of this proposition is straightforward seeing that the function f , which is part of the function defined in Proposition 8, is continuous and monotonic over positive or negative intervals. If I z or I y contains both positive and negative values, then it is simply split into one negative and one positive intervals.

Algorithm 2 describes the global polygonalization algorithm for a 3D curve taking into account all the elements presented so far. This algorithm requires an ordered set of voxels denoted V = {V 1 , . . . V n } as input parameter. This curve may be open or closed. The output of this algorithm is a set of ordered discrete segments s k together with a corresponding real segment r k which standard digitization is s k .

Let us give some precisions and explanations on this algorithm.

In the first line, two variables are set: i is used to track all the voxels of the curve and k counts the number of segments found by the algorithm. On line 2, the first extremity of the first real segment is chosen: this real point v 1 has to belong to the voxel V 1 .

The "while" loop that begins in line 3 and ends in line 29 is the global tracking of the voxels of the curve. Inside this loop, the discrete segment s k is initialized with the current voxel V i (line 4). The extremity of the first real segment r k is set to the current fixed real point v k (line 5). The last initializations are those of the domains: as illustrated in Figure 11, those domains are first lines and then segments after the first voxel addition (lines 6 to 8).

While those domains are not empty, the parameters are compatible and the end of the curve is not reached, the following voxel of the curve is added to the current discrete segment s k (line 12). The projections of the new voxel induce some reductions of the three domains (line 13). Next, the image of the two intervals I z and I y is computed and intersected with I x in order to check the compatibility of the domains (lines 14 to 19). When one condition of the while loop (line 10) is no more fulfilled, if V i is not the last voxel of the curve, then the voxel V i-1 is defined as the last voxel of the segment s k (line 22). One solution l k is chosen for the 3D line and the real point v k+1 (line 26) is both r k last point and r k+1 first point (lines 25 to 27).

Finally, the recognition process starts all over from line 3 with V i-1 as first voxel of the next segment s k+1 and v k+1 as first extremity of the real segment r k+1 .

In the case of closed curves, the same trick as for the 2D case can be used for the end of the recognition process. Indeed, if the discrete curve is closed, Algorithm 2 Polygonalization of a 3D 6-connected curve Polygonalization 3Dcurve(ordered set of voxels V )

1: i ← 1, k ← 1 2: Choose a real point v 1 belonging to V 1 . 3: while (i ≤ n) do 4: s k ← {V i } 5: r 1 k ← v k 6: Dom(π x (s k )) = P yz ∩ C -1 y (v k ) = I x 7: Dom(π y (s k )) = P xz ∩ C -1 x (v k ) = I y 8: Dom(π z (s k )) = P xy ∩ C -1 x (v k ) = I z 9:
compatible ← true i ← i+1

12:

s k ← s k ∪ {V i } 13:
Compute the reductions of the three intervals according to the constraints related to π x (V i ), π y (V i ) and π z (V i ). Compute f (I z , I y ).

16:

if (I x ∩ f (I z , I y ) = ∅) then 17:
compatible ← f alse s k ← s k -{V i } and retrieve I x , I y and I z before the addition of V i .

23:

i ← i -1

24:

end if

25:

Choose a solution l k in the domains.

26:

Choose a real point v k+1 such that v k+1 ∈ l k and v k+1 ∈ V i .

27:

r 2 k = v k+1 28:
k ← k + 1 29: end while the Euclidean curve reconstructed should also be a closed curve. This can be achieved fixing not one extremity but the two extremities of the maybe last real segment and check if its discretization contains the set of last pixels (see section 3.3.2 for details).

Finally, let us study the complexity of this algorithm. All the operations on lines 2, 4 to 9 and 11 to 12 can be done in constant time. The reductions of the intervals on line 13 are done either in constant time with a simple line/line intersection computation, or in logarithmic time if the same trick as in [START_REF] Vittone | Caractérisation et reconnaissance de droites et de plans en géométrie discrète[END_REF] is used. In practice, this second method is often the fastest mean to compute the intersection of two rational lines. Operations on lines 14 to 19 can also be done in constant time since f (I z , I y ) is an interval. In line 25, the choice of the solution is also done in constant time taking the middle point of the solution interval. In line 26, the new real point is computed as the middle of the segment l k ∩V i (V i is the last voxel of the discrete segment) which is done in constant time. At last, the global complexity of this algorithm is O(n log(n)) where n is the length of the curve.

Results

Algorithm 2 has been implemented in C using the multi-precision arithmetic library called GMP [START_REF]GMP library[END_REF]. Using this library enables to work with rational numbers of any precision without rounding any value. The Figures 12 and13 present some results on synthetic curves. Note that the curve of Figure 13 could theoretically be decomposed into 3 discrete segments only. Indeed, the choice of a fixed real point extremity for each Cartesian segment (used to ensure the reversibility) constraints also the discrete segments. In Figure 13 (a), the first real point chosen is the point of coordinates (x 0 , y 0 , z 0 ) which are the coordinates of the first voxel center, whereas in Figure 13 (a), the first real point chosen is the point of coordinates (x 0 + 1 4 , y 0 -1 4 , z 0 + 1 4 ). Even if the number of segments computed is the same for the two cases here (it could be different), note that the extremities of the segments are different.

Conclusions and future work

In this paper we described a framework to find a polygonal curve from a discrete curve with an invertible method in dimension two and three. In 2D, the method proposed in [START_REF] Breton | Towards an invertible euclidean reconstruction of a discrete object[END_REF] has been improved and leads to a new algorithm to vectorize a discrete curve without post-treatment patches. We use the Vittone's algorithm for line recognition and force the vertices of the reconstructed polygon to be inside the curve setting one extremity of each real segment to a chosen point before the recognition process. This ensures the reversibility of the reconstruction.

Then, we extended this algorithm for 3D curves. As in 2D, one extremity of each segment of the reconstructed curve is set before the recognition process. This process is composed of three 2D discrete segments recognitions which are done simultaneously. Moreover, a new constraint, called compatibility constraint, ensures that those three recognitions can lead to a solution 3D line.

The results for the polygonalization of 2D curves are good, but we can notice that the number of real segments obtained is greater than what we got with the algorithm presented in [START_REF] Breton | Towards an invertible euclidean reconstruction of a discrete object[END_REF]. Nevertheless this could be improve by relaxing the constraint of a fix extremity for each real segment while keeping the reversibility property: for instance, we may force each real segment to go through one part of the previous segment which is inside the curve instead of one given point. This should give better results.

In 3D, we presented first results for this problem. We can consider implementing the same improvements as in the 2D case in order to reduce the number of segments found. Moreover, in the context of surface polygonalization, this algorithm may be adapted to polygonalize the border of a discrete plane (coplanar curve). From a discrete surface segmentation into discrete plane pieces, such an algorithm would give an analytical modeling of a discrete surface.

Fig. 1 .

 1 Fig. 1. (a) A standard line of parameters (2, -5, 0); (b) A standard plane of parameters (2, 3, 5, 0). A 3D standard line is a 6-connected discrete line which projections are 2D standard lines. The exact analytical definition together with an example are given in Section 4.

Fig. 2 .

 2 Fig. 2. Illustration of the links between the Cartesian space and the parameter space.

Fig. 3 .

 3 Fig. 3. The 5 possible shapes of the solution set.

  Vittone's algorithm recognizes a discrete straight line 0 ≤ axby + µ < b (with b ≥ a ≥ 0) as a line with parameters ( a b , µ b ). If we have a Euclidean line axby + µ = 0 (with b ≥ a ≥ 0) the 8-connected discretization of this line provides the following discrete analytical line: -b 2 ≤ axby + µ < b 2 . It is a centered solution. Vittone's algorithm would therefore recognize the discretization of the straight line ax-by +µ = 0 as a straight line of parameter ( a b , µ b + 1 2

  settings 0 ≤ axby + µ < b b=b-a

4

 4

Fig. 4 .

 4 Fig. 4. (a) A set of pixels that has been recognized as a discrete line segment. (b) The domain of the pixels and the Euclidean line corresponding to p ′ k in the parameter space. (c) The solution π ′ k in the parameter space is a Euclidean line (in (a)) going through p ′ k .

Fig. 5 .

 5 Fig. 5. Same comments as in Figure 4 except that we look for a Euclidean line joining both points p ′ k and p ′ 1 . In (b) and (c) we can see how this translates in the parameter space.

Fig. 6 .

 6 Fig. 6. Result of the polygonalization algorithm: on the left, a discrete curve and a zoom on the woman hand; on the right, the polygonal curve computed by Algorithm 1 and the same hand detail.

Fig. 7 .

 7 Fig. 7. Illustration of the parameter spaces and operators definitions : the image of the point p 0 ∈ P x (left) by the operator C x is the plane x + β 0 z + γ 0 = 0. This equation is either the equation of a plane in the 3D space (0xyz) (upper right) or the equation of a line is the 2D space (0xz) (lower right).

Definition 5 (Fig. 8 .

 58 Fig. 8. Illustration of a discrete and a straight line and their projections.

Fig. 9 .Proposition 6 Fig. 10 .

 9610 Fig. 9. Representation of the projections' domains of S in the three parameter spaces P x , P y and P z .Proposition 6 Let L C be a straight line in the Cartesian space of direction (a, b, c), with b = 0 and c = 0, and let L P be the image of L C in the parameter space P x (L P = C -1 (L C )). Then the point L P ∩ (α = 0) maps to the line p y (L C ), and the point L P ∩ (β = 0) maps to the line p z (L C ).

Proposition 8

 8 Let L C be a Cartesian line, and let m 1 = (α 1 , β 1 ) (resp. m 2 = (α 2 , β 2 )) be the point in P xy (resp. P xz ) associated to p z (L C ) (resp. p y (L C )).

Fig. 11 .

 11 Fig. 11. Representation of what happens in the parameter space when a point v is fixed: on the right, the domains we consider are now straight line segments.

10 :

 10 while (I x = ∅ and I y = ∅ and I z = ∅ and compatible = true and i ≤ n) do 11:

14 :

 14 if (I x = ∅ and I y = ∅ and I = ∅) then 15:

  if (I x = ∅ or I y = ∅ or I z = ∅ or compatible = f alse) then 22:

Fig. 12 .Fig. 13 .

 1213 Fig. 12. Result of the polygonalization algorithm for 3D 6-connected curves.The algorithm decomposes the curves into discrete and Cartesian segments, one Cartesian segment for one discrete segment. The polygonal line computed is represented by a polygonal dark line, and its digitization is exactly the input discrete curve. The darkest voxels are the extremities of the discrete segments recognized. Remark that each dark voxel contains one extremity of a Cartesian segment.