Eric Marchand 
  
François Chaumette 
  
Fitting 3D Models on Central Catadioptric Images

Increasing the field of view of camera is an important issue practical in robot vision. One solution is to consider catadioptric camera that allows a 360 o field of view. In this paper we propose a 3D model tracking algorithm that allows a fast and reliable tracking of 3D objects within central catadioptric images. The proposed approach relies on the virtual visual servoing approach. All the modeling aspects have been reconsidered to consider the projection model. Results show the method to be robust and efficient.

I. INTRODUCTION

Increasing the field of view of camera is an important issue in robot vision. Omnidirectional cameras have been recently widely studied in the literature. One solution to build such systems is to combine a mirror with a classical camera. Such systems are called catadioptric cameras. When a single projection center is sufficient to describe the projection in the image plane, these systems are referred as central catadioptric cameras.

Theoretical aspects of central catadioptric images are now well known [START_REF] Baker | A theory of single-viewpoint catadioptric image formation[END_REF], [START_REF] Geyer | Catadioptric projective geometry[END_REF], and a lot of works have been done in the structure from motion area [START_REF] Geyer | Mirrors in motion: epipolar geometry and motion estimation[END_REF] and more recently in visual servoing [START_REF] Mezouar | Central catadioptric visual servoing from 3d straight lines[END_REF]. Nevertheless few researches have been done to solve visual tracking problem which is a fundamental issue for the development of visual systems that consider such images. In [START_REF] Mei | Homography-based tracking for central catadioptric cameras[END_REF] an SSD-based tracker allows the tracking of planar structure using an efficient second order minimization method. Closer to our problem is [START_REF] Barreto | Visual servoing/tracking using central catadioptric images[END_REF] where a model-based tracker based on a global non-linear minimization is presented. Although the goal is very similar, the modeling of the cost function as well as the minimization issue that we consider in this paper are different. In this paper we present a 3D model-based tracking, relying on the virtual visual servoing framework [START_REF] Comport | Robust model-based tracking for robot vision[END_REF], [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF]. Assuming small camera motions between two frames, tracking is reduced to a 3D camera localization problem.

When dealing with 3D camera localization or pose computation, most of the approaches proposed in the literature rely on a 3D registration issue. Considering perspective projection full-scale non-linear optimization techniques (e.g., [START_REF] Lowe | Fitting parameterized three-dimensional models to images[END_REF], [START_REF] Drummond | Real-time visual tracking of complex structures[END_REF], [START_REF] Comport | Robust model-based tracking for robot vision[END_REF], [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF]) which consists of minimizing the error between the observation and the forward-projection of the model have proved to be very efficient. In this case, minimization is handled using numerical iterative algorithms such as Newton-Raphson or Levenberg-Marquardt [START_REF] Lowe | Fitting parameterized three-dimensional models to images[END_REF]. These 2D-3D registration techniques rely on the use of a 3D model of the tracked objects. Considering catadioptric camera, Authors are with INRIA, IRISA, Lagadic, F-35000 Rennes, France ; email marchand@irisa.fr similar approaches can be considered [START_REF] Barreto | Visual servoing/tracking using central catadioptric images[END_REF]. Nevertheless, a new projection model has to be considered which implies the definition of new visual features and of new Jacobians to be used in the minimization approach.

In this paper, we first recall the virtual visual servoing approach in Section II. We then present in Section III the visual features used in the minimization process and we determine the analytical form of the corresponding Jacobian. In Section IV we describe the low level image processing method used in our tracker. Finally, experimental results are presented in Section V.

II. OVERVIEW OF THE TRACKING ALGORITHM

The fundamental principle of the virtual visual servoing approach [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF] is to define the pose computation problem as the dual problem of 2D visual servoing [START_REF] Espiau | A new approach to visual servoing in robotics[END_REF], [START_REF] Hutchinson | A tutorial on visual servo control[END_REF]. In visual servoing, the goal is to move a camera in order to observe an object at a given position in the image.

To illustrate the principle, consider the case of an object with various 3D features o P (for instance, o P are the 3D coordinates of some object points in the object frame). A virtual camera is defined, with same intrinsic parameters as the real camera and whose pose in the object frame is defined by the homogeneous matrix c M o . The approach consists of estimating the real pose by minimizing the error ∆ between the observed data s * (the position of a set of features in the image) and the position s of the same features computed by forward-projection according to the current pose:

∆ = k i=1 pr ξ ( c M o , o P i ) -s * i 2 , (1) 
where pr ξ () is the projection model according to the camera intrinsic parameters ξ and where k is the number of considered features. Usually pr ξ () is the perspective projection model [START_REF] Lowe | Fitting parameterized three-dimensional models to images[END_REF], [START_REF] Drummond | Real-time visual tracking of complex structures[END_REF], [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF] but any kind of projection model such as a projection model suitable for catadioptric cameras, can be considered as will be shown in the next section. We will see that all the modeling issue of the tracking have thus to be rewritten. The virtual camera initially at c0 M o is then moved using a visual servoing control law in order to minimise the error ∆. At convergence, the virtual camera reaches the pose c * M o which corresponds the real camera's pose (see Figure 1).

An important assumption is to consider that s * is computed from the image with sufficient precision. When outliers are present in the measures, a robust estimation is required. M-estimators can be considered as a general form of maximum likelihood estimators [START_REF] Huber | Robust Statistics[END_REF]. Many functions have been proposed in the literature which allow uncertain measures to be less likely considered and in some cases completely rejected. In other words, the objective function is modified to reduce the sensitivity to outliers. The robust optimisation problem is then given by

∆ R = k i=1 ρ s i (r) -s * i , (2) 
where ρ(u) is a robust function [START_REF] Huber | Robust Statistics[END_REF] that grows subquadratically and is monotonically non-decreasing with increasing |u|. Iteratively Re-weighted Least Squares (IRLS) is a common method of applying the M-estimator. It converts the M-estimation problem into an equivalent weighted leastsquares problem. Thus, the error to be regulated to 0 is defined as e = D(s(r)s * ),

where D is a diagonal weighting matrix given by D = diag(w 1 , . . . , w k ). Each element of D is a weight which is given to specify the confidence in each feature location. The computation of weights w i is described in [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF].

A simple control law that allows to move a virtual camera can be designed to try to ensure an exponential decoupled decrease of e. It is given by [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF]:

v = -λ( D L s ) + D s(r) -s * , ( 4 
)
where v is the virtual camera velocity (v = (v, ω) where v is the instantaneous linear velocity) and ω is the instantaneous angular camera velocity, L s is called the interaction matrix and links the motion of the feature in the image to the camera velocity (ṡ = L s v) and λ is a gain that tunes the convergence rate. The choice of s (and thus of L s ) is a key point of this algorithm and is now described in details.

III. MODELING ISSUES FOR CATADIOPTRIC CAMERAS

A. Projection models

A unified projection model for central panoramic systems has been proposed by Geyer and Daniilidis [START_REF] Geyer | A unifying theory for central panoramic systems and practical implications[END_REF]. According to this model such cameras can be modeled by a first projection on a sphere with coordinates (0, 0, ξ) followed by a perspective projection on the image plane (see Figure 2). Such a model can be defined using parameter ξ which depends intrinsically of the mirror parameters used in the catadioptric camera.

Assuming that sensor intrinsic parameters are known, point X = (X, Y, Z) projects in the image plane as x = (x, y, 1) such that:

x = f (X) with x = X Z+ξ √ X 2 +Y 2 +Z 2 y = Y Z+ξ √ X 2 +Y 2 +Z 2
(5) 2. Generic projection model as defined in [START_REF] Geyer | A unifying theory for central panoramic systems and practical implications[END_REF] 

ξ Image plane x = (x, y, 1) X = (X, Y, Z) M z x C Fig.

B. Choice of the visual features

In this paper we assume that the 3D model of the object is made of 3D lines. Therefore, following our previous work [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF] we consider as visual features the distance between a point p (extracted from the image, see Section IV) and the projection of this line, (that is an ellipse e(r)) for a given pose . Therefore in our case vector s(r) will defined by:

s(r) =     . . . s i (r) . . .     with s i (r) = d a (x, e(r)) (6) 
where d a (.) defines the algebraic distance between point x and the projection e(r) of the 3D line (see Section III-E).

In this paper we have chosen to represent a 3D line as the intersection of two planes. We will show how this line projects in the image according to the projection model [START_REF] Comport | Robust model-based tracking for robot vision[END_REF] and how to compute the interaction matrix related to the projection of a 3D line, and then to the considered distance. In [START_REF] Barreto | Visual servoing/tracking using central catadioptric images[END_REF] a distance between a point projected on the normal of the contour is considered (as in [START_REF] Drummond | Real-time visual tracking of complex structures[END_REF] in the perspective case). This leads to a different modeling of the visual feature and of its interaction matrix.

Let us note that other representations of 3D lines exist and can also be considered in this framework. For example Andreff et al. [START_REF] Andreff | Visual servoing from lines[END_REF] considered a Pluckerian representation of 3D lines and computed the related interaction matrix for classical perspective cameras. This work has been extended in [START_REF] Mezouar | Central catadioptric visual servoing from 3d straight lines[END_REF] to the case of a catadioptric camera. This Pluckerian representation of 3D line can be easily used as an alternative to the method presented in this paper.

C. Projection of a 3D straight line

A 3D straight line can be represented as the intersection of two 3D planes given by:

P 1 : A 1 X + B 1 Y + C 1 (Z -ξ) = 0 (7) P 2 : A 2 X + B 2 Y + C 2 Z + D 2 = 0
where the 3D parameters are constrained by:

   A 2 1 + B 2 1 + C 2 1 = 1 A 2 2 + B 2 2 + C 2 2 = 1 A 1 A 2 + B 1 B 2 + C 1 C 2 = 0 (8)
so that the two planes with unit normals

N 1 = (A 1 , B 1 , C 1 ) and N 2 = (A 2 , B 2 , C 2 ) are orthogonal.
According to the projection model defined by ( 5) the projection of a straight line is nothing but the perspective projection of a circle defined as the intersection between plane P 1 and the sphere S centered in (0, 0, ξ) with radius 1 (see Figure 2). This 3D circle is characterized by the following system:

X 2 + Y 2 + (Z -ξ) 2 = 1 (9) A 1 X + B 1 Y + C 1 (Z -ξ) = 0 (10) 
Considering [START_REF] Espiau | Closed-loop recursive estimation of 3d features for a mobile vision system[END_REF] we directly obtain from [START_REF] Espiau | A new approach to visual servoing in robotics[END_REF] the equation of the perspective projection of this circle in the image plane (that is an ellipse) defined by:

K 0 x 2 + K 1 y 2 + 2K 2 xy + 2K 3 x + 2K 4 y + K 5 = 0 (11) with      K 0 = 1 + A 2 1 C 2 1 ξ 2 (ξ 2 -1) K 1 = 1 + B 2 1 C 2 1 ξ 2 (ξ 2 -1) K 2 = A1B1 C 2 1 ξ 2 (ξ 2 -1) K 3 = -A1 C1ξ 2 K 4 = -B1 C1ξ 2 K 5 = -1 ξ 2
(12) To obtain a minimal representation, let us define a i = K i /K 5 , i = 0 . . . 4, we have then an ellipse defined by:

Q(x, y) = a 0 x 2 + a 1 y 2 + 2a 2 xy + 2a 3 x + 2a 4 y + 1 (13) with      a 0 = A 2 1 C 2 1 (1 -ξ 2 ) -ξ 2 a 1 = B 2 1 C 2 1 (1 -ξ 2 ) -ξ 2 a 2 = A1B1 C 2 1 (1 -ξ 2 ) a 3 = A1 C1 a 4 = B1 C1 (14) 

D. Interaction matrix

Let us now compute the interaction matrix related to parameters a i . It is possible to rewrite a 0 , a 1 and a 2 using a 3 and a

4    a 0 = a 2 3 (1 -ξ 2 ) -ξ 2 a 1 = a 2 4 (1 -ξ 2 ) -ξ 2 a 2 = a 3 a 4 (1 -ξ 2 ) (15) Leading to    ȧ0 = 2a 3 (1 -ξ 2 ) ȧ3 ȧ1 = 2a 4 (1 -ξ 2 ) ȧ4 ȧ2 = (1 -ξ 2 )(a 3 ȧ4 + a 4 ȧ3 ) (16) 
To obtain the interaction matrix L a i for a i , i = 0 . . . 4, we thus need to determine L a3 L a4 . Considering equation [START_REF] Huber | Robust Statistics[END_REF] we have

ȧ3 = Ȧ1 C1 -A1 C 2 1 Ċ1 = Ȧ1 C1 -a 3 Ċ1 C1 ȧ4 = Ḃ1 C1 -B1 C 2 1 Ċ1 = Ḃ1 C1 -a 4 Ċ1 C1 ( 17 
)
It is then necessary to use Ṅ1 = ( Ȧ1 , Ḃ1 , Ċ1 ). The time variation Ṅ1 is given by [START_REF] Espiau | Closed-loop recursive estimation of 3d features for a mobile vision system[END_REF][1]:

Ṅ1 = - 1 D 2 N 1 N 2 v -N 1 × ω (18) 
leading to:

L A1 = -A 1 A -B 1 A -C 1 A 0 -C 1 B1 L B1 = -A 1 B -B 1 B -C 1 B C 1 0 -A 1 L C1 = -A 1 C -B 1 C -C 1 C -B 1 A 1 0 (19) with A = - A 2 D 2 , B = - B 2 D 2 , C = - C 2 D 2
Plugging ( 19) in (17) then using ( 14), we easily obtain:

L a3 = αa 3 αa 4 α a 3 a 4 -1 -a 2 3 a 4 ( 20 
)
L a4 = βa 3 βa 4 β 1 + a 2 4 -a 3 a 4 -a 3 (21) 
with α = A -Ca 3 and β = B -Ca 4 . Considering ( 20) and ( 21) in ( 16) allows us to deduce the interaction related to a 0 , a 1 and a 2 since we have of course:

L a0 = 2a 3 (1 -ξ 2 )L a3 (22) L a1 = 2a 4 (1 -ξ 2 )L a4 (23) L a2 = (1 -ξ 2 )(a 3 L a4 + a 4 L a3 ) (24)

E. Distance between a point and the projection of a 3D line

There exist various ways to define the distance between a point and an ellipse (ie, here the projection of a 3D line) [START_REF] Rosin | Analyzing error of fit functions for ellipses[END_REF], [START_REF] Fitzgibbon | A buyer's guide to conic fitting[END_REF]. In this paper we decide to use the algebraic distance (see Figure 4). Considering a point (x, y), the algebraic distance between this point and an ellipse of equation ( 13) is given by: 

d a = Q(x, y) (25) 
L da =       x 2 y 2 2xy 2x 2y             L a0 L a1 L a2 L a3 L a4      

IV. LOW LEVEL IMAGE PROCESSING

When dealing with low-level image processing, the contours are sampled at a regular distance. At these sample points a 1 dimensional search is performed to the normal of the contour for corresponding edges. An oriented gradient mask [START_REF] Bouthemy | A maximum likelihood framework for determining moving edges[END_REF] is used to detect the presence of a similar contour. One of the advantages of this method is that it only searches for edges which are aligned in the same direction as the parent contour. An array of 180 masks is generated off-line which is indexed according to the contour angle. This is therefore implemented with convolution efficiency, and leads to real-time performance.

More precisely, the process consists of searching for the corresponding point p t+1 in image I t+1 for each point p t (see Figure 5). A 1D search interval {Q j , j ∈ [-J, J]} is determined in the direction δ of the normal to the contour. For each point P t i in the list L t , and for every integer position Q j , we compute a criterion corresponding to the square root of a log-likelihood ratio ζ j [START_REF] Bouthemy | A maximum likelihood framework for determining moving edges[END_REF]. The latter is nothing but the absolute sum of the convolution values, computed at p t and Q j respectively in images I t and I t+1 , using a pre-determined mask M δ function of the orientation of the contour. Then the new position p t+1 is given by:

Q * j = arg max j∈[-J,J] ζ j with ζ j =| I t ν(p t ) * M δ + I t+1 ν(Qj ) * M δ | (26) ν(.
) is the neighborhood of the considered pixel. In this paper the neighborhood is limited to a 7 × 7 pixel mask. That is a trade-off made between real-time performance and mask size. Likewise there is a trade-off to be made between the search distance and real-time performance while considering the maximum inter-frame movement of the object. Currently, this search distance has been set to 10 pixels. This low level search produces a list of k points which are used to calculate the distances defined in (25) and used in [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF] and [START_REF] Geyer | Mirrors in motion: epipolar geometry and motion estimation[END_REF].

V. EXPERIMENTAL RESULTS

The algorithm was tested on real data acquired using two different catadioptric cameras. Two objects were also considered: a box and a set of two plinths. Real images and sensor calibration used to obtain our experimental results are courtesy of the Lasmea. The whole system has been implemented using the ViSP software [START_REF] Marchand | ViSP for visual servoing: a generic software platform with a wide class of robot control skills[END_REF]. Computation time is 100ms for each frame using a 2.6 Ghz PC.

A. Tracking handheld box

In the first experiments a box has been considered. Although the object is quite simple, the object moves very fast which implies very large inter-frame motion as can be seen on the video (see Section V-C). Let us note that despite faces appearance and disappearance (eg, on Figure 6a and 6d or on Figure 7a and7b), partial occlusions of the object and large contour ambiguity (see Figure 7d) tracking is achieved successfully. The robust minimization based on the M-estimation considered in section II is one of the key of the success of the tracking for these difficult sequences. Without robust estimation (as considered in [START_REF] Barreto | Visual servoing/tracking using central catadioptric images[END_REF], that if D = I in equation ( 4 

B. Tracking plinths

In this experiment the considered object is composed by two orthogonal plinths and is made with 7 segments. In the image sequence, the sensor, initially static is then handled and moves in various directions (with translations and rotations) and is then finally put back down. Let us note that this is a very long sequence with more than 800 images. The goal is to use such a tracker within a visual servoing system such as [START_REF] Mezouar | Central catadioptric visual servoing from 3d straight lines[END_REF]. Indeed robust structure tracking is usually one of the key of the success of such navigation systems. In this experiment, despite large and fast motion of the camera tracking is successfully achieved along all the sequence. Figure 8 algorithm is displayed while red points correspond to the point x (found using the ME) used in equation [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF]. Yellow parts of ellipses correspond to the projection of the 3D model (e(r)) for the final estimated pose. Figure 9 shows height images of the full sequence with the forward projection of the 3D model. 

C. Videos

Videos are available from the demo section in the Lagadic web site http://www.irisa.fr/lagadic.

VI. CONCLUSIONS

We have presented in this paper a 3D model-based algorithm suitable for central catadioptric images. It has proved to be fast and reliable and can therefore be considered in various robotics applications such as visual servoing. Let us note that the parametrization used for the line is an interesting alternative to the one proposed in [START_REF] Mezouar | Central catadioptric visual servoing from 3d straight lines[END_REF] and can also be used for visual servoing purpose. 

Fig. 1 .

 1 Fig. 1. Virtual visual servoing principle

Fig. 3 .

 3 Fig. 3. Algebraic distance between a point and an ellipse

Fig. 4 .

 4 Fig. 4. Algebraic distance between a point and an ellipse It is then possible to compute the interaction matrix L d related to this distance. Considering the time derivative of (25): ḋa = ȧ0 x 2 + ȧ1 y 2 + 2 ȧ2 xy + 2 ȧ3 x + 2 ȧ4 y we obtain immediately:

Fig. 5 .

 5 Fig. 5. Determining points position in the next image using the moving edge algorithm: (a) calculating the normal at sample points, (b) sampling along the normal and searching new similar contour.

Fig. 6 .

 6 Fig. 6. Tracking a moving box: results on a first sequence

Fig. 7 .

 7 Fig. 7. Tracking a moving box: results on a second sequence

Fig. 8 .

 8 Fig. 8. Tracking plinths: red points correspond to the points p extracted by the ME algorithm along the normals to the previous 3D model projection (in blue). Forward projection of the model for a given pose is displayed in yellow.

Fig. 9 .

 9 Fig. 9. Tracking a plinth. The plinth is static while the sensor is moving.

http://www.irisa.fr/lagadic

ACKNOWLEDGMENTS

The authors wish to thank Hicham Hadj-Abdelkader and Maxime Lhuillier from Lasmea who provided respectively the image sequences and the sensor calibration. They also wish to thank Youcef Mezouar and, again, Hicham Hadj-Abdelkader for fruitful discussions.