Caroline Nadeau

Alexandre Krupa

Improving ultrasound intensity-based visual servoing: tracking and positioning tasks with 2D and bi-plane probes

Keywords: Visual servoing, intensity-based control, ultrasound, bi-plane sensor

Real time and non invasive, the ultrasound imaging modality can easily be used in minimally invasive surgery or needle insertion procedures to visualize an organ or a tumor to reach. However the manual stabilization of the ultrasound image while the organ moves with patient breathing or heart beating can be very tricky. In this paper, we present an intensity-based approach to control both in-plane and outof-plane motions of an ultrasound probe held by a robotic arm in order to reach and follow one organ cross section. Two methods are proposed to improve the accuracy of this intensity-based approach, by estimating on-line the 3D image gradient required in the control law and by considering a bi-plane sensor. Robotic experiments are performed with two different ultrasound sensors on a realistic abdominal phantom and validate this visual servoing approach.

I. INTRODUCTION

The visual servoing techniques consist in integrating the visual feedback provided by a vision sensor in order to control the movements of a dynamic system with a closedloop approach [START_REF] Espiau | A new approach to visual servoing in robotics[END_REF]. Historically, a large majority of the visual servoing approaches proposed in the literature consider a monocular camera to control a robotic system. More particularly in the medical field, the first image-guided robotic systems have been introduced in the context of laparoscopic surgery, using as vision sensor an endoscopic camera [START_REF] Uecker | Automated instrument tracking in robotically assisted laparoscopic surgery[END_REF]. Artificial markers fixed on a surgical tool were used as visual features to control the endoscope in order to keep the instrument centered in the image.

In the particular case of ultrasound (US) visual servoing, the image formation principle and the geometry of the vision sensor which gives information only in its image plane are far different from the ones of a camera. Consequently, the challenges in the visual control differ. While light changes and depth estimation are no more an issue with the US modality, the major challenges consist in the processing of the images and in the control of the out-of-plane motions of the sensor. First works in US visual servoing focus on the control of the in-plane motions of the probe. In [START_REF] Abolmaesumi | Image-guided control of a robot for medical ultrasound[END_REF], the three in-plane degrees of freedom (dof) of a robotic system are controlled to maintain the section of an artery centered in the US image during a manual out-of-plane translation of the US probe. In [START_REF] Hong | A motion adaptable needle placement instrument based on tumor specific ultrasonic image segmentation[END_REF], two dof of a needle-insertion robot are controlled by visual servoing to perform a percutaneous cholecystostomy C. Nadeau is with Université de Rennes I, IRISA, INRIA Rennes-Bretagne Atlantique, France. A. Krupa is with INRIA Rennes-Bretagne Atlantique, IRISA, France cnadeau@irisa.fr, alexandre.krupa@inria.fr while compensating involuntary patient motions. The target and the needle are automatically segmented in the US images and their respective poses are used to control the robot.

Some authors have proposed solutions to control the outof-plane motions of the probe by visual servoing. In [START_REF] Novotny | Real-time visual servoing of a robot using three-dimensional ultrasound[END_REF], an actuated surgical tool is guided to a desired target under 3D US imaging. Both instrument and target are segmented and tracked in the 3D image and the pose error is then used to control 4 dof of the instrument. However 3D US devices currently provide small and low-quality volumes at a low acquisition rate which limits their use in real-time robotic applications. In another work [START_REF] Lee | Ultrasound-based visual servoing system for lithotripsy[END_REF], the three translations of a XYZ stage robot equipped with two US probes and a HIFU transducer are controlled to follow a kidney stone while compensating for physiological motions during a lithotripsy procedure. With a 2D US probe, Vitrani et al. [START_REF] Vitrani | Robust ultrasoundbased visual servoing for beating heart intracardiac surgery[END_REF] proposed to use the two image points corresponding to the intersection of a surgical forceps with the probe plane as visual features to control 4 dof of the tool inside a beating heart in an eyeto-hand configuration. In relation with this work, the authors of [START_REF] Sauvee | US image based visual servoing of a surgical instrument through non-linear model predictive control[END_REF] developed a predictive control scheme to keep the forceps visible in the US image. More recently, different approaches have been proposed to control the 6 dof of the probe with six geometric features built from 2D moments extracted from a single US image [START_REF] Mebarki | 2D ultrasound probe complete guidance by visual servoing using image moments[END_REF] or three orthogonal images [START_REF] Nadeau | A multi-plane approach for ultrasound visual servoing : application to a registration task[END_REF]. However the moments computation requires a segmentation step whose efficiency is dependent on the organ shape and which is not robust to organ topology changes.

To overcome the image processing limitation, recent works have proposed to use directly as visual features the intensity value of the image pixels. In [START_REF] Krupa | Real time motion stabilization with B-mode ultrasound using image speckle information and visual servoing[END_REF], an approach based on the speckle correlation observed in successive US images is detailed. However, a step of learning of the decorrelation curves is required to control the out-of-plane motions of the probe and this local approach is only fitted for tracking tasks. In [START_REF] Nadeau | Intensity-based direct visual servoing of an ultrasound probe[END_REF], we proposed an intensity-based approach to control the 6 dof of a 2D probe, where the interaction matrix relies on the estimation of the 3D image gradient. In this previous work, this estimation was made possible once at the initial pose of the probe, which restrained the method to tracking applications. In this paper we present a new way to estimate on-line the 3D image gradient, which allows the computation of the current interaction matrix of the system and the extension of the intensity-based approach to positioning tasks. The medical applications we consider are needle insertion procedures or medical diagnoses where the visual servoing can be used to reach and stabilize the view of an organ or tumor while compensating for physiological motions of the patient. The second contribution of this work is the modelling of the US intensity-based approach for dealing with a bi-plane probe.

The structure of our paper is as follows. We initially detail the intensity-based approach with the on-line estimation of the interaction matrix. Validation and improvements with respect to the use of the initial interaction matrix for positioning tasks are demonstrated in simulation. In Section III, we propose to extend this intensity-based approach to take into account the visual information extracted from both orthogonal images provided by an US bi-plane device. Experiments of tracking and positioning tasks with 2D and bi-plane US sensors are then presented in Section IV. Finally, perspectives and concluding remarks are given in Section V.

II. INTENSITY-BASED CONTROL

An image-based visual servoing control scheme consists in minimizing the error between a current set of visual features s and a desired one s * . In the case of a number p of visual features exceeding the number m of dof of the controlled system, a combination matrix C of size m × p and full rank m is introduced to define the task function to minimize:

e(t) = C (s(t) -s * ) . (1)
In an eye-in-hand system, the instantaneous velocity applied to the imaging sensor v p is computed from the task function and the interaction matrix L s . More particularly, the classical visual servoing control law described in [START_REF] Espiau | A new approach to visual servoing in robotics[END_REF] is based on an exponential decrease of the task function (ė = -λ e, λ > 0) and presents the best behavior with C = L s + , which yields to:

v p = -λ L s + L s -1 L s + (s(t) -s * ) = -λ L s + (s(t) -s * ) (2)

A. Intensity features

With visual servoing techniques, the efficiency of the control law is highly dependent of the choice of the visual features. Typically, geometric features, which can be extracted by various methods based on image similarity measure, contour segmentation or image processing (threshold and morphological operators), are well-adapted to control the different motions of a robotic system. However, the processing of the US images is often more complex than with other imaging modalities due to the noise, called speckle, generated by the propagation of the US waves in the soft tissues. In this paper, we propose to avoid any US image processing step by choosing as visual feature the image itself. In this case, the visual features vector s corresponds to the intensity values of the pixels contained in a region of interest (ROI) of the US image:

s = {I(1, 1), ..., I(u, v), ..., I(M, N)} , (3)
where M and N are respectively the width and the height of the ROI and where I(u, v) represents the intensity of the pixel of coordinates (u, v) in the US image. We detail in the following the computation of the interaction matrix L s that links the variation of these image features to in-plane and out-of-plane motions of the US probe in order to control the six dof of the device.

B. Modelling of the interaction

Given R p (x p , y p , z p ) the frame attached to the US probe where (x p , y p ) defines the image plane and z p corresponds to the elevation axis, the coordinates x = (x, y, z) of an image point in this frame are such as:

  x y z   =   s x (u -u 0 ) s y (v -v 0 ) 0   , (4)
with (s x , s y ) the image pixel size and (u 0 , v 0 ) the pixel coordinates of the image center, and where z = 0 since the considered visual features belong to the US image plane.

The interaction matrix estimation relies on the constancy of the US wave reflection by a physical 3D point. Given such a 3D point at the position x at the time t that moves to the new position x + dx at the time t + dt and considering that the US reflection is converted in an intensity value in a B-mode US image, the US reflection conservation yields to the following intensity conservation equation:

I(x + dx,t + dt) -I(x,t) = 0, (5)
that can be expanded in the form of a first order Taylor series:

∂ I ∂ x dx + ∂ I ∂ y dy + ∂ I ∂ z dz + ∂ I ∂t dt = 0. (6)
Then, the time variation of each pixel intensity I(u, v) can be expressed as a function of the corresponding 3D point motion :

İ(u, v) = -∇I (u,v) ẋ, (7)
with ∇I (u,v) = (∇I x ∇I y ∇I z ) the 3D image gradient associated to the pixel (u, v).

According to the kinematics fundamental relationship, the velocity of the 3D point ẋ = ( ẋ, ẏ, ż) can be expressed as a function of the in-plane and out-of-plane components of the probe motion v p , through the interaction matrix L x :

ẋ = L x v p , L x =   -1 0 0 0 -z y 0 -1 0 z 0 -x 0 0 -1 -y x 0   . ( 8
)
From ( 7) and ( 8), the interaction matrix L I(u,v) of size 1×6 associated to the visual feature I(u, v) is written as:

L I(u,v) = -∇I x -∇I y -∇I z -y∇I z x∇I z y∇I x -x∇I y , (9)
and the complete interaction matrix L s , defined as ṡ = L s v p , is built by stacking the M × N matrices L I(u,v) .

C. On-line gradient estimation

To control the six dof of the US probe, the variation of the visual features is related to the in-plane and out-of-plane motions of the sensor through the interaction matrix. More particularly the out-of-plane component of the 3D image gradient ∇I z is essential to compute this interaction matrix.

In [START_REF] Nadeau | Intensity-based direct visual servoing of an ultrasound probe[END_REF], we used additional parallel images to compute the 3D gradient with a 3D filter applied on this set of images. However with a conventional 2D US probe mounted on a robotic arm, a small back and forth translational motion along the elevation direction was required to acquire these images. In this case, the interaction matrix was estimated once at the desired pose of the probe and used without being updated to perform a tracking task.

However, for positioning tasks the interaction matrix has to be estimated on-line, and we propose here a new approach to estimate the gradient component ∇I z from the past images acquired during the positioning task. This estimation is based on the Taylor expansion of the pixel intensity:

I(x + dx) = I(x) + ∇I x dx + ∇I y dy + ∇I z dz. ( 10
)
This equation is rewritten as:

Y = Φ θ , (11)
with:

   Y = I(x + dx) -(I(x) + ∇I x dx + ∇I y dy) Φ = dz θ = ∇I z , ( 12
)
where θ is the parameter of the system to estimate and where Y and Φ can be measured at each iteration. The pixel motion dx = (dx, dy, dz) due to the probe motion is computed with the robot odometry while the intensity values of this pixel and the in-plane gradient components (∇I x , ∇I y ) can be directly measured in the image.

To take into account the past values of the measures and then increase the accuracy of this estimation, we propose to solve the equation ( 11) by using the recursive least squares (RLS) algorithm. The RLS algorithm consists in computing the estimate of the parameter θ[k] that minimizes the cost function J( θ[k] ), where J( θ[k] ) corresponds to the squared sum of the prediction errors

E [k] = Y [k] -Φ [k] θ[k] ,
weighted by a forgetting factor 0 < β ≤ 1 used to lower the influence of past data:

J( θ[k] ) = k ∑ j=k 0 β k-j Y [k] -Φ [k] θ[k] 2 . ( 13
)
The parameter θ[k] that minimizes this convex function sets to zero its derivative and is recursively defined as follows [START_REF] Kreisselmeier | Stabilized least-squares type adaptive identifiers[END_REF]:

θ[k] = θ[k-1] + δ (e(k))F [k] Φ [k] Y [k] -Φ [k] θ[k-1] , (14)
where the factor F is also recursively defined as:

F [k] -1 = β F [k-1] -1 + δ (e(k))Φ [k] Φ [k] + ε 0 . ( 15
)
The term ε 0 is known as a stabilization term and ensures that the factor

F [k] remains invertible even if the signal Φ [k]
is not sufficiently excited, which corresponds to a pure inplane motion of the probe. Moreover, a dead zone is added to avoid the update of the parameter estimation when the euclidean norm of the prediction error e(k)

= ||E [k] || is low.
The function δ (e(k)) that characterizes the dead zone is chosen as defined in [START_REF] Middleton | Design issues in adaptive control[END_REF]:

δ (e(k)) =    1 - δ 0 e(k) i f e(k) ≥ δ 0 0 else , ( 16
)
where δ 0 is an estimation of the maximum value of the prediction error due to the noise and disturbance.

D. Simulation validation

To validate our approach, we use a software simulator that we have developed to reconstruct and display a dense volume from a set of parallel images. In addition to this display functionality, the simulator allows the control of a virtual 2D probe that generates an US image by cubic interpolation process. An US complete volume of the right kidney of a realistic abdominal phantom is loaded in the simulator. This volume is created from a set of 335 parallel images of size 250 × 250, the voxel size is 0.6 × 0.6 × 0.3mm 3 .

We simulate a positioning task, using the simulation environment to obtain a ground truth of the evolution of the pose error of the US probe. We position the probe on the kidney volume and we consider the corresponding US scan as the desired one. Then the probe is moved away to a new pose where the observed organ section is considered as the initial image. The results of the positioning task are given in Fig. 1. A visual error is defined as the euclidean norm of the error between the current and desired visual features vectors C = (s-s * ) ⊤ (s-s * ). The difference images between the current and the desired US scans are displayed for the initial and final probe poses. The uniform gray color of this difference image after the convergence of the algorithm shows the success of the positioning task since the final image perfectly coincides with the desired one. During the five first iterations of the algorithm an openloop out-of-plane translation is applied to the probe in order to initialize the estimation of the gradient component ∇I z .

Then, during the convergence process, a new estimate of this parameter is computed each time enough motion along the z-axis is observed. The evolution of the value of this estimate for one pixel is shown in Fig. 2. These results can be compared to the ones obtained with no update of the interaction matrix. We consider the same initial and desired poses of the probe and the 3D image gradient is computed from a set of five parallel images acquired at the initial probe pose. The interaction matrix is then estimated and used without being updated during the convergence process. The Fig. 3 shows the behavior of the control law, which fails to converge to the desired pose. III. BI-PLANE APPROACH Among the new generation of US transducers, the XYPA2.5 device commercialized by the Vermon company is a 2.5 MHz bi-plane US probe that provides two orthogonal images at a frame rate of 20Hz (see Fig. 4). With such a frame rate, this sensor can therefore be used for imageguided control of a robotic system, where the additional information provided by the orthogonal image plane can be integrated to improve the performance of the control. The second contribution of this work is then the modelling of the interaction to take into account the visual information extracted from both orthogonal images of a bi-plane probe.

A. Modelling of the interaction

In each image plane US i , a set of image pixels s i is extracted from the ROI and the new features vector is s = (s 0 , s 1 ). The variation of these image features ṡi is related to the instantaneous velocity v US i of the corresponding image plane through the interaction matrix previously described: Both image planes being rigidly attached to the probe, their velocities v US i can be expressed in function of the instantaneous velocity of the probe v p :

ṡ0 = L s 0 v US 0 ṡ1 = L s 1 v US 1 . (17)
∀i ∈ {0, 1} , v US i = US i W p v p (18)
with:

US i W p = i R p i t p × i R p 0 3 i R p , (19)
where i t p and i R p are the translation vector and the rotation matrix of the probe frame F p expressed in the image frame F US i . For the control of the probe, we consider that F US 0 coincides with F p , hence:

   US 0 W p = I 6 US 1 W p = R y 0 3 0 3 R y R y =   0 0 1 0 1 0 -1 0 0   . ( 20
)
The interaction matrix associated to the features vector s extracted from both image planes of the bi-plane probe is then:

L s =   [ -∇I x -∇I y -∇I z -y∇I z x∇I z y∇I x -x∇I y ] ∀(u,v)∈US 0 [ ∇I z -∇I y -∇I x -y∇I x +x∇I y x∇I z -y∇I z ] ∀(u,v)∈US 1   (21)

B. Simulation validation and discussion

For this simulation the same desired pose of the probe as for the 2D approach is considered and the desired features are extracted from both orthogonal images of a virtual biplane probe. A different initial pose is then chosen, farther than in the previous simulation and from which the oneplane algorithm falls in a local minimum. The results of the bi-plane positioning task are given in Fig. 5. In this case, no open-loop out-of-plane translation is required to initialize the estimation of the components ∇I z associated to each pixel of both images. At the first iteration of the algorithm, these out-of-plane gradient components are set to zero, as is shown on the probe pose evolution (see Fig. 5(a)) by the null initial rotational velocity of the probe around its y-axis. The five other dof of the probe are controlled from the in-plane gradient components, generating thus an outof-plane motion of the probe that allows the estimation of the parameters ∇I z with the method proposed in Section II-C. During the positioning process, these parameters are then updated each time enough out-of-plane motion is applied to the corresponding pixel until convergence.

IV. EXPERIMENTAL RESULTS

With the intensity-based approach we propose to provide assistance to a surgeon by reaching and following a desired organ cross section and thus stabilizing the view of the US probe. In this work, we validate both positioning and tracking tasks. First we implement the on-line estimation of the 3D image gradient on a convex 2-5 MHz 2D US probe to perform a tracking task while the phantom undergoes fast motions. Then, in a second part, we consider a bi-plane 2-5 MHz US probe to perform local positioning tasks.

A. Experimental setup

Experiments are performed with an anthropomorphic robotic arm equipped with an US transducer and a force sensor (see Fig. 6). We combine the visual control with a force control that guarantees a constant force of 1N applied to the phantom and servoes the translational motion along the y-axis of the probe frame. As we chose for safety reasons to give priority to the force control over the vision control, the latter can fail to converge to the desired image since the ytranslational velocity component due to the image control is not applied to the probe. To deal with this issue, we propose to separate the dof controlled by the force control from the others. The five velocity components corresponding to the translations along the x and z axis of the probe and to the three rotations are applied physically to the US device while the last component corresponding to the translation along the y axis is virtually applied to a window containing the ROI and included in the US image.

B. Tracking with a 2D probe

We position the 2D US probe on the abdominal phantom and we define a ROI in the observed US B-scan that contributes to the visual control (see Fig. 7(a)). The force/vision servo process is launched after a small automatic back and forth out-of-plane translation used to initialize the estimation of the 3D image gradient. Then we manually apply various large and fast translational and rotational motions to the phantom. The dynamics of these disturbances exceeds the one of the control law in order to create important delays in the tracking and assess the ability of the control to overcome these delays. The tracking results are shown in Fig. 7 and in the attached video.

The interaction matrix, initialized at the desired pose, is updated during the tracking task from the image gradient extracted directly from the current image (for the 2D components ∇I x and ∇I y ) and estimated from the image measures as described in Section II-C (for its out-of-plane component ∇I z ). Despite the important disturbances applied to the phantom, which generate some delays in the tracking, the US probe converges to the desired image.

C. Positioning task with a bi-plane probe

In this second experiment, the bi-plane probe is mounted on the robotic arm and positioned on the abdominal phantom to acquire both desired images. From a different initial pose of the probe, the bi-plane intensity-based approach is then used to reach the desired images. Contrarily to experiments performed with a 2D probe, the image gradient is not initialized at the desired probe pose. The in-plane components are computed directly from the initial images and the out-of-plane component is set to zero at this pose, before being estimated when the probe moves. The results of the positioning task are given in Fig. 8.

V. CONCLUSION AND PERSPECTIVES

This paper presented new methods to improve the accuracy of the US intensity-based approach by computing the current interaction matrix in the control law and by considering a bi-plane device. We have implemented and validated these methods in simulation environment but also on a robotic system with two different US sensors. In further works, we will propose a qualitative and quantitative comparison of the 2D and bi-plane approaches. Finally we plan to propose a full autonomous process of robotic assistance for a medical gesture, where both in-plane and out-of-plane motions of an US probe are controlled to reach and follow one or several desired views of an organ while compensating for physiological motions of a patient. This way, all along the medical gesture, the surgeon could have a display of a stabilized view of a chosen organ or tumor and could automatically navigate from one view to another.

Fig. 1 .

 1 Fig. 1. Positioning task using visual features extracted from a single US image. The visual convergence from the initial US view (a) to the final one (d) is shown by the corresponding difference image with the desired view (respectively (b) and (e)). The visual (c) and pose (f) convergence is displayed during the positioning task.

Fig. 2 .

 2 Fig. 2. Estimation of the ∇I z parameter. (a) Evolution of this value for one pixel. (b) Number of parameters updated throughout the positioning task.

Fig. 3 .

 3 Fig.3. Failure of the positioning task using the initial interaction matrix.

Fig. 4 .

 4 Fig. 4. The XYPA2.5 bi-plane device (Vermon) (left) provides two orthogonal US images. The probe frame F p is aligned with the frame of the image plane US 0 and the plane US 1 corresponds to a rotation of 90 • around the y-axis (right).

Fig. 5 .

 5 Fig. 5. Positioning task with a bi-plane US probe. The convergence is validated in term of pose (a) and visual (b) errors. (c) and (e) are both images of the bi-plane probe at its initial (top) and final (bottom) pose. (d) and (f) are the corresponding difference images with the desired US image. The contribution of the gradient component ∇I x computed from the image plane US 1 in the out-of-plane motion of the US probe allows the on-line estimation of the out-ofplane gradient without initialization step. On the other hand, the extraction of visual features from two orthogonal images increases the convergence domain of the visual servoing.

Fig. 6 .

 6 Fig. 6. Experimental setup. A six dof anthropomorphic robotic system (a) is equipped with a force sensor and a bi-plane or 2D US probe which provides images of a realistic abdominal phantom (b).

Fig. 7 .

 7 Fig. 7. Tracking of an abdominal cross-section (a) with a 2D probe. The tracking delay of the 2D probe is observed through the evolution of the visual error (b) and through the US image (c) and difference image (d) corresponding to the maximum visual error (at t = 33s). Nevertheless, at the end of the disturbance, the probe reaches its desired pose as is shown by the final US image (e) and difference image (f).

ACKNOWLEDGMENT

The authors acknowledge the support of the ANR project USComp of the French National Research Agency.


Antoine Petit

Eric Marchand

Keyvan Kanani

Vision-based Space Autonomous Rendezvous : A Case Study

For a space rendezvous mission, autonomy imposes stringent performance requirements regarding navigation. For the final phase, a vision-based navigation can be a solution. A 3D model-based tracking algorithm has been studied and tested on a mock-up of a telecommunication satellite, using a 6-DOF robotic arm, with satisfactory results, in terms of precision of the pose estimation and computational costs. Quantitative tests in open loop have been carried out to show the robustness of the algorithm to relative inter frame motions chaser/target, orientation variations and illumination conditions. The tracking algorithm has also been successfully implemented in a closed loop chain for visual servoing.

I. INTRODUCTION

A. Context : an autonomous space rendezvous mission

A space rendezvous consists in the approach of a chaser spacecraft towards a target spacecraft, from detection of the target (if necessary) until docking on the target [START_REF] Woffinden | Navigating the road to autonomous orbital rendezvous[END_REF], [START_REF] Fehse | Automated rendezvous and docking of spacecraft[END_REF]. It can be achieved manually or autonomously. A high level of autonomy can be preferred for several reasons (ground control cannot be used because of large communication delays, safety reasons...). This maneuver has various applications, from space station supplying, spacecraft refueling, grasping or repairing. The first rendezvous was done manually in 1966, between Gemini 8, commanded by Neil Armstrong, and Agena target vehicle. More recent examples include the European Space Agency (ESA) Automated Transfer Vehicle (ATV) which autonomously docks onto the International Space Station (ISS), and projects such as the ESA Geostationary Servicing Vehicle (GSV) or the RObotic GEostationary orbit Restorer (ROGER), both tasked to capture, inspect, assist or re-orbit satellites in trouble. The following main phases make up a rendezvous mission (Fig. 1):

1) the launch of the chaser (several days), which reaches the target orbital plan (A, B), in absolute navigation, and detection of the target. 2) the intermediate phase (a few hours until a day), in which the chaser approaches the target orbit, in absolute navigation. It consists in orbital transfers such as Hohmann transfers (B to E). 3) the terminal phase (a few orbit periods) consists in getting closer to the target down to a secure distance, in relative navigation, with successive impulsive orbital maneuvers. Depending on the mission (ATV, Soyouz/Progress, DART), these maneuvers are handled differently [START_REF] Woffinden | Navigating the road to autonomous orbital rendezvous[END_REF]. 4) a forced translation is performed during a few minutes to dock onto the target.

B. Scope of the study

In this work, which extends a former study performed by Astrium known as High-Integrity, Autonomous, Multi-Range Rendezvous and Docking (HARVD) [START_REF] Tingdahl | Harvd image processing algorithms detailed design[END_REF], the target is a telecommunication satellite (Amazonas-2). Since such a satellite can be tracked with good accuracy from the ground, the chaser can be driven on the target orbit and quite close to this target (distance < 50km). For the autonomous terminal phase, the secure distance is set to 1000m [START_REF] Tingdahl | Harvd image processing algorithms detailed design[END_REF]. Along with the forced translation, an alignment phase, or fly-around phase (Fig. 2), with the target docking port is introduced. We only focus on the final approach phase, the most critical one, with precise requirements regarding lateral, longitudinal and angular alignment and computational costs, and whose navigation is proposed to be achieved by computer vision (along with inertial sensors), using a monocular camera mounted on the chaser, for a distance ranging from 1000m to contact. To address this task, some computer vision approaches have already been implemented for proximity operations in a rendezvous mission. On the ATV, a videometer processes a laser beam reflected by retroreflectors installed on the ISS, which is the first operational vision system for spacecraft navigation [START_REF] Pinard | Accurate and autonomous navigation for the ATV[END_REF]. Using cameras, [START_REF] Blais | Integration of a tracking laser range camera with the photogrammetry based space vision system[END_REF], [START_REF] Space | Dart demonstrator to test future autonomous rendezvous technologies in orbit[END_REF], [START_REF] Woffinden | Navigating the road to autonomous orbital rendezvous[END_REF] propose to rely on easy to detect and track markers installed on the target. Other techniques under study in a space context deal with natural features with stereo vision methods to recover the pose, given a 3D model of the object [START_REF] Dionnet | Robust stereo tracking for space robotic applications[END_REF], [START_REF] Jasiobedzki | Pose determination and tracking for autonomous satellite capture[END_REF]. With a monocular camera, a feature-matching computer vision approach had been selected in [START_REF] Tingdahl | Harvd image processing algorithms detailed design[END_REF]. It consists in the extraction of invariant features in the image that are matched to a database built from preliminary learning sessions. It is however too computationally intensive to be used during the whole mission at the required frequency (a few Hz). A critical sensitivity to distance, illumination, relative orientation, and occlusions of the target had been observed [START_REF] Tingdahl | Harvd image processing algorithms detailed design[END_REF]. In this paper we apply a markerless 3D model-based tracking algorithm to the problem of satellite tracking using a monocular camera [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF]. It has already been implemented in an ESA aerospace applications [START_REF] Dionnet | Robust stereo tracking for space robotic applications[END_REF] and in an aeronautics context [START_REF] Coutard | Visual detection and 3d model-based tracking for landing on aircraft carrier[END_REF].

The tracking algorithm is recalled in Section II. In Section III, we present a visual servoing control loop to automatically perform a simulated rendezvous between a robot and a mockup of Amazonas-2. Results regarding tracking performances and robustness in open loop and regarding visual servoing are exposed in Section IV.

II. 3D MODEL-BASED TRACKING

Our problem is restricted to model-based tracking, using a 3D model of the target. This model is made up of lines, which provides a good invariance to pose and illumination changes, and robustness to some image noise or blur. The purpose is to compute the pose (camera, and thus chaser, position and orientation with respect to the target) which provides the best alignment between edges of the projected model and edges extracted in the image [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]. The approach considered in this paper is described in [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF]. Given a new image, the 3D model of the scene or the target is projected in the image according to the estimated previous camera pose r (see Fig. 3). Each projected line L i (r) of the model is then sampled along its orientation in the image with a sample step S. Then from each sample point p i,j a 1D search along the normal of the projected edge is performed, for a range R, to find a corresponding point p i,j in the image. As in [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF], p i,j is of maximum likelihood with regard to p i,j .

In order to compute the new pose, the distances between points p i,j and the projected lines L i are minimized with respect to the following criteria [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF] : where d ⊥ (L i (r), p i,j ) is the distance between a point p i,j and the corresponding line L i (r) projected in the image from a pose r. Here, ρ is a robust estimator, which reduces the sensitivity to outliers. This is a non-linear minimization process with respect to the pose parameters r. The minimization process follows the Virtual Visual Servoing framework [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF] similar to a Gauss-Newton approach. The tuning of the range R, sample step S, and the maximal number of iterations K in minimization process are discussed in Section IV, regarding estimation precision and computational costs.

∆ = i j ρ(d ⊥ (L i (r), p i,j )) (1)

III. TRACKING FOR VISUAL SERVOING

In order to servo the robot, we propose to use the 3Dmodel based tracking algorithm within a 2 1/2 D visual servoing control loop.

Visual servoing consists in using data provided by a vision sensor for controlling the motions of a dynamic system [START_REF] Chaumette | Visual servo control, Part I: Basic approaches[END_REF]. Classically, to achieve a visual servoing task, a set of visual features s has to be selected from the image enabling the control of the desired degrees of freedom. The goal is to minimize the error between the current values of visual features s extracted from the current image and their desired values s * . For this purpose, techniques [START_REF] Chaumette | Visual servo control, Part I: Basic approaches[END_REF] depend on the features s used : they can be 2D points directly extracted from the image, for Image-based Visual Servoing (IBVS) or 3D parameters recovered thanks to image measurements like pose computation for Position-based Visual Servoing (PBVS). Here we apply a hybrid solution, 2 1/2 D visual servoing approach [START_REF] Chaumette | Visual servo control, Part II: Advanced approaches[END_REF], [START_REF] Chaumette | 2 1/2 D visual servoing: a possible solution to improve image-based and position-based visual servoings[END_REF], which avoids the shortcomings of the two basic approaches, by combining features in 2D and 3D, in order to decouple position and rotational movements, with a simpler interaction matrix, and with a better stability than IBVS or PBVS :

s = x y θu z t T (2)
where x and y are the metric coordinates in the image of a point of the object, here the center of the mock-up, θu z is the third coordinate of the θu vector, which represents the rotation the camera has to perform to reach the desired pose, and t is the translation vector the camera has to perform to reach the desired pose, expressed in the desired camera frame. θu z and t have thus to be regulated to 0. We need to minimize the error e = ss * where features s are recovered thanks to model-based tracking. A kinematic controller, which is convenient for most of systems, is then designed to servo the camera. A proportional control scheme is defined, to make the error exponentially decreases, leading to the following control law:

v c = -λ L s + (s -s * ) (3)
with L s + the estimate of the pseudo-inverse of L s , the interaction matrix associated to the visual features. The paper described in [START_REF] Chaumette | Visual servo control, Part I: Basic approaches[END_REF] details how this matrix can be computed. It can here be estimated thanks to the parameters of the pose computed by the model-based tracking algorithm.

IV. RESULTS

A. Experimental conditions

To implement the tracking algorithm on a vision based rendezvous context, Astrium provided a complete 3D-model and a real reduced (1/50) mock-up of Amazonas-2, a telecom satellite built from the Eurostar-3000 platform, and similar to the one used for HARVD experiments [START_REF] Tingdahl | Harvd image processing algorithms detailed design[END_REF]. Amazonas-2 was launched in 2009 for Spanish company Hipsasat to cover the American (especially South America) position. It is located on a Geostationary Orbit.

The provided 3D model is too complex to deal with realtime applications. Thus we have considerably simplified the model, keeping the most significant geometrical features. Besides, depending on the relative size of the mock-up in the image, the relevant information to be considered varies. Indeed, regarding the central module of the satellite, the contours are not very precise and regular, due to the insulating film on most of the surface of this module and on the four circular antennas. Since the film reflects light, the irregularities are enhanced by the sun luminosity and make it hard for the algorithm to identify edges. For long distances (Fig. 6(a)), the information provided by the central module is not useful, nor significant, compared to the information given by the solar panels, with sharper edges. However, for shorter distances (Fig. 6(d)), details of the central module have to be included and the model has to be refined, to perform the approach properly. The solution implemented here is to switch between 3D models in order to use the model with the most relevant information with regard to the distance camera/mock-up.

Using a robot to simulate the rendezvous, with a camera mounted on its end-effector, enables to have regular and quite realistic movements (let us however note that the specific dynamic of spacecraft is not considered in this paper). Besides, the position of the robot is known very precisely (10 -4 m), providing us with ground truth. It enables to exactly replicate test procedures and so to carry out systematic tests to quantify the tracking performances in open loop. The robot we have used (Fig. 4) has 6 degrees of freedom and 640×480 images are processed. Sun illumination can also be simulated by spot lights located around the scene. As the mock-up size is 1/50 of the satellite and as we use for tests a different camera from the one that would be on-board (in this case, a camera with a 5 deg. FoV, along x and y, is supposed to be mounted on the chaser spacecraft), scaling factors need to be derived. Fig. 5 shows the equivalence between the real and experimental contexts with regard to the z-coordinate of the target center in the camera frame, so that the target has the same relative size in the image in both cases. The dynamics along x and y are also affected as relative motions chaser/target in the real context scaled by a factor close to 1/50 induce in the experimental context the same inter frame motions of the target in the image. Typical forced translations for rendezvous scenarios specified in [START_REF] Tingdahl | Harvd image processing algorithms detailed design[END_REF] are performed with velocities ranging from 0.05m/s to 0.10m/s, and with a camera frame rate set to 1Hz in these cases. dz real can be approximated to 0.0024 as we restrict to distances ranging from 56m to 800m, due to the robot joint limits.

B. Results on a nominal scenario

In this scenario, for the first translation, the camera is moving with a constant velocity (0.015m/s) along its optical axis, which targets the center of the mock-up, with a 15Hz frame rate (the equivalent velocity would be 0.40m/s in a space context, so larger than usual translation velocities), and with a distance chaser/target ranging from 1.56m to 0.62m, what corresponds to a 642m-255m range in the space context. Then the robot performs a fly-around phase to realign to the docking port of the mock-up and performs its final translation with a constant velocity (still 0.015m/s), until a distance chaser/target of 0.135m, which is the shortest achievable distance because of joint limits between the robot and the mock-up, and which is equivalent to 56m in the space context. Here, illumination conditions are favorable, leading to distinctive edges on the mock-up. Two tracking configurations, C1 and C2, are implemented, and are specified by the type of 3D model(s) used, their switching distance d s if it is the case, and the tracking key parameters defined in Section II (range R, step S, and the maximum number K of iterations in the control law), (see Tab. I). These values have been empirically determined.

Performances are satisfactory for both cases (see Fig. 6 and Fig. 10), with low positioning errors along x, y and z. In the equivalent space context, these estimation errors along x and y would be less than 0.08% of the distance chaser/target, except for C2 from 82m (90s on Fig. 10) to the center of the target, with around 0.3% of the relative distance, along x. They would always be less than 0.5% along z. Low angular misalignment is also observed (less than 5 deg., except around y for C2 at the far end of the estimation). In this scenario, solutions C1 and C2 appear to have quite similar performances, except that C1 is more stable than C2, due to its higher level of detail and its smaller sample step S. C. Robustness tests 1) Robustness to translation motions: As the target may not be initially centered in the image or because of uncertainties in the chaser and target motions, it has been relevant to test the sensitivity of the algorithm to translation motions of the robot along x and y in the camera frame, with relative orientation alignment. The purpose is to determine 5.5 6.0 6.6 6.8 8.1 6.0 5.8 5.4 /x (pix/frame) Target motion lim.

3.8 5.4 5.8 6.0 7.9 5.9 5.8 6.0 /y (pix/frame) Rel. motion lim. /x 9.4 17.0 33.2 50.1 13.8 16.9 29.6 39.8

(real) (cm/frame) Rel. motion lim. /y 8.7 20.4 39.5 59.0 18.0 22.3 39.5 59.0 (real) (cm/frame) if at different distances chaser/target the tracking is properly performed or not, given the inter frame motion constraints. Tab. II gathers the results, with limits in terms of the motion of the target in the image and of the relative motion chaser/target in the real space context. The tracking is less robust to inter frame translation motions of the target in the image at short distances and for motions along y. C2 appears to be more robust than C1 for distances inferior to d s as its range R is larger. Nevertheless, limits in terms of inter frame relative motions in the space context show that the tracking algorithm is more than suitable for a rendezvous approach with slow and smooth motions, as defined in [START_REF] Tingdahl | Harvd image processing algorithms detailed design[END_REF].

2) Robustness to orientation variations: Tests have shown that with C1 or C2 configuration, tracking can be performed, whatever the initial relative orientation. Some problems may occur when the solar panels tend to appear as single edges in the image (see Fig. 7(a)). This singularity leads to some ambiguities on the contours of these panels when the chaser performs its first translation and fly-around phases (Fig. 7(b), Fig. 7(c)). Fig. 11 shows performances of this approach. We observe greater but still acceptable position errors. In the equivalent space context, the error along z would be below 0.7% of the relative distance for C1, 1% for C2 during the first translation (0s -7s) and fly-around phases (7s -20s) (Fig. 11), then decreasing to remain under 0.36% until the end (at 56m from the target in a real space context), in both configurations. Along x and y, errors would always be below 0.25% of the relative distance for C1 and C2. The orientation error around z reaches 13 deg. for C1 during the fly-around phase, but then gets much lower. More generally, tracking is limited to rotational motions of the chaser in the camera frame implying inter frame target motions in the image similar to the limits defined in Tab. II.

3) Robustness to illumination changes: A quite extreme scenario has been experimented, with brutal changes (Fig. 8). Fig. 12 shows performances of two solutions, different from C1 and C2. Low luminosity triggers the disappearance of the panels. With the lack of this key geometrical information, the four circular antennas become necessary for the tracking to be performed. The information being concentrated in the middle of the image, we can observe some errors for rotations around x and y, and of position for z at the be-

D. Computational costs

In our application, the ability to switch between different models and so to adapt them to the conditions (illumination, distance) is a relevant issue to reduce computational costs. The simple 3D models used in C2 have thus been designed in this sense, what has lead to lower computational costs than for C1 and its single detailed model, with more lines. For parameters K, R and S defined in Section II and which tune the estimation convergence rate, C2 is a trade-off to reduce costs while preserving performances as compared to C1. With an Intel Core 2 Duo processor, the mean execution time per frame for a typical approach is 45ms for C1, 21ms for C2, with a 18ms standard deviation for C1, 4ms for C2.

E. Visual servoing results

We have performed an autonomous rendezvous in a closed loop chain by servoing the robot close to the mock-up, with angular alignment, using the 2D 1/2 visual servoing technique presented in Section III. As exposed in Section I-B, the maneuver has been divided into three phases: a first translation is achieved to drive the target into the center of the image, a fly-around phase to align to the docking port axis of the target, and a final translation until almost docking (0.13m or 55m). The servoing performed on the mock-up has successfully achieved the intended goal, for two illumination conditions (Fig. 9

V. CONCLUSIONS AND FUTURE WORKS

A real-time model-based based approach has been applied on the final approach of an autonomous rendezvous mission, in order to track and automatically approach the target spacecraft. The tests performed in open loop on a mock-up have shown promising performances in terms of pose computation, regarding motions, distance, orientation or illumination constraints. The estimation precision and the robustness of the algorithm would make it suitable for equivalent realistic rendezvous final approaches. This method is also adapted to perform an autonomous rendezvous only based on visual servoing, with successful results.

As the tracking method only deals with simple line models, with a simple hidden-face algorithm, future works would aim at processing a complete, complex surface 3D model, thanks to hardware graphics. Evaluating and propagating tracking uncertainty using statistical methods and improving initialization techniques would also be key issues.

Fig. 1 :

 1 Fig. 1: Phases of a typical rendezvous

Fig. 2 :

 2 Fig.2: Final approach of a rendezvous mission

Fig. 3 :

 3 Fig. 3: Moving Edge principle

Fig. 4 :

 4 Fig. 4: The experimental setup. It includes Amazonas 2 mock-up, camera mounted on the robot end-effector and a light.

  experimental context -z coord.

Fig. 5 :

 5 Fig. 5: Equivalence experimental/real contexts. The scaling factor dz exp

Fig. 6 :

 6 Fig. 6: Tracking for a nominal scenario, under good illumination, with C2. On (a), the chaser starts its first translation, then its flyaround phase on (b). During the second translation, the 3D model and tracking parameters change (c), until d = 0.14m (d).

Fig. 7 :

 7 Fig. 7: Potential problem when tracking with singular orientation, when the solar panels of the target appear as single edges (a).

Fig. 8 :

 8 Fig. 8: Tracking with 3D models switching at illumination changes.

  (a)-(c) and (d)-(f) and provided video).

Fig. 9 :

 9 Fig. 9: Complete approach within a visual servoing experiment under strong and low luminosity.

TABLE I :

 I Tracking configurations

		C1	C2
	Model	Single model	2 models
	ds (m)	No	0.43
	R (pixels)	6	6 (d > ds), 11 (d < ds)
	S (pixels)	5	10 (d > ds), 8 (d < ds)
	K	30	5

TABLE II :

 II Robustness to motions

		C1	C2
	Distance (exp) (m)	0.3 0.5 0.9 1.3 0.3 0.5 0.9 1.3
	Distance (real) (m)	124 206 370 535 124 206 370 535
	Target motion lim.

Fig. 12: Pose errors, for a standard approach with six brutal illumination changes, at times t = 7s, 45s, 52s, 57s, 74s and 88s, using a configuration with two 3D models switching at illumination changes (Multi-model), and a configuration with a single 3D model during the whole approach.


Antoine Petit

Guillaume Caron
email: guillaume.caron@irisa.fr

Hideaki Uchiyama
email: hideaki.uchiyama@irisa.fr

Eric Marchand
email: eric.marchand@irisa.fr

Evaluation of Model based Tracking with TrakMark Dataset

Keywords: I.4.8 [IMAGE PROCESSING AND COMPUTER VISION]: Scene Analysis-Tracking; H.5.1 [ INFORMATION IN-TERFACES AND PRESENTATION (e.g., HCI)]: Multimedia Information Systems-Artificial, augmented, and virtual realities

We benchmark two tracking methods developed in the INRIA Lagadic team with a TrakMark dataset. Since these methods are based on a 3D model based approach, we selected a dataset named "Conference Venue Package 01" that includes a 3D textured model of a scene. For the evaluation. we compute the error of 3D rotation and translation with the ground truth transformation matrix. Through these evaluations, we confirmed that the provided dataset was suitable for quantitative evaluations of 3D visual tracking.

INTRODUCTION

This paper is related to the benchmarking process of tracking algorithm suitable for augmented reality applications. Although realtime tracking and registration methods received much interest in the last few years, this is still a key feature and, unfortunately, one of the bottleneck for the development of augmented reality based applications. Benchmarking such algorithms is then a necessary process to assess the efficiency of the proposed methods.

Whereas benchmarking keypoint detectors and trackers received much interest in the literature [START_REF] Mikolajczyk | A performance evaluation of local descriptors[END_REF], few tentatives have proposed to benchmark and allowed a fair comparison of tracking algorithms. Among the few existing proposals, one can consider the Metaio benchmark [START_REF] Lieberknecht | Benchmarking template-based tracking algorithms[END_REF] and the TrakMark benchmark. The former one, which is very demanding, mainly focuses on the evaluation of 2D template-based tracking algorithms. Considered approaches are keypoint matching methods (such as FERNS [START_REF] Ozuysal | Fast keypoint recognition using random ferns[END_REF], SIFT [START_REF] Lowe | Distinctive image features from scale-invariant keypoints[END_REF] or SURF [START_REF] Bay | Speeded-up robust features (SURF)[END_REF]) or tracking approach based the minimization of the SSD (such as the ESM [START_REF] Benhimane | Homography-based 2D visual tracking and servoing[END_REF]). Other approach which maximized the mutual information shared by two images [START_REF] Dame | Accurate real-time tracking using mutual information[END_REF] can also be considered. The latter allows evaluation on both real scene and real image, and can also be considered model-based tracking approaches. In both cases, ground-truth are available which allows fair comparisons and qualitative evaluation of proposed approaches.

In this paper, we propose to benchmark two model-based tracking approaches developed in the INRIA Lagadic team [1]. The former is an extension to model-based 3D pose estimation [START_REF] Dame | A unified direct approach for visual servoing and visual tracking using mutual information[END_REF] of mutual-information based tracker proposed in [START_REF] Dame | Accurate real-time tracking using mutual information[END_REF]. The second one, derived from [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF][START_REF] Wuest | Tracking of industrial objects by using CAD models[END_REF], is an edge-based registration process than take advantage of GPU rendering capability. Both approaches are based on the virtual visual servoing framework proposed in [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF]. Virtual visual servoing is a formulation of a full-scale non-linear minimization based on the visual servoing scheme [START_REF] Chaumette | Visual servo control, Part I: Basic approaches[END_REF]. It established a duality between vision-based robot control and various computer vision problem such as pose estimation, calibration, motion estimation, etc. Since these are model-based tracking method, we consider the TrakMark benchmark with the same sequence for which groundtruth is available. Evaluation criterion are also discussed.

DATASET

In the TrakMark datasets, many kinds of image sequences are provided. We evaluated our tracking methods with one of the datasets named "Conference Venue Package 01". This dataset includes the 3D textured model of the ISMAR2009 conference venue generated by Ishkawa et al. [START_REF] Ishikawa | In-situ 3D indoor modeler with a camera and self-contained sensors[END_REF] and the computer-generated images with several camera motions. Especially, we selected the image sequence captured with parallel translation, panning and tilting of a camera because it was the most challenging motion of this TrakMark dataset for the evaluation.

METHODS

In this section, we provide the short descriptions of our tracking methods. The details will be presented in the future conferences and journals.

Mutual information based Tracking

Model based pose estimation using vision has been tackled using various feature types. We propose here to introduce the mutual information feature for this problem. Mutual information (MI) is defined thanks to the entropy H of two images I and I * and their joint entropy:

MI(I, I * ) = H(I) + H(I * ) -H(I, I * ) (1)
MI has been used for registration works [START_REF] Studholme | An overlap invariant entropy measure of 3D medical image alignment[END_REF][START_REF] Viola | Alignment by maximization of mutual information[END_REF] and more recently to track planes in image sequences [START_REF] Dame | Accurate real-time tracking using mutual information[END_REF]. In the latter work, a planar region is detected or selected in the first image of a sequence and defines the reference template of the plane to track. Then, for the following image of the sequence, a homography is incrementally computed to warp the template in the current image so that the mutual information shared by both reference and current image regions is maximal. This optimization process is done for each new image and the initial guess is the optimal warp for the previous image, since in tracking processes, the motion in images is assumed to be small for a short period of time. This feature has shown to be robust to noise, specular reflections and even to different modalities between the reference image and the current one.

We consider here an extension of [START_REF] Dame | Accurate real-time tracking using mutual information[END_REF] to the case of non planar model based pose estimation and tracking. This extension adapts the use of MI over SL(3) presented in [START_REF] Dame | Accurate real-time tracking using mutual information[END_REF] to the SE(3) space. It means that the parameters space is the full six 3D pose parameters (three translations and three rotations) whereas it was eight parameters to estimate in the previous work: the relative pose between current and reference cameras and the plane parameters, all up to scale.

We, hence, need to reformulate the cost function. goal is to perform the registration of the model with respect to the image and it can be formulated as the maximization of the mutual information shared between the input image I * and the projection of the model M . If γ are the intrinsic camera parameters (focal length and principal point) and r its pose, the pose estimation problem can be written as:

r = arg max r MI ( I * , I γ (M , r) ) . ( 2
)
Image I γ (M , r) is resulting from the projection of the model M at given pose r. From a first order Taylor expansion of the mutual information function at the current pose r, the link between the variation of the mutual information feature and the pose variation is expressed. The increment to apply to the pose is then obtained using a Newton's optimization like method.

To solve this function, a textured 3D model of the object to track is necessary and it has to be projected for each camera pose r. To generate images of the 3D model, we used OpenGL as a 3D renderer and more particularly the Ogre3D library [START_REF]Ogre3D, open source 3D graphics engine[END_REF]. OpenGL allows not only to generate photometric images but also depth images. More precisely, we obtain an image where each pixel contains the Z coordinate of the 3D point projected in this pixel. This is particularly interesting since the Z of each visible point appears in the Jacobian linking mutual information and pose variations.

The algorithm presented in Figure 2 sums up all the processes of the mutual information based pose estimation and tracking approach.

Depth and Texture Edge based Tracking

The second approach considers the use of a complete surface model, which can be textured or untextured. The method aims at realigning edges generated from the rendered model and edges extracted from the image by minimizing a criterion related to the error between these two kinds of edges. For this purpose, as in [START_REF] Wuest | Tracking of industrial objects by using CAD models[END_REF], at each acquired image, the model is rendered using hardwareaccelerated 3D graphics (GPU), through Ogre3D library, with respect to the pose computed for the previous image. The challenge is then to sort out visible and prominent edges from the rendered scene. We propose to determine these edges by combining information provided by both the depth buffer, and the rendered textures of the scene. From the depth buffer (see Fig. 3(a)), which corresponds to the depth values of the scene according to the camera viewing direction at each pixel point, we can extract changes or discontinuities that suit the visible geometrical properties of the scene. This is performed by applying a second order differential operator to the depth values. We also extract edges from the textures (Fig. 3(b)) by processing a classical Canny edge algorithm. Given the edge map of the complete scene (Fig. 3(c)), we can compute the 3D world coordinates of the edge points in the scene and so generate a 3D edge model, as we know the pose used for the projection and the depth of these points thanks to the z-buffer. As dealing with the whole edge map can be computationally heavy, we can sample it along i and j coordinates of the image in order to keep a reasonable number of these edge measurement points. From these points we need to find their corresponding edges in the image. In a similar manner to [START_REF] Wuest | Tracking of industrial objects by using CAD models[END_REF] and [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF], we perform a 1D search along the normal of the underlying edge, whose orientation is retrieved thanks to Sobel filters computed on the grey level image of the normal map of the scene. As a matching edge point in the image, we choose the gradient maximum along the scan line (Fig. 3(d)). Then, we compute the camera pose r which minimizes the errors between the projected edge measurement points p i (r) and the corresponding edge points p ′ i in the image. As the error we choose the distance between the projected 3D line L p i (r) underlying p i (r) and p ′ i . The criteria to minimize can be expressed as :

S = ∑ i ρ(d ⊥ (L p i (r), p ′ i )) (3)
where d ⊥ (L p i (r), p ′ i ) is the distance between a point p ′ i and the corresponding line L p i (r), and ρ is a robust estimator used to reject outliers. The optimization technique is then similar to the virtual visual servoing framework described in [START_REF] Comport | Realtime markerless tracking for augmented reality: the virtual visual servoing framework[END_REF].

EVALUATION CRITERIA

In [START_REF] Lieberknecht | Benchmarking template-based tracking algorithms[END_REF], the evaluation criterion was based on the reprojection error of four points located on the diagonal lines of a reference image. The ratio of tracked images against all input images was computed such that the number of tracked images was incremented when the reprojection error was less than a threshold. In the TrakMark datasets, some of them provide the ground truth of a 3 × 4 transformation matrix from the world coordinate system to the camera coordinate system for all images. Therefore, we can evaluate the accuracy of rotation and translation components of a camera pose in 3D space.

The provided transformation matrix for each image is composed of a 3 × 3 rotation matrix R and a 3 × 1 translation vector t as [R|t]. We compute a rotation matrix R c and a translation vector t c for each image and compared them with the ground truths as follows. For the rotation, the difference matrix R d with the ground truth rotation matrix R g is first computed:

R d = R g R T c . ( 4
)
The difference matrix is then decomposed into an axis and angle of rotation with Rodrigues' rotation formula [START_REF] Huynh | Metrics for 3D rotations: Comparison and analysis[END_REF]. Because we use the angle as a difference of two rotation matrices, we compute the angle such that

θ R d = arccos ( tr (R d ) -1 2
) . ( 5
)
For the translation, the 2-norm distance (Euclidean distance) d with the ground truth translation vector t g is computed such that

d = t g -t c . ( 6
)
Note that Huynh has reported that the best evaluation of 3D rotations was to use unit quaternions in [START_REF] Huynh | Metrics for 3D rotations: Comparison and analysis[END_REF]. However, we selected an angle component in the representation of rotation with an axis and angle as a criterion because it has more geometrical meaning.

EVALUATION PROCEDURE

The whole procedure of our evaluation is as follows.

1. Download dataset images, intrinsic parameters and the ground truth of extrinsic parameters for each dataset image 2. Convert the format of the 3D model into Ogre 3D [START_REF]Ogre3D, open source 3D graphics engine[END_REF] for our codes 3. Set the ground truth of an initial image in the dataset to our tracking methods 4. Compute camera poses for the rest of the images 5. Compute the difference between our result and the ground truth with the evaluation criteria.

RESULTS

We evaluated our two tracking methods with the criteria and discussed the results individually.

Mutual Information based Tracking

The nonlinear process of this method maximizes the mutual information, shared by the desired and current images, optimizing the pose of the virtual camera. It is not very clear to give mutual information results. Thus, we chose to display image differences between reference images from the dataset and images obtained at optimal poses computed thanks to our method. Image differences should be grey when both images are identical. Figure 4 shows some difference images at different locations along the trajectory. Since our virtual images are not rendered with the same 3D engine as the one used by TrakMark to obtain datasets, some rendering properties (texture filtering, colors, etc) may be different. That is why image differences are not perfectly grey at convergence. One can note that despite this problem, the registration succeeds. This is due to the mutual information feature which is robust to such issues, whereas more classical registration cost function, such as the sum of squared differences, are not. After having evaluated results qualitatively in images, the evaluation of estimations is done quantitatively in 3D. Figure 5 shows the estimated trajectory which is extremely close to the ground truth. This nearly perfect superimposition of both trajectories was clearly waited from the qualitative evaluation part (Fig. 4) since if desired and current images are perfectly aligned, it means the pose of the virtual camera is quasi identical to the desired one.

Estimations are also evaluated on the translational and on the rotational parts of each pose (Fig. 6).

Depth and Texture Edge based Tracking

For this method, we propose the comparison of the results of two approaches, one for which the generation of the edge map from the rendered model only relies on depth discontinuities, which are related to the geometrical properties of the scene, and the other for which this generation process is based on both texture and depth discontinuities, as presented in Section 3.2. Fig. 7 and Fig. 8 qualitatively show the performances of both approaches. The yellow points represent the edge measurement points reprojected with respect to estimated camera pose. We observe that relying only on depth discontinuities is not really suitable for this very textured scene, as ambiguities between these depth edges and the textures edges in the image occur. This leads to some local minima, making the tracking fail. When including the texture information, the tracking is properly performed throughout the whole sequence as seen on Fig. 8. For this approach, Fig. 9 shows the estimated trajectory which is very close to the ground truth and Fig. 10 shows the errors in terms of translation and rotation. This method is also quite computationally efficient as we reach a 13 frames per second for the first approach and 7 fps for the second.

CONCLUSION

We evaluated two visual tracking methods developed in the INRIA Lagadic team with a TrakMark dataset. One was an extension of 2D template based plane tracking with mutual information [START_REF] Dame | Accurate real-time tracking using mutual information[END_REF] to nonplanar object tracking, and the other was based on the minimization of reprojection error of 3D edges extracted from both texture and depth images. Since the ground truth provided in the dataset was the transformation matrix from the world coordinate system to the camera coordinate system, we use the error of 3D rotation and translation as evaluation criteria. Through the evaluations, we confirmed that the dataset was applicable to benchmarking visual tracking for both augmented reality and visual servoing.

Figure 1 :

 1 Figure 1: Example of images. The scene of "Conference Venue Package 01" was captured at the ISMAR2009 conference venue. The computer-generated images and the 3D textured model are provided.

Figure 2 :

 2 Figure 2: Synopsis of the mutual information pose estimation algorithm. The process loops until the mutual information between I and I * is stable.

Figure 3 :

 3 Figure 3: On (a) and (b) are represented the z-buffer and the textures of the rendered 3D model using Ogre3D, from which the edge map is generated (c). This edge map is then sampled to extract measurement points, reprojected on the current image and from them a 1D search along the edge normal is performed to find a matching edge point in the acquired image (d).

Figure 4 :

 4 Figure 4: Qualitative evaluation of MI based pose estimation convergence. Difference between some reference images and images generated at optimal estimated pose.

Figure 5 :

 5 Figure 5: Estimated trajectory using our mutual information based pose estimation (red) superimposed over the ground truth (blue).

Figure 7 :

 7 Figure 7: Tracking relying on depth edges. Until (c), the tracking is properly performed, but then the ambiguities between geometrical and textures edges appears, and lead to local minima (d).

Figure 8 :

 8 Figure 8: Tracking relying on depth and texture edges. The tracking is properly performed throughout the sequence.

Figure 9 :

 9 Figure 9: Estimated trajectory using the depth and texture edge based pose estimation (red) superimposed over the ground truth (blue).

Figure 10 :

 10 Figure 10: Estimation errors in (a) position and in (b) orientation over all the sequence, with respect to the ground truth, for the depth and texture edge based pose estimation.