T Sentenac
email: sentenac@laas.fr

F Bugarin

B Ducarouge
email: benoit.ducarouge@laas.fr

M Devy

Automated thermal D 3 reconstruction based on a robot equipped with uncalibrated infrared stereovision cameras

Keywords:

In many industrial sectors, Non Destructive Testing (NDT) methods are used for the thermomechanical analysis of parts in assemblies of engines or reactors or for the control of metal forming processes. This article suggests an automated multi-view approach for the thermal D 3 reconstruction required in order to compute D 3 surface temperature models. This approach is based only on infrared cameras mounted on a Cartesian robot.

The low resolution of these cameras associated to a lack of texture to infrared images require to use a global approach based first on an uncalibrated rectification and then on the simultaneous execution, in a single step, of the dense D 3 reconstruction and of an extended self-calibration. The uncalibrated rectification is based on an optimization process under constraints which calculates the homographies without prior calculation of the Fundamental Matrix and which minimizes the projective deformations between the initial images and the rectified ones.

The extended self-calibration estimates both the parameters of virtual cameras that could provide the rectified images directly, and the parameters of the robot. It is based on two criteria evaluated according to the noise level of the infrared images. This global approach is validated through the reconstruction of a hot object against a reference reconstruction acquired by a D 3 scanner.

Introduction

This article addresses the problem of a fully automated D 3 thermal reconstruction [START_REF] Maldague | Theory and Practice of Infrared Technology for Nondestructive Testing[END_REF][START_REF] Maynadier | One-shot measurement of thermal and kinematic fields: infrared image correlation (IRIC)[END_REF][START_REF] Vidas | Heatwave: a handheld 3D thermography system for energy auditing[END_REF] from sensors embedded on a robotic system. Such a method can be suitable for performing diagnostics on mechanical assemblies, such as nuclear reactors [START_REF] Reichle | Defining the infrared systems for ITER[END_REF], or for improving energy efficiency in building construction [START_REF] Ham | Epar: energy performance augmented reality models for identification of building energy performance deviations between actual measurements and simulation results[END_REF] or for monitoring forming processes [START_REF] Jeswiet | Asymmetric single point incremental forming of sheet metal[END_REF]. The first sub-problem that arises is to define the system architecture using heterogeneous sensors. The second sub-problem is selecting the dense D 3 reconstruction methods with overlaid thermal data. The third sub-problem is to make the inspection task totally automatic with a self-extended calibration (i.e. a calibration without a specific target, covering all the geometric parameters of the robotsensor system, including the intrinsic and extrinsic sensor parameters and the robot parameters).

The most conventional architecture is based on a D 3 laser scanner and infrared cameras mounted on a robot [START_REF] Gonzalez-Aguilera | Novel approach to 3D thermography and energy efficiency evaluation[END_REF][START_REF] Borrmann | A mobile robot based system for fully automated thermal 3D mapping[END_REF][START_REF] Nagatani | Three-dimensional thermography mapping for mobile rescue robots[END_REF]. To overcome the significant cost of the D 3 laser scanner, several articles [START_REF] Ham | Epar: energy performance augmented reality models for identification of building energy performance deviations between actual measurements and simulation results[END_REF][START_REF] Zalud | Fusion of thermal imaging and CCD camera-based data for stereovision visual telepresence[END_REF][START_REF] Vidas | 3D thermal mapping of building interiors using an rgb-d and thermal camera[END_REF] have suggested an architecture using only cameras. Inexpensive and readily available digital visible cameras (CCD camera, color camera, Kinect, etc.) give images processed by a D 3 modeler, while infrared cameras provide the thermal data mapped on the D 3 model. One successful system is the HeatWave system [START_REF] Vidas | Heatwave: a handheld 3D thermography system for energy auditing[END_REF][START_REF] Moghadam | Heatwave: the next generation of thermography devices[END_REF], i.e. a hand-operated device consisting of rigidly attached infrared and color cameras. These multisensory architectures face the difficulty of fusing D 3 data provided by a D 3 modeler and temperature data acquired by infrared cameras into a common coordinate frame. A joint geometric calibration of heterogeneous sensors [START_REF] Prakash | Stereo techniques for 3D mapping of object surface temperatures[END_REF] must be performed, which requires finding a pattern that is completely visible by both the D 3 sensors and the IR cameras. This could be a tricky task, because these heterogeneous sensors have different spectral sensitivities, spatial resolutions and fields of view. The ideal architecture is then only based on infrared cameras for a direct D 3 thermal reconstruction. Assuming that thermal methods already described in [START_REF] Sentenac | Bi-color near infrared thermoreflectometry: a method for true temperature field measurement[END_REF][START_REF] Sentenac | Noise effect on the interpolation equation for near infrared thermography[END_REF][START_REF] Gilblas | Quantitative temperature field measurements on a non-gray multi-materials scene by thermoreflectometry[END_REF] are not in the scope of this paper, the challenge is then to provide a dense image-based D 3 reconstruction [START_REF] Vidas | A mask-based approach for the geometric calibration of thermal-infrared cameras[END_REF][START_REF] Yu | Camera calibration of thermal-infrared stereo vision system[END_REF] with infrared cameras.

Several image-based D 3 reconstruction algorithms have been pro- posed using visible cameras. The first step, image registration, requires the detection and matching of features between images. Many feature detectors (e.g., Harris, SIFT, SURF, FAST, ORB, …) automatically and correctly extract and match interest points on infrared images. Then two classes of methods can simultaneously build a sparse D 3 model and the camera trajectory using only these matched interest points. The Robotics community has developed several Vision-based Simultaneous Localization And Mapping (VSLAM) techniques [START_REF] Durrant-Whyte | Simultaneous localization and mapping: part i[END_REF][START_REF] Rusu | 3D robotic mapping: the simultaneous localization and mapping problem with six degrees of freedom[END_REF], taking advantage of other proprioceptive data acquired from the robot (odometry, IMU…), but assuming generally that the intrinsic camera parameters are known. The Vision community has proposed Structure from Motion (SfM) [START_REF] Hartley | Multiple View Geometry in Computer Vision[END_REF] approaches (Bundler, OpenMVG…) in order to recover from an image sequence both the D 3 environment structure and the camera Motion; extrinsic and intrinsic camera parameters can be estimated simultaneously when computing the D 3 point positions. The recovered parameters should be consistent with the reprojection error (i.e., the sum of distances between the projections of each set of D 3 corresponding feature points and its corresponding image features). This minimization problem can be formulated as a non-linear least squares problem and solved from a Bundle Adjustment (BA) algorithm [START_REF] Lourakis | Sba: a software package for generic sparse bundle adjustment[END_REF]. Exploiting VSLAM or SfM methods, an accurate and dense D 3 model could be incrementally and gradually built and refined with, typically, a sequence of one thousand images, either from an Iterative Closest Point (ICP) algorithm (stereovision) or by a Multi-View Stereo (MVS) [START_REF] Furukawa | Accurate, dense, and robust multiview stereopsis[END_REF] technique (monocular vision).

The paper proposes an automated thermal D 3 reconstruction based on an architecture composed of a Cartesian robot equipped only with uncalibrated infrared cameras. The architecture requires a coupled method that deals simultaneously with a multi-view D 3 thermal re- construction and a self-extended geometric calibration. An infrared stereo vision rig provides a compensation to the lower spatial resolution of infrared images. It also improves the number of reconstructed points and thus the density of the D 3 model. Moreover, it gives an initial guess for the D 3 position of every point. For the first step of the method, a reasonable amount of stable and tractable matched points is obtained through a specific method for infrared images based on the phase congruency model [START_REF] Hajebi | Dense surface from infrared stereo[END_REF][START_REF] Hajebi | Structure from infrared stereo images[END_REF] which is combined with more classical feature detectors. With few and low-textured infrared images, the result would be limited to a sparse D 3 reconstruction. Next, the simultaneous reconstruction and self-calibration with uncalibrated infrared cameras is solved by the minimization of a cost function which integrates all the geometrical calibration parameters for both the cameras and the robot. Estimations variables are: four intrinsic parameters for each camera, six for the relative position and orientation between the two cameras and six for the rotational and the translational components of the Euclidean transformation. This latter transformation is named hand-eye, between hand (robot gripper) and eye (camera). The total number of parameters is twenty if it is assumed that geometrical distortions due to the lens are corrected beforehand. This assumption is a good trade-off between the accuracy and the computation time. Indeed, it decreases the accuracy, but it avoids additional degrees of freedom and high non-linearities in the geometric model which are consuming in computation times. For the targeted application in a robotic context, a real-time processing of the calibration, compatible with the speed of the robot, is preferred even if the accuracy is not optimal. Finally, these parameters have to be estimated on-line from features extracted and matched from images acquired on the unknown object from two or more positions of the robot. These positions are known accurately in the robot reference frame from its forward kinematics.

A projective rectification method is first introduced for improving the matching problem between pixels on the left and right images, limiting the search space from the whole image to only a line. Another objective is to reduce the number of parameters of the cost function. The projective rectification method, applied to uncalibrated stereovision infrared cameras, requires a specific algorithm. Two homographies, applied to rectify the left and right images, have to be estimated in a single step to cope with the problem of noisy and lowtextured infrared images. A new cost function is proposed which is minimized under geometric constraints by a non-linear optimization process. These constraints are defined to keep the structure and the skewness of the images, which are essential to preserving their geometries.

The set of points matched from rectified images is then used to perform simultaneously the extended self-calibration and the multiview infrared D 3 reconstruction. It is based on the minimization of two different functions depending on the observed noise level in infrared images. The first objective function is based on the minimization of the reprojection errors. When the noise level increases, a second objective function is expressed in the projective space which takes advantage of the epipolar constraint between two images. This second objective function is expressed using the intrinsic camera parameters and the essential matrix between two robot positions which depends itself on the hand-eye parameters and the robot motion. Fig. 1 summarizes the suggested coupled method which merges an extended self-calibration and a multi-view infrared D 3 reconstruction applied on rectified images.

The paper is organized as follows. Section 2 briefly describes existing literature on rectification and places our suggested method with uncalibrated infrared images in this context. The method is then described and the results are compared to those in the literature. Section 3 outlines the formulation of the multi-view infrared D 3 reconstruction simultaneous to the extended self-calibration. The method is evaluated on synthetic data. Finally, Section 4 summarizes the results of the fully automated D 3 reconstruction performed from multiple views acquired by an uncalibrated infrared stereo rig mounted on a Cartesian robot. The whole approach is evaluated comparing the D 3 model of a re- constructed object with a reference CAD model.

Suggested rectification method of uncalibrated cameras

The suggested rectification method, applied to uncalibrated stereovision cameras, takes the advantage of calculating the homographies, projective transformations applied to rectify images, without a previous calibration and in a single step to cope with the problem of noisy and low-textured infrared images. These homographies are then calculated with only one non-linear optimization process under geometric constraints.

The section begins with a short description of the background for calculating the homographies and works related to the rectification problem. The suggested projective rectification method is then introduced and the non-linear optimization process under geometric constraints is detailed. These constraints are defined to keep the structure and the skewness of the images, which are essential to preserving their geometries. This property is essential for the self-calibration introduced in the next section. Finally, the method is evaluated by comparison with conventional methods that only work on visible images to prove its effectiveness even under these conditions.

Background and works related to the rectification problem

The rectification process reduces the two-dimensional matching problem on stereo images to a one-dimension matching problem. Using epipolar geometry, it consists of aligning the epipolar lines to make them parallel to the horizontal axis of the image. For uncalibrated cameras, knowledge of the epipolar geometry is packaged in the Fundamental matrix. The problem is then the computation of two projective transformations (homographies) from the Fundamental matrix to align the epipolar lines parallel to the horizontal image axis. The following paragraphs detail the computation of homographies and discuss the approach used to rectify uncalibrated images.

Background to epipolar geometry

Epipolar geometry defines the geometry between a pair of images I and I′ from two stereoscopic cameras or two different locations of a mobile camera. Let Q,1 a D 3 point, be simultaneously seen by two pinhole cameras in D 3 space. Let c and ′ c the optical centers of these two cameras. Let q and (resp. ′ q ) the projections of Q through c and ′ c in images I (resp. and I′). Q c , and ′ c define an epipolar plan in D 3 space, denoted by P . The left epipolar line ′ l q (resp. right epipolar line ′ l q ) in I (resp. I′) is defined by the intersection of P and I (resp. P and I′). By geometric construction, q has to be on the right epipolar line ′ l q (and resp. ′ q on the left epipolar line ′ l q ). This constraint is the epipolar constraint: for a given point I ∈ q , its corresponding point

I ′ ∈ ′ q
lies on its epipolar line ′ l q , i.e. ′ ′ = ⊤ q l 0 q . Similarly q lies on the epipolar line ′ l q , i.e. = ⊤ ′ q l 0 q . Because the relationships between retinal coordinates of corresponding points ′ q q ( , ) and their epipolar lines ′ ′ ( ) l l , q q are projective linear, the epipolar constraint can be rewritten as follows:

′ = ⊤ q F q 0. ( 1
)
where

F is Fundamental Matrix ( M ∈ F ( ) 3
). It encapsulates the projective motion between two uncalibrated perspective cameras. The epipoles e and ′ e are points which satisfy the following equation:

= ′ = ⊤ Fe F e 0 . 2 (2)
The epipole e (resp. ′ e ) is the intersection of all the epipolar lines in- cluded in I (resp. I′). Moreover, Eq. ( 2) implies that the rank of F is lower or equal to two. F is then defined up to a scale factor. It theo- retically depends upon seven independent parameters. More details about epipolar geometry and the fundamental matrix are provided in the book [START_REF] Hartley | Multiple View Geometry in Computer Vision[END_REF].

Computation of homographies

Using epipolar geometry, the rectification provides corresponding epipolar lines which are collinear and parallel with the x-axis. The rectified fundamental matrix F 0 is then expressed as follows:

= ⎡ ⎣ ⎢ - ⎤ ⎦ ⎥ F 0 0 0 0 0 1 0 1 0
I I ′ ∈ × ′ q q ( , ) to I I ∼ ∼ ′ ∈ × ′ ∼ ∼ q q ( , )
as:

∼ = ∼ ′ = ′ ′ q Hq q Hq ; (3)
It follows from Eq. (1) that:

∼ ′ ∼ = ′ ′ = ⊤ ⊤ ⊤ = q F q q H F Hq 0. F 0 0 (4)
Finally, H and ′ H are compatible homographies only if they satisfy the following equation:

= ′ ⊤ F H F H. 0 (5)
The practical computation of H and ′ H is then achieved by solving [START_REF] Ham | Epar: energy performance augmented reality models for identification of building energy performance deviations between actual measurements and simulation results[END_REF]. However, due to the first row of F 0 , the pair of homographies is not unique. This remaining degree of freedom can introduce an undesirable distortion to the rectified images. Methods suggested in the literature to reduce distortions can be divided into two classes. The first class assumes that F is fixed and the second class F is implicitly recomputed. The first class is divided into two steps. A preliminary accurate estimation of F and a calculation of ′ H H ( , ) are first introduced from Eq. ( 5). The distortions are then corrected by applying symmetrical [START_REF] Devernay | Vision stéréscopique et propriétés différentielles des surfaces[END_REF][START_REF] Loop | Computing Rectifying Homographies for Stereo Vision[END_REF] or independent [START_REF] Trivedi | Estimation of stereo and motion parameters using a variational principle[END_REF][START_REF] Hartley | Theorie and practice of projective rectification[END_REF][START_REF] Gluckman | Rectifying transformations that minimize resampling effects[END_REF][START_REF] Mallon | Projective rectification from the fondamental matrix[END_REF][START_REF] Laganière | Projective rectification of image triplets[END_REF]

∑ ′ + ′ ∼ ∼ ′ = ⊤ d d q Fq F q q min ( , ) ( , ). k
N k k k k H H , 1 (6)
where

′ = … q q ( , ) k N 1
is a set of N matched points between the images I and I′.

′ ∼ d q Fq ( , ) k k (resp. ′ ∼⊤ d F q q ( , ) k k
) is the distance from the point ′ q k (resp. q k ) to the epipolar line ∼ Fq k (resp. ′ ∼⊤ F q k ). However, this mini- mization alone is not enough to overcome the distortions problem which can be fixed thanks to the specific parametrization of the homographies H and ′ H [START_REF] Monasse | Three-step image rectification[END_REF][START_REF] Zilly | Three-step image rectification[END_REF]. These parametrizations go hand in hand with assumptions on the parameters of F. For instance, in articles [START_REF] Monasse | Three-step image rectification[END_REF][START_REF] Zilly | Three-step image rectification[END_REF], the authors assume that the principal point is centered and the aspect ratio is equal to one (so that the skew is equal to zero).

Formulation of the uncalibrated rectification

The suggested formulation of the uncalibrated rectification is consistent with the second class of rectification methods, which works in only one step to take as fully as possible into account the low number of points matched in infrared images (images with little texture and a weak spatial resolution). Iteration of many steps could propagate and amplify the uncertainty due to low numbers of matched points. The suggested method then consists in computing the Fundamental matrix F and the rectification homographies H and ′ H in a single step without prior assumptions on the geometry and on the parametrization of H and ′ H . The contribution is then to achieve the calculation of homographies by applying the non-linear objective function of Eq. ( 6) under constraints to minimize the geometrical distortions induced by H and ′ H . Introduced by [START_REF] Mallon | Projective rectification from the fondamental matrix[END_REF] and commonly used in the literature [START_REF] Monasse | Three-step image rectification[END_REF][START_REF] Zilly | Three-step image rectification[END_REF], criteria based on the aspect ratio and the orthogonality of the images are applied after the rectification process to control loss or pixels creation. The main idea of the suggested computation of the homographies is to include these classical criteria as constraints of minimization process in order to keep the structure and the skewness of the image. The constraint space is then based on the aspect and size ratio and the orthogonality of the image which are invariant to affine transformations [START_REF] Loop | Computing Rectifying Homographies for Stereo Vision[END_REF].

Definition of the non-linear objective function

The distance d from a point ′ q to its corresponding epipolar line ∼ Fq of Eq. ( 6) is defined in2 as follows:

′ = ′ ∼ ∼ ∼ ⊤ d π q Fq q Fq Fq ( , ) | | ‖ ( )‖ 2 2 (7) with → π x x x x x : ( , , ) ( , ) 1 2 3
1 2 the canonical projection. Conversely, the distance for a point q to its corresponding epipolar line ′ Fq can also be defined.

Eq. ( 6) is then rewritten symmetrically on both images as follows:

∑ ′ + ′ ∼ ∼ ∼ ∈ = ⊤ ⊤ π π q Fq Fq F q min | | ‖ ( )‖ ‖ ( )‖ P C k N k k k k 1 2 2 2 2 2 (8)
where

= ′ ∈ × ≃ P H H ( , ) 9 9
1 8 is the vector of the nine parameters of H and the nine parameters of ′ H . This vector is minimized on a constraint space ⊂ C 18 to find a pair of homographies reducing the geometric distortions. The notation of only one constraint space, denoted by C, actually unifies and merges heterogeneous space constraints (i.e. with different units) applied either on I or on I′. C is then the intersection of several constraint spaces which are detailed in the next subsection.

Space constraints definition

The need to minimize distortion requires keeping the image structure and the skewness. Hence, each pixel of the original image should map to a single pixel in the rectified images (pixel creation and loss should be minimal). In the article [START_REF] Hartley | Theorie and practice of projective rectification[END_REF] the creation or the loss of pixels is modeled by the change in local area of a patch around the point before and after the rectification. This change can be determined by the numerical properties of the Jacobian homography. These numerical properties can be controlled by matrix operators based on the determinant [START_REF] Gluckman | Rectifying transformations that minimize resampling effects[END_REF] or on singular values [START_REF] Mallon | Projective rectification from the fondamental matrix[END_REF]. Thus, minimizing these operators involves reducing the distortion in the whole image and their implementation then becomes very time-consuming in terms of computation. Our approach is to define a set of numerical stable constraints C which only enforces the physical properties of the image. The image distortions are then defined as modifications of the image structure: aspect ratio, size (width and height) and orthogonality of the image.

Aspect ratio criterion. The constraints of the image structure can be first quantified by the invariance of the ratio of image diagonals computed from the aspect ratio between original and rectified images. Ideally, the aspect ratio is equal to 1 (see Fig. 2(a)). The four image corners are sufficient to compute this criterion. Let p p p p { , , , }

‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ 1 32 2 4 2 1 3 2 2 4 2 (9) = ∈ ∈ -∊ + ∊ C P E P { | ( ) [ 1 ,1 ]} a a 18 2 (10)
where ∊ is the tolerant variation of the aspect ratio. Orthogonality criterion. The angles invariance between original and rectified images can be modeled by the orthogonality which is ideally equal to °90 (see Fig. 2(b)). Article [START_REF] Loop | Computing Rectifying Homographies for Stereo Vision[END_REF] considers that applying con- straints on the four side mid-points preserves the orthogonality of the image. With the same previous notations, orthogonality criterion E o and the orthogonality constraint C o are then defined by:

→ → ′ E P γ γ : (acos( ), acos( )) o 18 2 (11)
= ∈ ∈ °- °+ C P E P θ θ { | ( ) [ 9 0 ,9 0 ]} o o 18 2 (12)
where ′ γ γ

, are the angle defined by:

̃̃̃̃̃̃̃= - + - - --- - - γ
q q q q q q q q q q q q

‖ ‖ ‖ ‖ ‖( ) ( )‖ 2‖ ‖‖ ‖ 1 3 2 2 4 2 1 3 2 4 2 1 3 2 4 ( 13
)
̃̃̃̃̃̃̃′ = ′-′ + ′-′ -′-′ -′-′ ′-′ ′-′ γ
q q q q q q q q q q q q

‖ ‖ ‖ ‖ ‖( ) ( )‖ 2‖ ‖‖ ‖ 1 3 2 2 4 2 1 3 2 4 2 1 3 2 4 ( 14
)
and θ is the tolerant variation of the orthogonality. Size ratio criteria. The invariance of the width and the height of images between original and rectified images can be only computed by a set of points defined by the side mid-points of the horizontal side for image I (resp. I′), q q , 1 3 (resp. ′ ′ q q , 1 3 ) and vertical image side, q q , 2 4

(resp. ′ ′ q q , 2 4 ). In the rectified images, the corresponding points are defined as: ̃̃̃q q q q { , , , }

4 and ̃̃̃′ ′ ′ ′ q q q q { , , , }

4 . The size ratio E s and the size ratio constraint C s are then defined by: → E :

s 18 4
(15)

̃̃̃̃̃̃̃→ = ∈ ∈ - + ′ ′ ′ ′ ′ ′ ′ ′ - - - - - - - - ( ) P C P E P δ δ , , , { | ( ) [ 1 ,1 ]}
s s q q q q q q q q q q q q q q q q

‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ ‖ 18 4
1 32

1 32 1 3 2 1 3 2 2 42 2 42 2 4 2 2 4 2 (16)
where δ the tolerant variation of the image size. Ideally, E P ( )

s
is equal to 1 (see Fig. 2(c)).

Finally, the constraints space C is the intersection of the previous constraints, i.e.

= ∩ ∩ C C C C a s o .
Although C merges heterogeneous constraints, the optimization method used to solve Eq. ( 8) considers each constraints separately and then deals with constraints of different units. The minimization of Eq. ( 8) is performed by a local optimization method which is the Sequential Quadratic Programming algorithm (S.Q.P.) [START_REF] Lawrence | A computationally efficient feasible sequential quadratic programming algorithm[END_REF][START_REF] Lawrence | User's Guide for cfsqp version 2[END_REF][START_REF] Panier | On combining feasibility, descent and superlinear convergence in inequality constrained optimization[END_REF]. A good initial estimate P 0 is then necessary to guar- antee the convergence to a right solution. The choice of the constraints space C affects the ease of finding this initial estimate. Therefore, the nine initial parameters of the homographies have to satisfy the constraints C. The identity homographies are the simplest initial parameter vector that satisfies the constraints. Indeed, the image without bias (distortion) that satisfies the constraints is the original image. Moreover, the identity has the property (that we verified experimentally on all our image tests) of being situated in a basin of attraction of a minimum of physically practicable premises.

Validation experiment

The suggested rectification method was evaluated on images from the Mallon's test set [START_REF] Mallon | Projective rectification from the fondamental matrix[END_REF] which are taken by the same camera under a fixed lens configuration. Each set (Arch, Boxes and Drive) consists of two RGG images with × 640 480 pixels resolution. The lens distortion has been removed and the ground truth is not available. Results carried out by the suggested rectification method were compared to Monasse's method [START_REF] Monasse | Three-step image rectification[END_REF], which also works in a single step. This article [START_REF] Monasse | Three-step image rectification[END_REF] also compares the results of other conventional methods implemented by Loop and Zhang [START_REF] Loop | Computing Rectifying Homographies for Stereo Vision[END_REF], Hartley [START_REF] Hartley | Theorie and practice of projective rectification[END_REF] and Mallon [START_REF] Mallon | Projective rectification from the fondamental matrix[END_REF].

Evaluation criteria

The rectification performance is concerned with quantifying the rectification error E r . For each point of one image, it is the distance between its corresponding point and epipolar line in the other image. It is represented by the mean (mean) and the standard deviation (σ ). Then, the other criteria are the aspect ratio E a and the orthogonality E o . Ideally, the aspect ratio E a must be 1 and the orthogonality E o must be °90 .

Results

Table 1 gives the relative performance on the three set of images (Arch, Boxes and Drive) from the Mallon's dataset and the comparison with Monasse's method.

The values of the orthogonality E o (respectively the aspect ratio E a ) criterion are similar for both methods and close to the reference value of °90 (respectively of 1). The image structure and the skewness are well preserved and remain invariant. The rectified images can be applied in a geometric calibration process. The suggested method provides a smaller rectification error E r (mean and standard deviation) than the reference method (notably for the Drive set). The method is based on no prior assumptions on the parametrization of the homographies H and ′ H . No geometry is assumed and all parameters of homographies can be estimated with all possible values. However, at the same time the estimation is constrained to preserve the structure of the image and the generality of the solution. Finally, it should be emphasized that the identity provides a good initial estimate and convergence to a relevant solution.

Suggested extended self stereovision hand-eye calibration

The extended self hand-eye calibration with one uncalibrated camera is described with ten parameters: four intrinsic parameters of the camera (focal length, optical center), and six ones for the rotational and the translational component of the Euclidean transformation between hand and eye. With a stereovision rig, ten parameters must be added: four intrinsic parameters for the second camera, and six for the relative position and orientation between the two cameras. The previous rectification step decreases the number of parameters to be estimated, from twenty to thirteen (six intrinsic parameters for cameras, six parameters for the rigid hand-eye transformation and one for the rigid rectified stereovision configuration). It is assumed that image distortions are neglected. The estimation problem of the high number of parameters of the extended self stereovision hand-eye calibration is highly relevant in the presence of noise in the matching points between each camera and each displacement.

A first paragraph deals with the background and related works of hand-eye calibration. The second paragraph introduces the suggested extended self stereovision hand-eye calibration which is a simultaneous calibration of hand-eye, camera-intrinsic and stereo parameters. In the third paragraph, the problem is formulated as two non-linear optimization processes in relation to the noise on the matching points. In the last paragraph, experimental results on synthetic data show the influence of the noise on the accuracy of the D 3 reconstruction.

Background and works related to hand-eye calibration

Hand-eye calibration is the computation of the unknown transformation, named X between the camera frame and the frame of the hand of the robot. The hand-eye problem is solved knowing the displacement of the camera, named A, in a reference frame and the displacement of the hand of the robot, named B, between two positions i and j.

Let R j c (resp. R ′ j c ) be the coordinates system associated to c (resp. ′ c ) at the j th position and R w , the world coordinates system associated to the initial position of the robot. The hand-eye calibration then consists in estimating the parameters of the homogeneous Euclidean transformation, X, between the robot (hand) and the camera (eye) coordinate frame. The main assumption is that the single camera is rigidly coupled to the robot. As shown in Fig. 3, the robot performs discrete displacements from the position i to the position j, which are encoded through the euclidean transformation, B ij . The displacements are also measured by the camera c as a transformation A ij in its own coordinate frame. The parameters of hand-eye transformation, X, are then calculated ac- cording to the graph of coordinate frames of Fig. 3 and by solving the following equation: Fig. 3. Graph of transformations. Transformation matrix B ij displaces the coordinates system of the robot, transformation matrix A ij displaces the co- ordinates of the camera c, hand-eye transformation X displaces the coordinates system of the robot (hand) on coordinates of the camera (eye).

= - A XB X ij ij 1 (17)
By referring rigid transformations by × 4 4 homogeneous matrices composed of a rotation and a translation R t ( , ), Eq. ( 17) is rewritten as follows:

= ⊤ R R R R A X B X ij ij (18) = - - ⊤ ( ) t R t R R R I t A XB X B X R X ij ij ij ij A (19)
Several methods exist to solve the hand-eye calibration problem with different optimization processes. A general taxonomy of the classical approaches is provided in article [START_REF] Shah | An overview of robot-sensor calibration methods for evaluation of perception systems[END_REF] and in the focus of our contribution, they can be split into two classes: non-unified and unified methods. Using the decomposition [START_REF] Yu | Camera calibration of thermal-infrared stereo vision system[END_REF], the non-unified methods separately estimate the translation t X and the rotation R X . Several para- metrizations can be used to estimate the rotation: the angle-axis parametrization of the group of rotations SO (3) [START_REF] Tsai | Real time versatile robotics hand-eye calibration using 3D machine vision[END_REF][START_REF] Shiu | Calibration of wrist-mounted robotic sensors by solving homogeneous transform equations of the form AX=XB[END_REF], quaternions [START_REF] Chou | Finding the position and orientation of a sensor on a robot manipulator using quaternions[END_REF] or dual quaternions [START_REF] Horaud | Hand-eye calibration[END_REF]. Moreover, with noisy inputs, according to investigations carried out in the article [START_REF] Chen | A screw motion approach to uniqueness analysis of head-eye geometry[END_REF], unified methods with a simultaneous estimation of both rotation and translation are highly relevant to balance the estimation of the rotation part of X. It has also been shown in the article [START_REF] Horaud | Hand-eye calibration[END_REF], that camera intrinsic parameters are not independent from hand-eye parameters. The extended unified methods thus include both hand-eye and camera parameters [START_REF] Zhao | A flexible method combining camera calibration and hand-eye calibration[END_REF] (e.g. camera pose, intrinsic parameters and lens distortion [START_REF] Malti | Hand-eye and radial distortion calibration for rigid endoscopes[END_REF][START_REF] Malti | Hand-eye calibration with epipolar constraints: application to endoscopy[END_REF]) which are simultaneously estimated. The optimization for solving Eq. ( 18) is linear for unified methods and non-linear for non-unified methods. According to prior knowledge of the external camera calibration, the conventional non-linear optimization process is based on the norm two:

-

A X XB ‖ ‖ ij ij 2 .
The most recent non-linear optimization processes use branch-and-bound [START_REF] Heller | Structure-from-motion based handeye calibration using ∞ l minimization[END_REF][START_REF] Seo | A branch-and-bound algorithm for globally optimal calibration of a camera-and-rotation-sensor system[END_REF][START_REF] Heller | Globally optimal hand-eye calibration using branch-and-bound[END_REF], Second Order Cone Programming [START_REF] Zhao | Hand-eye calibration using convex optimization[END_REF] or polynomial optimization [START_REF] Heller | Hand-eye and robot-world calibration by global polynomial optimization[END_REF] techniques. For extended methods, all parameters are mainly simultaneously estimated with the Structurefrom-Motion (SfM) technique. In the article [START_REF] Andreff | On-line hand-eye calibration[END_REF] camera pose parameters are recovered up to a scale factor and in articles [START_REF] Horaud | Hand-eye calibration[END_REF][START_REF] Daniilidis | The dual quaternion approach to hand-eye calibration[END_REF][START_REF] Schmidt | Calibration free hand-eye calibration: a structure from motion approach[END_REF] they are explicitly included. Branch-and-bound algorithms have recently been introduced to estimate camera poses with a ∞ L -norm formulation [START_REF] Seo | A branch-and-bound algorithm for globally optimal calibration of a camera-and-rotation-sensor system[END_REF] and with a hypothesis that assumes all the translations are equal to zero. For extended methods with several parameters [START_REF] Malti | Hand-eye and radial distortion calibration for rigid endoscopes[END_REF], the bundle adjustment technique is the most promising approach to retrieve all parameters simultaneously.

Extended self stereovision hand-eye calibration suggested for uncalibrated infrared cameras rig

For uncalibrated stereovision infrared cameras, the extended self stereovision hand-eye calibration suggested thus relies on a generalized extended unified method. Indeed, the suggested method performs simultaneously, and in a single step, the extended self-calibration and the multi-view infrared D 3 reconstruction. The first paragraph recalled the parameters of the calibration after the previous rectification step, which are estimated using a bundle adjustment technique. The second paragraph introduces two objective functions to solve the double problem of calibration and reconstruction. The first one aims at minimizing the reprojection errors (in pixels). The second one is expressed in the projective space, from the epipolar constraints between two images acquired by the same camera from two robot positions. These two functions make explicit all the parameters of the transformations (intrinsic camera parameters, hand-eye parameters and robot motion parameters). Finally, the third paragraph validates the extended self stereovision hand-eye calibration using simulated data generated with different noise levels. The last paragraph concludes with the best objective function according to the noise level.

Parameters of extended self stereovision hand-eye calibration

As illustrated in Fig. 3, the rigid body transformation between R j c and R ′ j

c is denoted by = T R t [ ]. Given a D 3 point R Q ∈ w w
, its projections q j and ′ q j in I j and I′ j after a displacement B j 0 of the robot from the initial position are:

Q = - q KXB X j j w 0 1 (20)
Q ′ = ′ - q KTXB X j j w 0 1 (21)
where the matrix T is reduced to a single translation over the horizontal axis for a rectified stereovision rig:

= × T I t [ ] 3 3
with = ⊤ t t ( 0 0) . t is the baseline of the rectified stereovision rig. K and ′ K are the matrices of intrinsic parameters which are given by:

= ⎛ ⎝ ⎜ ⎞ ⎠ ⎟ ′ = ⎛ ⎝ ⎜ ⎜ ′ ′ ′ ′ ⎞ ⎠ ⎟ ⎟ α u α v α u α v K K 0 0 0 0 1 ; 0 0 0 0 1 u v u v 0 0 0 0 (22) where u v , 0 0 (resp. ′ ′ u v , 0 0 ) represent the centers of images I (resp. I′), α α , u v (resp. ′ ′ α α , u v
) represent the number of pixels per millimeter for I (resp. I′). For an uncalibrated rectified stereovision rig, the two rec- tified images are coplanar, abscissa are collinear and the epipolar lines are aligned between the two images. Consequently, = ′ v v 0 0 and = ′ α α v v . Finally, projection functions are fully described with only thirteen parameters: six intrinsic parameters ,,,, )

′ ′ α α α u u v ( ,
u u v 0 0
0 , the baseline t and six for the hand-eye transformation (three for the rotation and three for the translation).

Discussion on the definition of the objective functions

Different approaches can be exploited to design the objective functions.

The first approach is inspired by the estimation pose problem in which the D 3 location of observed points is estimated in terms of re- projection errors. The approach operates the geometry of the stereo image sequences. The objective function minimizes the reprojection errors (see articles [START_REF] Wei | Active self-calibration of robotic eyes and hand-eye relationships with model identification[END_REF][START_REF] Jordt | Automatic high-precision self-calibration of camerarobot systems[END_REF]) between the measured points by the camera and points from the projection model (see Eq. ( 20)) in a recursive bundle adjustment. However, the objective function is symmetric to calculate the two distances in each image and the number of parameters involved is high. Instead of minimizing the reprojection error, it is more convenient to calculate the reconstruction error of D 3 points built from the stereo image sequence between two robot positions. Each D 3 point is expressed in the robot coordinate frame versus the intrinsic, the stereo and hand-eye parameters and the displacement of the robot. The estimation of intrinsic, stereo and hand-eye parameters is then performed in a single step starting from approximate initial guesses.

The second approach is based on the epipolar constraint geometry of the stereo image pairs and the stereo images sequences. It has some theoretical drawbacks compared to the bundle adjustment, in that it does not provide as much information and cannot achieve the same level of accuracy. However, the first benefit is a smaller parameter space. Moreover, methods based on epipolar geometry do not involve the D 3 position of the observed object point, i.e. the epipolar constraint decouples the extrinsic camera parameters from the D 3 structure of the observed object. Eq. ( 2) only depends on the D 2 coordinates image. Eq. ( 20) additionally requires the depth of each observed point. If the matching point is not accurate, the error is not propagated in the D 3 space. The epipolar constraint is then calculated between two images acquired by the same camera from two robot positions. This constraint is expressed using the Fundamental matrix, i.e. using the intrinsic camera parameters and the essential matrix which itself depends on the hand-eye parameters and the robot displacement. This second suggested objective function is thus expressed using the intrinsic camera parameters and the essential matrix between two robot positions. All the parameters are also estimated in a single step.

Based on these different approaches, our approach combines both the geometry of the stereo image sequences and the epipolar constraint geometry of the stereo image pairs with two objective functions being derived in relation to the noise level of the images. Given noisy image points, the minimization of these objective functions provides intrinsic, stereo and hand-eye parameters and starts from approximate initial guesses. The main contribution of the suggested extended self stereovision hand-eye calibration is the proposal of two objective functions which depend on the noise level on the infrared images. .

j j 0 . Consider a static D 3 point R Q ∈ k w
w from the initial position of the robot. It can be expressed in R j c using the following equation:

Q Q = XB . k j j k w (23)
However, the localization of Q k w remains the same for two successive displacements:

Q Q = + - - + -- B X B X , j k j j k j 1 1 1 1 1 1 ( 24
)
The aim is then the minimization of the D 3 locations (in mm) for the successive displacements:

Q Q ∑ ∑ - = - = + - - + -- B X B X min ‖ ‖ , j m k n j k j j k j X 1 1 1 1 1 1 1 1 1 2 (25)
with X the vector of the hand-eye transformation. Recall that all the B j matrices are known. However, the problem expressed in equation ( 25) is highly overdetermined. This equation can be rewritten by incorporating the intrinsic parameters of the two cameras, K and ′ K , and the transformation between both cameras, T.

The first step is the introduction of equality between ′ q q ( , ) k j k j , the projections of Q k w in the images taken at j th displacement, which is given by the following equations:

Q = q KXB k j j k w (26) Q ′ = ′ q KTXB . k j j k w ( 27
)
In the next step, Eq. ( 25) is rewritten so as to minimize the reprojection error (in pixels) as follows:

Q Q ∑ ∑ - + ′ -′ ′ = = q KXB q KTXB min ‖ ‖ ‖ ‖ . j m k n k j j k w k j j k w X K K T , , , 1 1 2 2
(28)

Note that solving Eq. ( 28) means knowing precisely the D 3 locations of the set Q ( ) k w k . However, without a calibration object and without knowing the D 3 points, the third step is the estimation of these D 3 lo- cations jointly to ′ X K K , , and T as follows:

Q Q Q ∑ ∑ - + ′ -′ ′ = = q KXB q KTXB min ‖ ‖ ‖ ‖ . j m k n k j j k w k j j k w X K K T , , , , 1 1 2 2
k w [START_REF] Hartley | Theorie and practice of projective rectification[END_REF] This formulation, known as bundle-adjustment, requires a initial estimate of Q = ( )

k w k n 1.
. . The fourth step is to recover these D 3 locations (unknown in our context) from the matched pixels by inverting Eq. [START_REF] Devernay | Vision stéréscopique et propriétés différentielles des surfaces[END_REF]. A D

3 location Q k w is then seen as Q [ ] k w j , the intersection of two rays Q ∼ k j and Q ′ ∼ k j
given by:

Q = - - ∼ ⊤ ⊤ - ⊤ ⊤ ⊤ R R K q R R t R t k j k j B X B X X B B 1 j j jj ( 30
)
Q ′ = ′ - ∼ ⊤ ⊤ ⊤ - ⊤ ⊤ ⊤ R R R K q R R R t k j k j B X T B X T T 1 j j (31) - - ⊤ ⊤ ⊤ R R t R t , B X X B B j j j
Assuming that the localization of Q k w remains the same for all positions of the rig, the errors ε k j of two successive triangulations (see Fig. 4) are minimized in the following formulation:

∑ ′ = ε min . k n k j X K K T , , , 1 (32)
with

Q Q = - + ε ‖[ ] [ ] ‖ k j k w j k
w j 1 2 (see Fig. 4).

The optimization method, noted (M 1 ), provides an estimated vector p of thirteen parameters (6 intrinsic parameters ( ′ K K , ), the baseline (T) and 6 parameters of the hand-eye transformation (X)).

Objective function from the epipolar constraint

These parameters can be also estimated in the second optimization problem by using the epipolar geometry. As explained in Section 2, q k j lies on the epipolar line ′ + l q k j 1 . The fundamental matrix F j associated to camera c and camera c moved by applying the displacement

+ A j j
, 1 is estimated by solving:

∑ + = + ⊤ ⊤ + π π q Fq F q F q min | | ‖ ( )‖ ‖ ( )‖ k n k j j k j j k j j k j F 1 1 2 2 1 2 j ( 33
)
In the second step, F j is decomposed as follows:

= -⊤ - F K E K j j 1 ( 34
)
where E j is the transformation between the retinal plane of c and the retinal plane of c moved by applying the displacement A j . Then, the cost function of Eq. ( 33) can be rewritten as follows:

∑ ∑

+ = - = + ⊤ -⊤ - -⊤ - -⊤ ⊤ - + π π q K E K q K E K q K E K q | | ‖ ( )‖ ‖ ( )‖ j m k n k j j k j j k j j k j 1 1 1 1 1 2 1 2 1 1 2 (35)
Note that the rigid transformation R t [ ] j j E E extracted from the essential matrix E j is strictly equal to + A j j , 1 :

= = + + - R t A XB X [ ] j j jj jj E E , 1 , 1 1 (36)
Then the rotation R j E and the translation t j E can be deduced from Eq. ( 18). The same idea can be applied to the fundamental matrix ′ F j Fig. 4. Definition of the error ε k j (in red) minimized in Eq. ( 32). This error is defined by the distance between Q [ ] k w j (the intersection of the raysin black -

Q ∼ k j and Q ′ ∼ k j ) and Q + [ ] k w j 1 (the intersection of the rays -in green -Q ∼ + k j 1 and Q ′ ∼ + k j 1 ).
(For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) associated to camera ′ c and camera ′ c moved by applying the dis- placement

+ A T j j , 1 : ′ = ′ ′ ′ -⊤ - F K E K j j 1 (37)
= ′ ′ + - R t TA T [ ] E j E j j j , 1 1 (38)
= + -- TXB X T . j j , 1 1 1 (39)
Finally, the following optimization problem is stated:

E E ∑ ∑ + ′ ′ = - = min . j m k n j k j k X K K T , , , 1 1 1 (40)
with:

E E = ′ = ′ ′ ′ ′ ′ ′ ′ ′ ′ ′ ′ ′ ′ + + + ⊤ -⊤ - -⊤ - -⊤ ⊤ - + + ⊤ -⊤ - -⊤ - -⊤ ⊤ - + j k π π j k π π q K EK q K E K q K E K q q K E K q K E K q K E K q | | ‖ ( ) ‖ ‖ ( )‖ | | ‖ ( ) ‖ ‖ ( )‖ k j j k j j k j j k j k j j k j j k j j k j 1 1 2 1 2 1 1 2 1 1 2 1 2 1 1 2
The method, noted (M 2 ), minimizes Eq. ( 40) and exhibits sufficient conditions to estimate the thirteen parameters of ′ X K K , , , and T (six for ′ K K

, , one for the baseline T and six for the hand-eye transformation X).

Validation of extended self stereovision hand-eye calibration

The experiments were carried out on synthetic images to evaluate and compare the accuracy and the robustness of the methods (M1) and (M2) against noisy images.

Generation of noisy synthetic images

For six given displacements

= … j B ( ), 1, , 6 j 0
(given by Table 2) between two successives images from the stereorig, for two given matrices K and ′ K of internal parameters (defined by Eq. ( 22) and Table 3) and for an hand-eye transformation X (given by Table 3), a set of D

3 points Q = … i ( ) , 1, ,100 i i
is projected into the two images as ′ = … q q ( , ) i i i 1, , 100 using Eq. [START_REF] Rusu | 3D robotic mapping: the simultaneous localization and mapping problem with six degrees of freedom[END_REF]. At each coordinates pixel, a centered Gaussian noise with zero mean and standard deviation σ q from 0 to 2 pixels is added. In order to have statistical evidence, the results are averaged over 100 trials. The averaged resulting noisy points ′ = … q q ( , ) i i i 1, , 100 are used to estimate the performance of the two proposed methods.

Performances of methods M 1 and M 2

The minimization of Eqs. ( 32) and ( 40) was performed by a BFGS algorithm [START_REF] Lewis | Non-smooth optimization via quasi-newton methods[END_REF] with the initial guess vector composed of the initial intrinsic and baseline parameters and the initial parameters of the hand-eye transformation which are tabulated in Table 3.

The criteria for assessing the performance of the two methods are the relative error

∊ = - x x x ‖ ‖/‖ ‖ x r
and the absolute error

∊ = - x x ‖ ‖ x a
between true vector x and its estimation x . These criteria are calculated after each simulation run. The quantities ∊ a r or ∊ x a denote the mean value over all simulation runs. Fig. 5 gathers the mean value of the relative error of intrinsic parameters, ∊ r K , and the absolute error of translation parameter ∊ t a .

With synthetic images with little noise ( < σ 0.6 q pixels), the fit on intrinsic parameters with Method M 1 and Method M 2 are equivalents. In contrast, when the level of noise increases ( > σ 0.6 q pixels), ∊ r K increases with method M 1 whereas it remains lower than 5% with Method M 2 . The same behavior can be observed for the translation parameter t. The value of ∊ t a remains lower than 0.4 mm with Method M 1 and it rises above 1.5 mm with Method M 2 when σ q increases from 0.6 to 2 pixels.

The mean value of the absolute error on translation and rotation hand-eye parameters, ∊ a tX and ∊ a RX , shown in Fig. 6, confirms the previous trend. When the standard deviation of noise, σ q , is low, Method M 1 provides slightly lower values of ∊ a tX and ∊ a RX . When the standard deviation of the noise, σ q , is higher than 0.6 pixels, Method M 2 performs significantly better than Method M 1 and provides a maximal value for ∊ a tX (resp. ∊ a RX ) of 4.5 mm (resp. of °0.2 ). Finally, it might seem that the two methods are equivalent in a context with a low standard deviation of noise (σ q ) point extraction (less than 0.6 pixel). The Method M 2 is globally less sensitive to the noise level. In our application, the standard deviation of noise remains less than 0.6 pixel. The two methods can thus be used. However, the Table 3 Intrinsic parameters of both rectified cameras, K and ′ K , the baseline t and the parameters of Hand-Eye transformation, = X R t [ ] X X , used for the generation of synthetic images (S.) and for the initial guess (I.) of the minimization process. Fig. 5. Evolution of the mean value of ∊ r K , the relative error of intrinsic parameters of camera c, and ∊ t a , the relative error of translation parameter versus the standard deviation of noise for fits to Eqs. [START_REF] Laganière | Projective rectification of image triplets[END_REF] and [START_REF] Shiu | Calibration of wrist-mounted robotic sensors by solving homogeneous transform equations of the form AX=XB[END_REF].

α u ′ α u α v u 0 ′ u 0 v 0 t (-) ( -) ( -) (pix.) (pix.) (pix.) (mm)
advantage of Method M 1 is that it calculates the D 3 points conjointly with the calibration parameters.

Conclusion

The first objective function (Method M 1 ) is very efficient with a low level of noise. When the noise level increases, it is necessary to introduce a second objective function in the projective space which takes advantage of the epipolar constraint between two images acquired by the same camera from two robot positions. However, Method M 1 is able to calculate the D 3 points conjointly with the calibration parameters.

Automated thermal D 3 reconstruction

The automated thermal D 3 reconstruction will be implemented to monitor incremental forming processes [START_REF] Jeswiet | Asymmetric single point incremental forming of sheet metal[END_REF]. Its principle is to locally and gradually deform a metal sheet using a hemispherical tip tool (small diameter compared to the dimensions of the sheet) until the desired shape. The trajectory of the tip tool is controlled by a numerical control machine. This process makes it possible to form sheets of geometries complex with an extremely simple and therefore inexpensive tooling. The automated multiview thermal D 3 reconstruction method suggested aims at controlling both the local thermal gradient under the tip tool and guarantying the shape of the part during the forming. The section aims at testing and evaluating the shape measurement of a reference part, called pyramid, which is a test part of the incremental forming process.

The first section details the automatic NDT architecture composed by uncalibrated infrared cameras mounted on the Cartesian robot. The method, explained in the two previous sections, is summarized in the second subsection. The last subsection provides the results the D 3 re- construction of the reference part provided by the suggested thermal D 3 reconstruction method. It also compares this D 3 reconstruction with D 3 reconstruction performed by a reference D 3 Digitizer, Konica Minolta Range, with ± 50 µm accuracy.

Automatic NDT architecture for the automated thermal D 3 reconstruction

The automatic NDT architecture is a Cartesian robot equipped only with stereovision bench composed of uncalibrated infrared cameras which is displayed in Fig. 7(a).

The Cartesian robot positions the end-effector (hand) at the desired position for the acquisitions in a working volume of × × 1000 1000 800 mm. The position of the end-effector is known and given by the calculation of euclidean transformations (translation and rotation, named B ij in the previous section) between the world and the end-effector frames. A Cartesian robot offers a high rigidity (and so a good reproducibility of motion) and a good accuracy (thanks to the partial decoupling of the axes). The stereovision infrared cameras, installed at the end-effector, are compact and dedicated to embedded applications, such as drones. Their advantages are size and very low weight (around 400 g). These cameras are equipped with matrix sensors of uncooled microbolometers. They operate in 8-14 µm spectral band and work at 25 fps. Their resolution is × 160 120 pixels with a pitch of 25 µm. The Noise Equivalent Temperature Difference (NETD) is ap- proximately 100 mK. The focal length is 11 mm and acquisition distance is approximately 50 cm.

The automated multiview thermal D 3 reconstruction is tested on the pyramid part (see Fig. 7(b)), manufactured by an incremental forming process. Its base is about thirty centimeters with sides around ten centimeters. Its height is twenty centimeters.

Overview of the process of automated thermal D 3 reconstruction

As shown in Fig. 8, the first step of the automated thermal D 3 re- construction is the acquisition of multiviews from the robot displacement, so as to extract interest points on the part of the object to be reconstructed. Considering the displacements of the end effector of the robot expressed by the matrix B j and + B j 1 , at the position j, the k Fig. 6. Evolution of the mean value of the error of hand-eye translation parameters, ∊ a t X (high figure) and hand-eye rotation parameters, ∊ a R X , (low figure) versus the noise parameter ν for fits to Eqs. [START_REF] Laganière | Projective rectification of image triplets[END_REF] and [START_REF] Shiu | Calibration of wrist-mounted robotic sensors by solving homogeneous transform equations of the form AX=XB[END_REF]. matched points on each of the pairs of images and on each image between the displacement j and + ′ j q q 1, ( , ) k k j are obtained with a combination of detectors.

Thus, the matched points on each of the pairs of images are the inputs to calculate the pair of rectification homographies of the stereovision bench which minimize the projective deformations according to the method described in Section 2. The two obtained homographies H j and + H j 1 are used to compute the two rectified images I j and I′ j . Next, the end effector is moved to the position + j 1 using + B j 1 . At this new position, a new step of rectification is performed. Four rectified images I I I ′ + , , j j j 1 and I′ + j 1 are then available. Edge detection is then performed on this set in order to obtain n quadruplets

′ ′ + + = …
q q q q ( , , , )

k j k j k j k j k n 1 1 1
of matched points. This set of quadruplets is used in the self-calibration and reconstruction step presented in Section 3 and to solve Eq. [START_REF] Laganière | Projective rectification of image triplets[END_REF]. Solving this problem provides the calibration parameters ′ K K T , , , the hand-eye parameters X and two sets Q Q

= … + = … ([ ] ) , ([ ] ) k w j k n k w j k n 1 1 1
of D 3 -points. Note that the criterion of Eq. ( 32) is the sum of euclidean distance

Q Q - + ‖[ ] [ ] ‖ k w j k w j 1 2 .
Consequently, the two sets of obtained D 3 -points are very close at the end of the optimization process. Nevertheless, a choice must be made between the two sets. We define Q

= … ( [ ] ) k w j k n 1 , the final set of D 3 points, by the mean between Q = … ([ ] ) k w j k n 1 and Q + = … ([ ] ) k w j k n 1 1 .
This process leads to the D 3 partial reconstruction of the observed object summarized in Fig. 8.

Compared D 3 reconstruction on the reference part

To assess the performances of the suggested approach, a D 3 re- construction comparison was carried out on the pyramid part presented in Fig. 7. The reference D 3 reconstruction was performed using the D 3 Digitizer, Konica Minolta Range. It measures the shape of the part by profilometry with ± 50 µm accuracy for a measurement distance of 60 cm.

The reference D 3 reconstruction was compared to the D 3 points cloud performed by the suggested automated thermal D 3 reconstruction with = m 15 displacements of the robot. These fifteen displacements only provide a partial D 3 reconstruction of the part.

Rectification of uncalibrated cameras

As presented in the flowchart of Fig. 8, the first step consists in matching N quadruplets ′ ′ + + = … q q q q ( , , , )

k j k j k j k j k N 1 1 1
of points on the stereo images acquired for the position j of the end-effector. The Harris's detector [START_REF] Harris | A combined corner and edge detector[END_REF], coupled with ZNCC (Zero mean Normalized Cross-Correlation) method, provide the highest number of uniformly distributed matched points. Fig. 9a shows an example of the number of matched points ( = N 110). These matched points are the inputs of the rectification method suggested in Section 2. Fig. 9b displays the original and rectified images calculated with this set of matched points.

The criteria used to qualify the rectification image of Fig. 9 have already been presented in Section 2.3 and are tabulated in Table 4. The mean (respectively the standard deviation) on the rectification error E r is very low: 0.09 pixels (respectively 0.03 pixels). The maximal orthogonality error °-max E (90 ) o (respectively the maximal aspect ratio error max E (1 ) a ) is 0.7°(respectively 0.02). These results, with an orthogonality error of less than one degree highlight a very small deformation of the image and will ensure the introduction of a small deformation in the next calibration step.

Extended self stereovision hand-eye calibration

As presented in the flowchart of below. Next, fifty-seven quadruplets of matched points are extracted (i.e. = n 57) using phase congruency model [START_REF] Hajebi | Dense surface from infrared stereo[END_REF][START_REF] Hajebi | Structure from infrared stereo images[END_REF] associated with ZNCC method (see Fig. 10).

The selected criteria to analyze the calibration-reconstruction quality are: d , the mean value of reconstruction and e , the mean value of the epipolar distance. They are tabulated in Table 5. The mean value of reconstruction error d remains lower than 1.32 mm. The mean value of the epipolar distances e is higher than one pixel at around 3.25 pixels. The extended hand-eye self-calibration is then accurate and consistent with the spatial resolution of the infrared cameras.

The mesh performed from the D 3 point cloud calculated during this step is displayed in Fig. 11 where an image projection of the part is also plotted.

Comparison of D

3 point cloud calculated by the suggested method and provided with the laser scanner (ground truth)

Fig. 12 provides the D 3 deviations performed with Geomagic®soft- ware between the previous mesh and the reference model (ground truth) built with the laser scanner.

On this figure, the grey color represents the reference model. The over-printed color is the D 3 deviation values between the D 3 re- construction performed by the method and the reference model. This colorbar provides the accuracy of the reconstruction for the different areas of the part. The left box indicates the colorbar associated to each D 3 deviations and displays the histogram. The two images differs from the range of colorbar. For the first image, the range of colorbar is calculated using the maximal deviation. The image is called full-scale. The second image is displayed with a range containing the most of the deviation values. This image is named restricted-scale. These different scales make full use of the histogram of D 3 deviations, in mm. The first image with a full-scale highlights the maximum values of the reconstruction deviation. The maximal deviation value reaches 12.8 mm and is located on the left field of the part. Most of D 3 devia- tions values are in the interval -[ 1.5 mm, 1.5 mm]. Thus, the part is mainly green. In the second image, the colorbar is then restricted to the interval -[ 1.5 mm, 1.5 mm] to analyze accurately deviation locations. The mean error on the top of the part is then equal to 1.7 mm with an standard deviation of 1.8 mm. The maximal deviation is recorded on the edges of the part where the area is not textured enough.

These deviations are comparable to the spatial resolution of the camera which is 1 mm (observation distance 50 mm, focal length 11 mm, pixel size 51 µm and a half correlation window of 4.5 pixel). This reconstruction performance on a plane with very low resolution cameras validates the D 3 automated thermal reconstruction with a multi- view approach performed with the help of a robot equipped with two uncalibrated infrared cameras.

Conclusions and future works

This article presented algorithms required for performing an automated thermal D 3 reconstruction of a part with a system composed of a stereovision rig of uncalibrated infrared cameras mounted on a six-axis Cartesian robot. This multi-view D 3 thermal reconstruction relies on a geometric extended self-calibration for estimating both intrinsic parameters for each camera, relative position and orientation between the two cameras and the Euclidean transformation between the robot and the camera reference frames.

A new projective rectification method was first introduced on one hand for improving the stereo matching problem between pixels on the left and right images reducing the search space from two dimensions to a single one, and on the other hand for decreasing the number of parameters of the cost function of the calibration problem. Because the suggested method must cope with noisy and low-textured infrared images, the two homographies needed for the rectification were calculated in a single step, from only a few matched points and without previous information on the Fundamental Matrix. The suggested rectification keeps the aspect ratio and the orthogonality of images, and consequently, limits the creation or loss of pixels. The structure and the skewness of the images are preserved with an error of less than 0.8% on both criteria. Before demonstrating its effectiveness on infrared images in the last section, this method was benchmarked against recent rectification methods on images classically used in the vision community.

Next, the set of matched points from rectified images was used to perform simultaneously, also in a single step, the extended self-calibration and the multi-view infrared D 3 reconstruction. The minimiza- tion problem of calibration was solved through two cost functions selected according to the noise level measured in infrared images. The first objective function aims at minimizing the reprojection errors (in pixels). When the noise level increases, a second objective function is expressed in the projective space, from the epipolar constraints between two images acquired by the same camera from two robot positions. These two functions make explicit the intrinsic camera parameters and the essential matrix between the robot positions, which itself depends on the hand-eye parameters and the robot motion. The tests with simulated images, with noise on the extracted points, confirmed that the method based on the reprojection error is suitable, with a noise less than 1 pixel. Above this value, the second method is more suitable. The reconstruction based on a rectified stereovision rig gives low registration errors with an approach in a single step. Finally, the automated thermal D 3 reconstruction was validated by comparing its results with those obtained with a D 3 scanner with an accuracy to 50 µm. An accurate D 3 model of a hot object with an ac- curacy around ± 1 mm was achieved with several views acquired by uncalibrated infrared cameras.

Fig. 1 .

 1 Fig. 1. Flowchart of the suggested method.

44

  be the corners of the images I and I′. In the rectified space, let ̃̃̃p p p p { , be the corners of rectified images. The aspect ratio E a and the aspect ratio constraint C a are then defined by:

Fig. 2 .

 2 Fig. 2. Criteria summary. Image (a) shows deformations due to an aspect ratio different to 1. Image (b) displays deformations due to an orthogonality criterion different to °90 . Image (c) displays deformations due to size ratio different to 1.

3. 3 .

 3 Definition of the two objective functions 3.3.1. Objective function from the geometry of the stereo image sequencesGiven a set of D3 points Q = ( ) in R w .It is assuming that all the points are seen by a stereo rig displaced using m rigid transformation = . For an easy reading of this section, we denote =

Fig. 7 .

 7 Fig. 7. (a) Thermal robotic set-up composed of a Cartesian robot and infrared cameras. (b) Reference part, called pyramid, used for testing the automated multiview thermal D 3 reconstruction.

Fig. 8 , 1 Fig. 8 .

 818 Fig. 8. Flowchart of the automated thermal D 3 reconstruction.

Fig. 9 .

 9 Fig. 9. Example of stereo infrared images of the pyramid part. (a) 110 Matched points between the stereo infrared images with Harris detector and ZNCC method. (b) Original (on the top) and rectified infrared images I I′ , j j (on the bottom).

Fig. 10 .

 10 Fig. 10. Example of = n 57 matched points on four rectified images.

Fig. 11 .

 11 Fig. 11. Mesh from the D 3 points cloud with a projection of the image part.

Fig. 12 .

 12 Fig. 12. D 3 deviations at different scales (full and restricted scales) between the reference model and the mesh emerging from the D 3 points cloud calculated by the suggested method. The box indicates the scale of each image in mm and the values of the displayed histogram.

  .

	infinity. Hence, H and	′ H transform the matched points
	0
	This process is accomplished by applying two homographies I and I′. These homographies map the epipoles e and ′ e to points at ′ H H on ,

Table 1

 1 Comparison of the suggested rectification method with Monasse's method of the article[START_REF] Monasse | Three-step image rectification[END_REF].

	Sample	Method	Orthogonality E o		Aspect Ratio E a		Rectification error E r
			′ H (°)	H (°)	H ′ (no unit)	H (no unit)	mean (pix.)	σ (pix.)
	Boxes	Proposed	89.44	90.66	0.9897	1.0105	0.11	0.06
		Monasse	89.60	89.63	0.9884	0.9892	0.1293	0.0887
	Arch	Proposed	89.67	89.88	0.9955	0.9983	0.17	0.16
		Monasse	89.80	90.05	0.9942	1.0014	0.2520	0.2349
	Drive	Proposed	90.61	89.32	1.0128	0.9876	0.55	0.44
		Monasse	89.95	90.00	0.9977	1.0001	0.7139	0.8253

Table 2

 2 Six displacements of the end-effector transformation,

	the generation of synthetic images.		B	j 0	=	R t [ B 0j B 0j ] , used for
	Displ.		rotation (°)		translation (mm)
	j	r j x B 0	r j y B 0	r j z B 0	t j x B 0			t j y B 0	t j z B 0
	1	0	180	0	0			0	500
	2	45	180	14	10			200	356
	3	-43	-167	22	-318			-158	424
	4	49	175	4	12			189	832
	5	60	167	34	-321			376	294
	6	-34	-163	-45	531			-269	756

Table 4

 4 Rectification performances on the pyramid part.

	Orthogonality E o	Aspect Ratio E a	Rectification error E r
	′ H (°)	H (°)	H ′ (no unit)	H (no unit)	mean (pix.)	σ (pix.)
	89.29	90.07	0.9882	0.9991	0.09	0.06

Table 5

 5 Results on extended hand-eye self-calibration with the reference part.

	Σ d	d	σ d	Σ e	e	σ e
	(mm)	(mm)	(mm)	(pix.)	(pix.)	(pix.)
	37.34	1.32	0.86	185.74	3.25	7.96

Note that to improve the readability of Sections

and

3, sans-serif font upper-case is used for 3D points (e.g. Q) while bold lower-case is used to denote D 2 points (e.g. q).