Rémi Boutteau
email: boutteau@esigelec.fr

Romain Rossi
email: rossi@esigelec.fr

Lei Qin
email: qin@esigelec.fr

Pierre Merriaux
email: merriaux@esigelec.fr

Xavier Savatier
email: savatier@esigelec.fr

A vision-based system for robot localization in large industrial environments

Keywords: vision-based localization, fisheye cameras, robust pose estimation, extrinsic calibration 1

In this paper, we propose a vision-based system to localize mobile robots in large industrial environments. Our contributions rely on the use of fisheye cameras to have a large field of view and the associated algorithms. We propose several calibration methods and evaluate them with a ground-truth obtained by a motion capture system. In these experiments, we also evaluate the influence of the parameters as the number of points used for calibration, or the influence of the accuracy of these points. Our system is then experimented in a real industrial environment, where we localize blue an ROV (Remotely Operated underwater Vehicle) in a basin. As shown is this paper, this system can also localize wheeled mobile robots in the same way.

4.0, they evolve to be more autonomous, interact with each other and collaborate with humans safely. The new technologies used in the Factory of the Future have brought out the concept of Cyber-Physical System (CPS), a system that interacts with its environment, extracts data, and processes them to control a physical process. In this context, Cyber-Physical Systems deployed in companies must be autonomous and embedded.

One of the main applications of mobile robotics in factories is logistics, through the use of AGV (Automatic Guided Vehicle), for example to move products or raw materials between different workplaces or storage spaces. These Cyber-Physical Systems can also perform inventory tasks, for example by moving in the factory while counting and locating objects with RFID tags [START_REF] Zhang | Mobile robot for retail inventory using rfid[END_REF].

To be able to navigate autonomously, a mobile robot must be able to localize itself accurately in its environment. In this purpose, several types of sensors can be used. The goal being to localize mobile robots in indoor environments, we can discard GPS-type technologies whose signals are not received in buildings. If we exclude proprioceptive sensors (e.g. odometers) that are generally not accurate enough for a reliable localization, two technologies are generally used: lidars and cameras.

Lidar (Light detection and ranging) sensors have the advantage of providing accurate and reliable data. Many localization systems work with these sensors [START_REF] Merriaux | Robust robot localization in a complex oil and gas industrial environment[END_REF][START_REF] Nobili | 16 channels velodyne versus planar lidars based perception system for large scale 2d-slam[END_REF], but there are some limitations. The first is its high price, especially if it is necessary to use 3D lidars. This is specially true if a whole fleet of robots must be localized (each robot will have to embed one or more lidars). The second limitation is their electrical consumption since lidars are active sensors.

Vision sensors is the second technology mainly used for localization and has undeniable advantages: the cameras are passive sensors, and are inexpensive. Many approaches have been developed in recent years for robot or vehicle localization based on vision, the most advanced solutions being [START_REF] Mur-Artal | Orb-slam2: an open-source slam system for monocular, stereo and rgb-d cameras[END_REF] and [START_REF] Wang | Stereo dso: Large-scale direct sparse visual odometry with stereo cameras[END_REF].

In this work, we are interested in a vision-based system for robot localization by instrumenting the factory rather than the robot itself. This makes it possible to locate a whole fleet of mobile robots without multiplying the number of cameras required. The application described in this paper is the remote control of an ROV (Remotely Operated underwater Vehicle) to recover industrial sludge in a silo, but this system can easily be used for the localization of AGVs in a factory. Our experimental results also illustrate this use case since we use wheeled mobile robots.

Our main contributions are:

the use of omnidirectional vision theory, to develop a system for robots localization in large environments using a reduced number of cameras, the development of several calibration methods, depending on the constraints of the environment, the evaluation of these calibration methods and the localization accuracy by studying the influence of several parameters, the use of our localization system in a real industrial context.

The remainder of the paper is organized as follows. Section 2 presents the overall architecture of the system we have developed. Section 3 describes the Unified Projection Model (UPM) that we used to model fisheye cameras. Section 4 is dedicated to the various calibration methods that we have developed according to the operational constraints of the application (prior knowledge of the building, etc.). Section 5 presents our method of detection and localization. Section 6 presents our experimental results, both quantitative in the laboratory with a ground truth and qualitative in our industrial application. We present the results with the different calibration methods and we study the influence of the different parameters on the localization accuracy. Finally, Section 7 summarizes this paper and includes some ideas for future work.

System overview

Our localization system is based on the tracking of a light pattern embedded on the mobile. In our work, we used a luminous pattern composed of 3 LEDs forming a "T" and spaced 300mm apart on one axis and 200mm apart on the second axis as shown in the Figure 1. This system can be embedded on a mobile robot or an AGV as illustrated by our experimental results and in Figure 8.

The industrial application of our work is localizing an ROV in a silo. This aim of our ROV is to collect industrial sludge in a silo using a delamination screw. To empty the silo, it is therefore necessary to make adjacent passages to remove the sludge to a depth of a few centimeters. To be able to do this, the operator must know accurately the position of the ROV because it does not have a visual feedback, the ROV being in a closed silo and therefore very dark.

The complete system we have developed for the ROV localization is described in Figure 2. One or two fisheye cameras (depending on the size of the silo) are positioned on the roof of the silo and observe the surface of the water. The ROV lies on the sludge and is totally immersed in the water. To be able to locate it, it is equipped with a buoy having three LEDs as shown in Figure 3. This buoy is attached to an articulated rod equipped with a sensor allowing to obtain its angle of inclination with respect to the ROV. A rangefinder provides us the distance between the roof of the silo and the surface of the water since the water level varies during the sludge extraction process.

The general idea is to detect the LED points in the image captured by the fish-eye camera and then to reconstruct their 3D positions. Concretely, the proposed system consists of the following components: calibration (intrinsic parameters, extrinsic parameters), LED detection and 3D reconstruction, localization and finally multiple camera fusion.

The Unified Projection Model (UPM)

Fisheye cameras present a very great interest in our application because we can cover a very large area with only one or two cameras. However, fisheye cameras do not obey the pinhole projection model commonly used in computer vision [START_REF] Hartley | Multiple View Geometry in Computer Vision[END_REF].

Fig. 2 System overview. The ROV lies on the sludge and is totally immersed in the water. To be able to locate it, it is equipped with a buoy having three LEDs. This buoy is attached to an articulated rod equipped with a sensor allowing to obtain its angle of inclination with respect to the ROV. A rangefinder provides us the distance between the roof of the silo and the surface of the water since the water level varies during the sludge extraction process.

Fig. 3 The ROV in a test basin similar to the silo. On the left when it is on the surface, on the right in the working situation.

Many models have been proposed to represent a fisheye camera. These models can be classified into 3 categories.

The first category of models directly uses the optical models that were used by the manufacturer to build the optics. The projections commonly used are then stereographic, equiangular, sine-law or equi-solid depending on the reference of the optics [START_REF] Sturm | Camera models and fundamental concepts used in geometric computer vision[END_REF].

The second category is based on the use of the pinhole model commonly used in computer vision. A 3D point is projected on the image plane and then distortions are added to the projected point to take into account the deformations induced by the fisheye. The mapping between the perspective projection and the fisheye point can be logarithmic [START_REF] Basu | Alternative models for fish-eye lenses[END_REF], polynomial [START_REF] Shah | Intrinsic parameter calibration procedure for a (high-distortion) fish-eye lens camera with distortion model and accuracy estimation[END_REF], or rational [START_REF] Brauer-Burchardt | A new algorithm to correct fish-eye-and strong wide-angle-lensdistortion from single images[END_REF][START_REF] Fitzgibbon | Simultaneous linear estimation of multiple view geometry and lens distortion[END_REF].

Finally, the third category is based on the use of the Unified Projection Model (UPM) initially used for catadioptric sensors [START_REF] Barreto | A unifying geometric representation for central projection systems[END_REF][START_REF] Geyer | A unifying theory for central panoramic systems and practical implications[END_REF][START_REF] Mei | Single view point omnidirectional camera calibration from planar grids[END_REF]. It has been demonstrated on numerous occasions that this model is also adapted to fisheye cameras [START_REF] Courbon | A generic fisheye camera model for robotic applications[END_REF][START_REF] Courbon | Evaluation of the unified model of the sphere for fisheye cameras in robotic applications[END_REF][START_REF] Ying | Can we consider central catadioptric cameras and fisheye cameras within a unified imaging model[END_REF].

The model used in this paper belongs to this third category. We used the model introduced by Mei and Rives [START_REF] Mei | Single view point omnidirectional camera calibration from planar grids[END_REF] and inspired from those of Barreto [START_REF] Barreto | A unifying geometric representation for central projection systems[END_REF] and Geyer [START_REF] Geyer | A unifying theory for central panoramic systems and practical implications[END_REF]. In this section, the equations for projecting a 3D point to obtain the coordinates of the corresponding pixel are presented. The inverse projection, called lifting, allowing to pass from one pixel of the image to the direction of the corresponding light ray, is also presented. In general, points are expressed in a reference frame different from the one of the camera, called world frame. The two reference frames are related by a rigid transform T cam world . If X world = [X Y Z 1] is the homogeneous coordinate vector of a point expressed in the world frame, and X cam Is the coordinate vector of this same point expressed in the camera frame, then we can write

X cam = T cam world X world . (1)
If C represents the coordinates of the camera center in the world coordinate system, and R is the matrix expressing the orientation of the camera frame, then Equation (1) can be written

X cam = R -RC 0 1×3 1 X world . (2)
There are several representations to express a rotation. In our model, we use quaternions. The rigid transform function, denoted by w, expressed in Equation (1), will be characterized by a vector of 7 elements (4 for the quaternion and 3 for the translation): V 1 = [q w q x q y q z t x t y t z ] .

Once the coordinates of the point are expressed in the camera frame X cam = [X c Y c Z c 1] , the latter is projected onto the unit sphere using the following equation:

X S =   X S Y S Z S   =     X c √ X 2 c +Y 2 c +Z 2 c Y c √ X 2 c +Y 2 c +Z 2 c Z c √ X 2 c +Y 2 c +Z 2 c     . (3)
This point undergoes a second projection on the normalized plane from a point situated at a distance ξ from the sphere center. The projection of X S onto the normalized plane, denoted by m = [x y] , is therefore obtained by: m

= x y = X S Z S +ξ Y S Z S +ξ . (4)
These two transformations (projection on the sphere and then projection onto the normalized plane) are grouped together in the h function which depends only on a parameter:

V 2 = [ξ ].
Distortions are then added to the point m by applying the distortion function d which involves five coefficients: three radial distortion coefficients (k 1 , k 2 and k 5 ) and two tangential distortion coefficients (k 3 and k 4 ). Let V 3 = [k 1 k 2 k 3 k 4 k 5 ] be the parameter vector containing these coefficients, ρ = x 2 + y 2 and m d = [x d y d ] the resulting point, its coordinates are obtained by:

m d = m + d(m), (5)
and consequently:

m d = x(1 + k 1 ρ 2 + k 2 ρ 4 + k 5 ρ 6 ) + 2k 3 xy + k 4 (ρ 2 + 2x 2 ) y(1 + k 1 ρ 2 + k 2 ρ 4 + k 5 ρ 6 ) + 2k 4 xy + k 3 (ρ 2 + 2y 2 ) . (6)
The final projection is obtained by the generalized projection matrix K. The latter contains 5 parameters: the generalized focal lengths γ u and γ v , the coordinates of the principal point u 0 and v 0 , and the skewness parameter α. Let k be this projection function and V 4 = [α γ u γ v u 0 v 0 ] the vector of parameters on which it depends, the projection p = [u v 1] of the point m d on the image plane is then given by:

p = K m d 1 , (7) p
=   u v 1   =   γ u γ u α u 0 0 γ v v 0 0 0 1   m d 1 . (8)
It is common to simplify the matrix K by setting α = 0, which corresponds to assuming that rows and columns of the pixels are orthogonal.

Let V = [V 1 V 2 V 3 V 4 ]
be the vector containing the 18 parameters of the model. The global projection function of a point X, denoted by p(V, X) is obtained by composition of the different functions presented previously:

p(V, X) = k • d • h • w(V, X). (9)
These steps are applied to obtain the projection onto the image plane of a 3D point knowing its coordinates in the 3D space.

The inverse projection, called lifting, consists in retro-projecting a point, i.e. knowing its pixel coordinates p pix = [u v 1] , determining the direction of the corresponding light beam.

It is necessary to start by re-projecting the p pix point on the normalized plane by performing an inverse perspective projection to obtain the point m = [x y 1] :

m = K -1 p pix , (10)
hence:

m =    u-u 0 γ u v-v 0 γ v 1    . (11)
To obtain the corresponding X s point on the sphere, we use the formula:

X S =      ξ + √ 1+(1-ξ 2 )(x 2 +y 2 ) x 2 +y 2 +1 x ξ + √ 1+(1-ξ 2 )(x 2 +y 2 ) x 2 +y 2 +1 y ξ + √ 1+(1-ξ 2 )(x 2 +y 2 ) x 2 +y 2 +1 -ξ      . ( 12
)
4 Calibration

Two kinds of calibration have to be taken into account in such an application. First, we have to estimate the parameters V 2 , V 3 and V 4 of the unified projection model. These parameters are inherent properties of a camera and this kind of calibration is called intrinsic calibration.

To do so, we used a calibration method described in [START_REF] Boutteau | An omnidirectional stereoscopic system for mobile robot navigation[END_REF][START_REF] Boutteau | Mobile Robots Navigation, chap. Chapter 1 : A 3D Omnidirectional Sensor For Mobile Robot Applications[END_REF]. Once V 2 , V 3 ,V 4 are estimated, we are about to calibrate V 1 . Indeed, V 1 encodes the geometric relation between the camera coordinate frame and the world coordinate frame. In this sense, the process for calibrating V 1 is usually called extrinsic calibration.

For notation convenience, we denote hereafter the position of a LED point in the world coordinate frame as X led world , the position of the camera in the world coordinate frame as X cam world , the position of the camera in the ROV coordinate frame as X cam rov and the pixel coordinates of a LED point in the image as p led pix . Similarly, we denote the orientation of the camera in the world coordinate frame as R cam world and the orientation of the camera in the ROV coordinate frame as R cam rov . The extrinsic parameters to be estimated, i.e. V 1 , can be decomposed into R cam world and T cam world . After calibration, to reconstruct the 3D position of a point by its projection in the image, we also need the equation of the plane π in addition to the projection parameters V. Depending on whether we have a set of precisely known 3D positions in the scene, we present in the following two methods for estimating R cam world , T cam world and π.

Extrinsic calibration using points with known 3D positions

Estimating R cam world and T cam world requires a set of 3D-2D point correspondences. If we have a set of N points whose precise 3D positions in the world coordinate frame are known, we can take a picture of the scene and get the pixel locations of these points in the image. In this way, for the ith point in the set, both its 3D position in the world coordinate frame X i world and the 2D pixel location of its projection in the image p i pix are available. The estimation of R cam world and T cam world consists in finding the parameters that minimize the sum of re-projection errors:

E(V 1 ) = 1 2 N ∑ i (p i pix -P(V 1 , X i world )) 2 . ( 13
)
The parameters to be estimated have six degrees of freedom (three for the translation and three for the rotation). As each correspondence leads to two equations (x-axis and y-axis), N must be equal to (or greater than) 3. In practice, at least one supplementary point correspondence is needed in order to assess the quality of the optimization result. We use the Levenberg-Marquardt algorithm [START_REF] Levenberg | A method for the solution of certain problems in least squares[END_REF] to solve this nonlinear least-squares optimization problem (Fig. 5). We can ask ourselves the impact of the accuracy of the 3D points coordinates (which are not always easy to determine) on the result of the extrinsic calibration. Similarly, the minimum number of points is 3, but does having more points have an impact on accuracy? All these questions are discussed in section 6.

Auto-Calibration by moving the ROV

In real environments, it can be difficult to measure 3D points very accurately. In our case, for example, it is difficult, if not impossible, to access the interior of the silo. Consequently, it may be difficult to obtain a set of precisely known 3D points for calibration. To circumvent this problem, we propose in this section an alternative calibration method by moving the ROV in the scene.

Thanks to the camera intrinsic parameters and the knowledge of the LEDs geometry, we can estimate the position of the ROV with respect to the camera. The method is exactly the same as the one presented in section 4.1 except that the known points are no longer those of the silo but the ones of the ROV LEDs. This problem is similar to what is known in literature as the PnP (Perspective-n-Point) problem [START_REF] Lepetit | Epnp: An accurate o(n) solution to the pnp problem[END_REF]. The methods proposed in the literature are not directly applicable, however, since they are intended for conventional perspective cameras; in our case, we minimize the function given in Equation ( 13) to recover the pose of the ROV with respect to the camera.

We therefore have a set of ROV positions. It is consequently possible to estimate the plane parameters on which the ROV moves in the camera frame: its normal n cam = [n x n y n z ] and the orthogonal distance d between the camera and the plane. We know that in the silo frame, this normal is vertical since the target floats on the water. We consequently have: n world = [0 0 1] . To obtain the orientation of the camera, we look for the rotation matrix R world cam which allows us to switch from the normal expressed in the camera frame to the one expressed in the silo frame:

  n x n y n z   =   r 11 r 12 r 13 r 21 r 22 r 23 r 31 r 32 r 33     0 0 1   . ( 14
)
Rotation around the water plane's normal (yaw) cannot be recovered. We thus only consider pitch and roll angles. The Denavit-Hartenberg [START_REF] Boutteau | Circular laser/camera-based attitude and altitude estimation: minimal and robust solutions[END_REF][START_REF] Hartenberg | A kinematic notation for lower pair mechanisms based on matrices[END_REF] parameterization of R world cam with these two angles leads to

  n x n y n z   =   cos(θ ) -sin(θ )cos(α) sin(θ )sin(α) sin(θ ) cos(θ )cos(α) -cos(θ )sin(α) 0 sin(α) cos(α)     0 0 1   . (15)
From [START_REF] Kalman | A new approach to linear filtering and prediction problems[END_REF], α (roll) and θ (pitch) can be easily recovered since

   α = arccos(n z ) θ = arcsin(n x ) √ 1-cos 2 (α) . ( 16
)
With this method, it is therefore possible to find the attitude (roll and pitch) and the altitude of the camera relative to the plane of the water. To find an absolute pose, it is possible to use one of the ROV poses (for example the first one) as the origin of the world frame. The ROV can for example be positioned at a corner of the silo to set the origin.

Detection and Localization

In this section we present the processing chain that allows us to locate the ROV in the silo from the images. The points corresponding to the LEDs are detected in the images, then their 3D coordinates are calculated. A series of tests is carried out to reject possible false detection. Finally, the location of the ROV can be estimated. At the end of this section we also present the method used when the environment is large and requires the use of several cameras.

LEDs detection in the image

In order to get p led pix , we need to detect the projections of the LEDs in the image. Here we do not use the background subtraction methods as in Section 4.2 due to speed considerations. Instead, we use a series of image processing techniques to detect the LEDs of the ROV in the image.

Given an image captured by the camera with the LEDs of the ROV in the scene, we first blur it to remove singular points. Then, we convert the image to binary by an intensity thresholding. Dark regions are rejected. Next, we apply a connected components labeling in the binary image to get a number of separated areas. In fact, each area represent a candidate LED point for further validation. We can reject the areas that are either too small or too large if we have some prior knowledge about the size of the LEDs in the image. Finally, we use the center location of the area to represent the pixel coordinates p pix of a candidate point.

It is often the case that after this series of processing there are still a large number of candidate points, among which we need to select the real LED points. To tackle this problem, we propose to first reconstruct the 3D positions of the points and then to detect the LEDs in the reconstructed plane according to their geometrical configurations.

3D Position Reconstruction

We now present the method to reconstruct the 3D position X world of a point X, which is supposed to lie on the plane π, from its projection p pix in the image.

First we lift p pix back to the unit sphere to obtain X s with equations ( 11) and ( 12). This point is then expressed in the world frame using the rigid transform T world cam obtained by the extrinsic calibration step. Let A be this point.

Obtaining the 3D coordinates of the LED point consists in computing the intersection of the line that connects X s and the unit sphere center (see Figure 4) with the plane π. Let B be the position of the unit sphere center in the world coordinate frame.

The line that joins A and B can be expressed by a Plücker matrix [START_REF] Hartley | Multiple View Geometry in Computer Vision[END_REF] as:

L = AB -BA . (17)
Finally, the reconstructed 3D position is the intersection of the line L and the plane π. We get X world = L • π. Note that the obtained X world is represented in homogeneous coordinates. We can get the original coordinates back by dividing X world with its last element.

LEDs detection in the Reconstructed Plane

Let S be the set of the candidate points detected in the image and N be the number of LEDs of the ROV. If size(S) is greater than N, we first reconstruct the 3D positions of the points in S using the method described in Section 5.2. Then, we reject the points that are out of the region of interest. In fact, we know a priori the area of the silo in the world coordinate frame. As the ROV can only appear in the silo, we reject from S all the points whose 3D positions are not within the silo area.

Then, we cluster the remaining 3D points into local groups by hierarchical clustering [START_REF] Maimon | Data Mining and Knowledge Discovery Handbook[END_REF]. The distance threshold used for clustering is determined adaptively according to the geometrical configuration of the LEDs. Specifically, the longest distance between two LEDs in our system is 316 mm (LED 2 to LED 1 or LED 2 to LED 3 in Figure 1). We use a threshold of 500 mm (slightly greater than 316 mm ) for clustering. After clustering, we remove from S the groups that have either too few or too many points. The groups with too few points consist of isolated light spots while the groups with too many points are reflections of light on the surface of the water.

Finally, we recognize the real set of the LEDs from the retained groups by verifying the internal geometrical characteristics, namely the distances between the points, the shape of the contour and the area. The distance constraint indicates that a plausible candidate should have some "support" points in S that are within the predefined distance tolerance range I dist from it. For example, in our ROV configuration as depicted in Figure 1, a valid candidate for LED 0 should have two support points at about 100 mm and a third support point at about 300 mm. Similarly, a valid candidate for the LED 2 should have two support points at about 316 mm and a third support point at about 300 mm. The points that have no valid supports from other points in S are removed. This process is repeated and S shrinks. After stabilization, if Fig. 6 The selection/rejection cascade to detect the LED points. size(S) is still greater than N, we decrease I dist by a small amount and restart the shrinkage iteration until I dist reaches the predefined minimum I dist min or size(S) is less than or equal to N.

For the shape constraint, we first build a reference contour representing the geometrical layout of the LEDs in the ROV coordinate frame. Then for each local group, the points in the group form a candidate contour. We employ the Hu Moments [START_REF] Hu | Visual pattern recognition by moment invariants[END_REF] for comparing the shape of a candidate contour with the shape of the reference contour. Contours with comparing scores exceeding the predefined tolerance range, I shape , are rejected. Their corresponding 2D points are removed from S. Similarly, for the area constraint we reject the candidates whose areas are out of the tolerance range, I area , from the reference area. These three kinds of geometrical verification are coherent and complementary. Combining them together can help increasing robustness for recognizing the LEDs. We summarize the point selection process in Figure 6.

ROV Localization

The points remaining in S have successfully survived the rejection/selection cascade. If size(S) is equal to N, we claim that the target is found. We can determine the order of the points by distance verification. The 3D position of the ROV can then be computed using the method described in Section 5.2.

Conversely, if size(S) is either greater or less than N, we claim that the target is not found in the test image. Further developments can use the motion history of the ROV to help distinguishing the ROV from similar objects when size(S) is greater than N.

Scalability: Multiple Cameras Fusion for Large Areas

In large areas, multiple cameras are needed to cover the entire terrain or to improve the localization accuracy. The purpose of fusion is twofold. First, the objective is to be able to estimate the pose of the ROV as soon as it is observable by one of the cameras. Second, when the ROV is observable by several cameras, the information is merged to optimize the estimation of its pose. The objective is therefore to obtain a more complete trajectory while avoiding jumps in the estimation of the ROV pose when switching from one camera to another. We treat this requirement as a state estimation and data fusion problem.

Several methods can be used for data fusion applied to localization: histogram filters, Kalman filters and particle filters [START_REF] Thrun | Probabilistic robotics[END_REF]. Histogram filters only allow a discrete estimation of the state space, which can introduce a loss of precision depending on the chosen discretization step. In addition, their complexity increases exponentially with respect to the sampling step and size of the environment. Kalman filters have the advantage of estimating a continuous state space and are relatively easy to implement. Their main disadvantages are that they do not allow the estimation of a multimodal distribution and that their efficiency is quadratic with respect to the state space to be estimated. When it is necessary to estimate both the pose of a robot and the positions of the landmarks, this can become problematic and we then turn to the particle filter. Particle filters allow the estimation of a continuous state space and a multimodal distribution. In our case, we want to merge the position from the two cameras without having to estimate the map (this is a localization problem and not a SLAM problem), so the Kalman filter is the best compromise between ease of implementation and efficiency. In addition, it has been shown in some studies, for example in [START_REF] Gutmann | An experimental comparison of localization methods[END_REF], that the particle filter is more robust but that the Kalman filter can be more accurate.

Since each camera measures the state of the ROV independently, we use the Kalman filter [START_REF] Kalman | A new approach to linear filtering and prediction problems[END_REF] to fuse the localization results of each camera and to estimate the state of the target. Specifically, at each time step, the Kalman filter predicts the state and each camera provides measurements to update this state. The state that gives the maximum a posteriori probability is used as the state of the target.

Experimental Evaluation

In this section, we evaluate the performance of the proposed localization system using different extrinsic calibration methods and different hypotheses on the input data. Quantitatively, we compute the position errors and the orientation errors with respect to the "ground truth" provided by a motion capture system (Vicon).

Experimental Settings

In order to have a practical evaluation of our algorithms, a dataset has been collected with a reliable ground-truth obtained by a motion capture system. The experiments have been conducted in a room equipped with a Vicon motion capture system composed of 20 Vicon T40S cameras. With such a system, we can assure a 6 DoF (Degrees of Freedom) localization of our system with a sub-millimetric accuracy as demonstrated in [START_REF] Manecy | X4-mag: A low-cost open-source micro-quadrotor and its linux-based controller[END_REF] and [START_REF] Merriaux | A study of vicon system positioning performance[END_REF] and at a high framerate (500fps). We implement the localization system in C++ with the OpenCV library. The validation environment in our laboratory is presented in Figure 7. The two blue cameras on top (about 3.60 meters from the ground) are the fish-eye cameras used for localization. The cameras on the wall are those of the Vicon system that provides the "ground truth" localization. The ROV with its LEDs on is emulated by a remotely operated Wifibot robot on the ground. A closer view of the robot with its 4 LEDs used in our laboratory for algorithm validation is shown in Figure 8. As a reminder, the real ROV in its environment is shown in Figure 3.

Our system takes less than 70 milliseconds for processing one frame when executing on a PC with an Intel Core i5-2520M CPU running at 2.50GHz and 8GB of memory. Firstly, we evaluate the calibration method using points with known 3D positions. We place some markers in the room and use the Vicon to get their 3D positions. To understand the influence of the number of points, we tested three configurations: one with 4 points, one with ten points and the last one with with 18 points. The resulting localization errors with respect to the results provided by the Vicon are summarized in Table 1. The estimated and real trajectories for the 4 and 18 points configurations are displayed in Figure 9. We see that as long as the 3D positions of the points are accurate and are geometrically scattered, using more points for calibration does not necessarily improve localization results. On the contrary, noise can be introduced using more points when calculating their pixel locations in the image, especially when the points are near the borders of the image.

Injecting errors to the 3D positions used for calibration

In order to assess the robustness of the localization system with respect to the accuracy of the 3D positions of the points used for calibration, we introduce errors to the positions of the points and recompute the position and orientation errors. Specifically, 50mm, 100mm and 200mm errors are injected to the position of each of the 4 points. The resulting localization errors are summarized in Table 2 and the estimated trajectories are displayed in Figure 10. We see a roughly linear relation between the introduced error and the resulting localization error. 3 and the trajectories are displayed in Figure 11. Note that for the "1 pose" case, we assume the equation of the plane to be horizontal, i.e its normal vector is assumed to be [0 0 1] , as we do not have enough points to estimate it. For the other configurations, π is estimated by the RANSAC procedure using the poses of the ROV.

Fusing the localization results of multiple cameras

Finally, we present the localization results obtained by fusing the results of multiple cameras and compare them with those obtained by a single camera. As shown in Figure 7, we have two fisheye cameras that can perform localization independently at the same time. We have recorded a sequence for evaluation. The localized trajectories are presented in Figure 12, where the advantages of the fusion algorithm can be observed qualitatively.

As we can see, the ROV trajectories estimated by camera 1 or camera 2 are not complete. Indeed, there may be a occultation of the LEDs at certain times, or the system may not be able to estimate the position of the ROV (for example if it is too far from the camera). In the real situation, i.e. for the ROV in its silo as shown in Figure 3, the tether of the ROV can mask the LEDs at certain times, for example when it turns around.

We can see in the figure that the trajectory reconstructed by fusion is much more complete than those estimated by a single camera. The holes in the trajectory are filled by fusion. It is also important to note that the trajectory remains smooth, i.e. there are no jumps when switching from one camera to another or when we use both. This shows that the entire system is correctly calibrated and the estimates are accurate.

Qualitative evaluation in real conditions

In the previous subsections, we quantitatively assessed the performance of our system under different conditions using a ground truth provided by a motion capture system. These experiments were carried out in our laboratory and it is questionable how robust the system is under real conditions.

The approach developed is very robust to illumination changes. Although the approach was initially developed for a very dark environment (inside the silo), it worked well in our laboratory but also in a test silo. This test silo is a basin in a shed with a transparent roof. The system is therefore subjected to large lighting variations depending on the time of day and weather conditions. As we can see in Figure 3, there are reflections of the sun directly in the basin and shadows. The system has been in operation for more than a year now without any failure. This is due to the use of high power LEDs that represent strong light sources directly visible to the camera. The filtering algorithms explained in the article also make it possible to reject potential false detections.

Conclusion

In this paper, we proposed a system and the associated algorithms to accurately localize a mobile in a large industrial environment, whether it is a mobile robot, an AGV, or in our case an ROV. This system relies on the use of a fisheye camera observing LEDs placed on the mobile. We proposed two methods for calibrating the system: the first one uses known points in the scene, the second one is an auto-calibration method. We have studied the accuracy of the system as a function of the calibration, and by varying several parameters as the number of points used, the introduction of errors in the coordinates of the points, etc. We have shown that the number of points used does not matter, but it is important that the points are known accurately. The auto-calibration method makes it possible to overcome this constraint. We have proved both quantitatively thanks to a motion capture system in the laboratory and qualitatively on the actual site that our system makes it possible to localize a mobile with a high accuracy (up to 22mm). In our experiments, we successfully localized a mobile robot in our laboratory and an ROV operating in an industrial environment.

Fig. 1

 1 Fig.1Geometric layout of the LEDs on the mobile. The luminous pattern is composed of 3 or 4 LEDs (LED 0 is not mandatory).

Fig. 4

 4 Fig.4The Unified Projection Model (UPM).

Fig. 5

 5 Fig.5Reprojection of the known points in the camera image. The red points are the selected (ground-truth) points, the green points are the reprojected points after extrinsic calibration. The error is less than one pixel per point.

Fig. 7

 7 Fig.7The validation environment in our laboratory. The two blue cameras on top (about 3.60 meters from the ground) are the fish-eye cameras used for localization.

Fig. 8

 8 Fig.8The model ROV with 4 LEDs used for validation in our laboratory. The ROV is remotely operated thanks to the Wifibot attached to its front.

Fig. 9

 9 Fig.9Estimated trajectories (in red and blue) using different numbers of points with known 3D positions for the calibration step. The ground-truth provided by the Vicon is plotted in green.

Fig. 10

 10 Fig. 10 Estimated trajectories using 4 points with noisy 3D positions for the calibration. The ground-truth is plotted in green.

Fig. 11

 11 Fig. 11 Estimated trajectories using different numbers of ROV poses for the calibration. The ground-truth is plotted in green.

FusionFig. 12

 12 Fig.12Fusing the results of multiple cameras. Left: the estimated trajectory obtained by Camera 1. Middle: the estimated trajectory obtained by Camera 2. Right: the trajectory obtained by fusing the results of both cameras. As we can see, this method fills the holes in the trajectories while keeping a smooth trajectory.

Table 1

 1 Localization errors using different numbers of points with known positions for calibration 6.2 Calibration using points with known 3D positions

	Calibration Configurations Position Error (mm) Orientation Error (degree)
	4 points	21.92 ± 11.42	0.54 ± 0.66
	10 points	23.73 ± 15.36	0.54 ± 0.68
	18 points	25.81 ± 19.06	0.54 ± 0.67

Table 2

 2 Localization errors using 4 points with noisy 3D positions for calibration.

	Calibration Configurations Position Error (mm) Orientation Error (degree)
	4 points + 50mm error	83.39 ± 64.13	0.77 ± 0.81
	4 points + 100mm error	143.66 ± 109.52	1.37 ± 1.32
	4 points + 200mm error	253.60 ± 185.04	4.06 ± 3.26

Table 3

 3 Localization errors using different numbers of ROV poses for calibration 6.4 Calibration by moving the ROV Next, we evaluate the localization performance of the calibration method by displacing the ROV. The localization errors are summarized in Table

	Calibration Configurations Position Error (mm) Orientation Error (degree)
	1 pose	306.68 ± 276.71	2.47 ± 1.64
	50 poses	58.57 ± 63.46	0.94 ± 0.99
	100 poses	67.71 ± 69.88	0.90 ± 0.98
	150 poses	49.07 ± 32.13	0.61 ± 0.76
	214 poses	31.21 ± 16.32	0.56 ± 0.75
	353 poses	28.14 ± 14.89	0.56 ± 0.75
	439 poses	27.62 ± 14.85	0.55 ± 0.74


J G Tartarin
email: author.tartarin@laas.fr

Non-Destructive Techniques For Evaluating The Reliability Of High Frequency Active Devices

SiGe and GaN technologies have achieved rapid development over the last two decades. High level of RF circuit integration on Si low cost substrates open the way for large development of SiGe HBTs, while needs for high power density make GaN HEMT a key technology for solid state power modules. As both of these technologies achieve very elevated frequencies, they become strong contenders to GaAs technologies. Then reliability studies are needed to improve the process at the lower technology readiness level scale, and to stabilize the technological process till the final qualification step. To make an efficient diagnostic on the causal origin of the physical root mechanisms involved during the application of a stress, a multi-tool approach is mandatory to secure the diagnostic. In this paper, case studies on SiGe HBT and GaN HEMT stressed devices are proposed through the cross-analysis of low frequency noise spectral densities, of electrical transient measurements, and of TCAD simulations.

Introduction

Over the last five decades, the modern era of telecommunications has developed considerably, beyond military and space applications, and with consumer applications (mobile telephony, wireless sensors networks or home automation systems). All these systems need to modulate a high frequency carrier (ranging between few MHz to more than one hundred GHz) with analog or digital data.

The spreading of the markets and related applications have been made possible due to the emergence of solid state active devices, from the development of silicon based transistors in the early 1960s to the various available solutions using GaAs, silicon(germanium) (Bi)CMOS, or also GaN technologies. First, it must be considered that technological development (both with top-down or bottom-up strategies) involves customers and foundries; this means that the final context of use of a technology can be very different according to the customer targeted markets (and subsequent performance, size, cost, environment, expected operating life, …). Therefore, strategies of development of each party must be combined as far as possible: this can become hard to meet those requirements (development of general purpose devices for foundries, requirement of specific performances and reliability tests for customers).

In section 2, development strategy of solid state devices towards higher frequency and elevated power are discussed, also considering the qualification final step. Considering the roadmap of a qualification, the first studies are focused on the development and on the improvement of demonstrator devices at low technology readiness levels (TRL 3 initial gate), then the implementation in the final product needs to pass through the critical threshold of TRL 5 main gate when reliability analysis usually starts. From a simple qualification procedure, first generic stresses are performed both on die or packaged devices. Then customers carry out specific tests according to the mission of the designed system. Anyway, considering the most stringent operating conditions, stresses involving d.c life tests, RF life tests, thermal cycling or also radiations tests are at least good indicators to evaluate the reliability of a technology.

In section 3, the outline of a procedure for the characterization and for data analysis is proposed: usually, specific characterization procedures are established, closely bind to the device under test specificities. Various experimental workbenches are available to achieve studies on virgin and stressed devices, and to improve a technology process to its qualification level. From the different sets of experimental tools, it can be distinguished destructive and non-destructive (non-invasive) techniques. If the first ones are able to reveal the structural defects (mechanical or chemical, in the surface or in the volume of the stacked layers), it is not possible to definitely conclude about their impact on the active device performance. On the other hand, non-invasive techniques are good markers of the effects of defects on electrical performances, but it is more speculative to determine their precise locations and physical nature.

Section 4 focuses on two case studies concerning HBT devices (SiGe) and the disruptive GaN HEMT technology. Cross-analyses between transient and harmonic electrical characterizations, as well as electrical low frequency noise (LFN) measurements are summarized. For each case study, TCAD models are used to consolidate the experimental results.

2) High Frequency applications: consequences of downscaling solid state devices

The mainstream of the developments of solid state technologies target the increase of the frequency carrier and the improvement of the maximum achievable power, according to radio-link budgets studies. If both sides of the emitter and receiver are of prime importance to make the telecommunication system efficient, the most critical part concerning the reliability of a system is the emitter module; elevated d.c and RF powers are achieved, together with thermal management at the device level. Therefore, the reliability studies bringing up to the qualification of a process are mainly oriented towards elevated current or voltage (direct or reverse) and elevated temperatures (with or without d.c or RF electrical signal Then, the main timeline that make up the development of a technological process is presented.

The demonstrators which are developed from the lower Technology Readiness Levels scale up to the validation in relevant environment at TRL 5. Then a new phase starts at product demonstration level TRL 6; it is quite unreasonable to consider an innovation process as linear, and iterations are needed to master a process both on the performances and on the reliability figures of merit. In keeping with the high speed device development objective, aggressive geometrical downsizing rules (vertical and lateral scaling of devices with optimized lengths, widths, heights) and new materials (lattice and thermal mismatches, spontaneous polarization or piezoelectricity) bring new paradigms the manufacturers have to manage with. It becomes particularly difficult when solutions to achieve good performances metrics are balanced with targeted reliability figures of merit (at least maximum junction temperature Tj and elevated mean time to failure MTTF). This cycle involves different partners from the lowest TRL scale (where academic institutes are largely associated during fundamental and technological levels studies) to the highest TRL scale concerning the industry management of the final portfolio (TRL 8 for technology qualification or initial market, and TRL 9 for final successful missions and market expansion). Fig. 1 gives a brief overview of frequency and power performances, with time to failure metrics for SiGe Heterojunction Bipolar Transistors (HBT) and GaN High Electron Mobility Transistors (HEMT); some of the best academic and commercial figures of merits are reported (transition or maximum oscillation frequency, power density, mean time to failure). The step between academic studies and end-of-line manufacturing is relative to the field of studies concerning solutions to heal and stabilize a process, and then to defining the security operating area (SOA) for customers.

Reliability studies and associated technological development concerns TRL scale mainly focused between TRL 5 to TRL 7. TRL scale is used for decision making in terms of technological investments; it is strongly correlated to the expected markets. The robustness of the development model depends on the market volume (placement versus other competitors on performance / price of a technological process) and the end-users variety. This last consideration makes it hard to define the decision maker about the generic stress conditions to

HBT SiGe

230/300 [START_REF] Nicolson | A low-voltage SiGe BiCMOS 77-GHz automotive radar chipset[END_REF] 505/720 [START_REF] Heinemann | SiGe HBT with ft/fmax of 505/720 GHz[END_REF] 6,5 mW/µm² [START_REF] Nicolson | A low-voltage SiGe BiCMOS 77-GHz automotive radar chipset[END_REF] @77 GHz 10 6 hours [START_REF] Kim | Considerations for forward active mode reliability in an advanced heterojunction bipolar transistor[END_REF] (Tj=125°C or Jc=10mA/µm²))

GaN HEMT 320/580 [START_REF] Micovic | Ka-band LNA MMIC's Realized in Fmax > 580 GHz GaN HEMT Technology[END_REF] (Lg=20 nm)

22 W/mm [START_REF] Ohki | An over 20-W/mm S-Band InAlGaN/GaN HEMT with SiC/diamond-bonded heat spreader[END_REF] (2-4 GHz) GaN HEMT 115/155 [START_REF] Moron | 12W, 30% PAE, 40 GHz power amplifier MMIC using a commercially available GaN/Si process[END_REF] (Lg=100 nm) 115/155 [START_REF] Floriot | GH25-10 : New qualified power GaN HEMT process from Technology to Product overview[END_REF] (Lg=250 nm)

>3.3 W/mm [START_REF] Moron | 12W, 30% PAE, 40 GHz power amplifier MMIC using a commercially available GaN/Si process[END_REF] @ 30 GHz >4,5 W/mm [START_REF] Floriot | GH25-10 : New qualified power GaN HEMT process from Technology to Product overview[END_REF] @ 20 GHz <10 6 hours [START_REF] Moron | 12W, 30% PAE, 40 GHz power amplifier MMIC using a commercially available GaN/Si process[END_REF] (Tj=200°C) 10 6 hours [START_REF] Floriot | GH25-10 : New qualified power GaN HEMT process from Technology to Product overview[END_REF] (Tj=200°C) 10 7 hours [START_REF] Ums | [END_REF] (Tj=225°C) develop for a given technology: automotive and ground applications criteria are less stringent than those of space and military applications, still considering the expected performances together with the environment of use. Anyway, a qualification strategy requires close cooperation between manufacturers and end-users customers.

As previously mentioned, the most critical aspects of reliability defining the SOA for active solid state devices are more or less closely related to current-voltage and thermal operating conditions for high power segment (transmitter). However, due to its small dimension, the active element of the receiver must sustain critical doses of radiation (space applications) or elevated RF jamming signals (radar). Technologies must be hardened in a representative environment (intrinsic, extrinsic related failures). The two main active devices used for high frequency applications are Heterojunction Bipolar Transistors (HBT or double-HBT) and Field Effect Transistors (FET, with various MESFET, HEMT, pHEMT or mHEMT declinations). This paper reports on intrinsic failures case studies (i.e. wear out of the device due to structural defects activated during stress); process variability or other low mastering technological events are not under the scope of this review. Two sets of reliability stresses can be distinguished for a qualification process:

-the first package of stresses concerns representative usual conditions under which the devices are supposed to be operated. Direct current life tests (HTRB, HTOL, IDQ-IGQ) can be operated at different constant biasing conditions, or by step variation (current or voltage). RF life tests (CW or pulsed with single or multi-tone signals) are used to assess the aptitude of the devices to sustain critical RF levels for high performance or rugged applications.

-the second group concerns specific extreme environment: cryogenic temperatures are used for very stable oscillators in space applications, whereas very high temperature can be largely out of the commercial applications standards (with usual stresses at 125°C instead of 85°C respectively for military and commercial specifications). Also radiation in space is much more critical than at ground level, and must be carefully considered. Even if those space and military applications represent socalled "niche" applications (i.e. usually concerning small market volumes), their impact is strongly strategic; it concerns defense, weather forecast (and resultant economy), and other commercial satellite services (open, commercial and public regulated services, safety of life, search and rescue). T-CAD models (t 0 -t final ) Electrical models (Linear, non-linear, noise)

3) Collective use of experimental setups and models for the assessment of active devices reliability

Experimental analysis

Model assisted feedback Fig. 2. Cross interpretation of non-destructive and destructive studies with electrical models for secured stress analysis. Characterizations are performed on virgin (or initial state) and stressed (or final state) devices, with possible intermediate measurements.

Results

Gathering experimental and theoretical results:

Feedback for technological process optimization and qualification

Electrical-dynamic models using T-CAD simulation facilities: -Coarse electrical effects and fine electrical manifestations (LFN, transient) : calibration of the model at t0 (virgin) -generalized model with sizing/biasing for technology feedback -Adjustment of the model after stress(es): correlation of the model evolution with characterization results

Step 6

T-CAD models Deep Level Transient Spectroscopy (I or C DLTS): only on I OUT (detection limitations for low I IN )

Evidence of e -or h + centers (E a + cross section) -correlation with LFN spectral analysis ?

Yes : active zone and/or contact resistance degrade Yes : traps in the volume or interface layers.

Comparison between Virgin/Stressed devices + analysis vs electrical or technological param.

No : no structural defect evolution in the No : no trap in the spectral window homog. materials Low Frequency Noise measurements : LFN noise current spectra of I IN and I OUT (*) DC biasing @ constant I OUT (variable V OUT ), @ constant V OUT (variable I OUT ) (cf. step 2 and step 3 for refining the DC conditions)

Evolution of 1/f ? Presence of GR centers ?

Electrical / thermal activation of lag effects (bulk effect/interface effects ?)

Pulsed measurements for discriminating input and output lag effects Methodology for device measurement @ different stress time -virgin device and stressed devices or -before, (during,) and after the application of the stress

Step 1

N Y

Step 2

Step 3

Step 4

Check DC plots (step1)

(*) Thermal characterization if needed
Step 5

Check DC plots (step1) Fig. 3. Methodology to reveal defects in a device, based on transient, harmonic electrical and noise measurements. TCAD models strengthens the experimental analysis outputs. Comparison between initial/final measurements allows to locate and to determine the nature of the defects. and their activation mode in the transistors is slightly simplified in Fig. 3 to match both for FET and HBT devices. In this presentation, the case studies are focused on the analysis of non destructive techniques and T-CAD models; non-destructive experimental setups allow the user to keep the devices (virgin or stressed) operational, and to apply a large variety of experimental procedures. Two distinct experimental setups can be differentiated:

-electrical measurements provide the user with functional data, and resulting parameters extracted from models. In addition to d.c characterizations, transient and frequency domain measurements are achieved (small or large signal, CW or pulsed, with or without temperature as a varying parameter).

-Low Frequency Noise spectra and Random Telegraph Noise reveals microscopic defects in the zones where the various currents flow (output current from the controlled source, input or leakage currents). Step 4 in Fig. 3 represent a crucial part of the survey, as the involved mechanisms that are detected in the noise spectra at the input and output of the device are very sensitive to the technological process. The noise signature is a fine marker of latent defects, some of them remain stable with time and stress, some are activated by the stress conditions.

The workflow and the subsequent interpretations can reveal the location, the nature and the operating mode of these defects by tuning adequate d.c biasing conditions, by varying the geometry or by performing measurements at different temperatures. The next section proposes some chosen case studies on SiGe HBT and on GaN HEMT devices, for which electrical and LFN measurements have been used in a complementary way to reveal the defects, and to propose potential technological solutions towards a stabilized technological process.

For very complex case studies, additional experimental setups and electrical models are developed to reduce (or remove) any speculative assumption. In Fig. 3, TCAD simulations are also used at the final step of the procedure. Model is calibrated on virgin device measurements, and is modified to account for the experienced degradation. In the next section, TCAD simulations are used as a complementary tool for providing clues on how devices can be degraded under stress application.

4) Reliability case studies on SiGe HBT and GaN HEMT

With the still evolving landscape of technologies, the emergence of SiGe BiCMOS integrated solutions in the early 1990s has shaped new architectures and new possibilities thanks to a high integration level and to low cost. However, GaAs HBTs and FETs were still positioned as frontline technologies for high performance and high frequency fields of applications. In parallel with the development of SiGe bipolar transistors, a second step was achieved in the 2000s with the advent of wide bandgap GaN technologies: thus critical areas concerning high frequency and high power applications (as well as robust low noise modules) have led to a redistribution of systems architectures. GaAs are still present in the roadmaps and portfolios of electronic devices manufacturers. But, as the benefits of GaAs technologies are more and more reduced to niche applications, most manufacturers change their strategy by centering their activities on SiGe or GaN solutions. The first case study gathers some reliability analysis on precursor devices of SiGe HBT (BiCMOS technologies), while the second one concerns some achievements on reliability of GaN HEMT devices.

Electrical and LFN measurements and models for reliability case studies on SiGe HBTs

Low frequency noise and random telegraph signal are sensitive metrics for determining defects in the active devices [START_REF] Labat | Link between low frequency noise and reliability of compound semiconductor HEMT and HBTs[END_REF][START_REF] Vandamme | Review of lowfrequency noise in bipolar tranzsistors over the last decade[END_REF][START_REF] Mukherjee | Random telegraph noise in SiGe HBTs: Reliability analysis close to SOA limits[END_REF], as it can reveal defects prior to any stress application, while most electrical characterization setups are not able to track those defects (except for DLTS measurement techniques). For long-term operation, the main reliability issue in SiGe HBTs is the cumulative degradation of the base current that occurs under combined high emitter current and high collector voltage stress. An illustration of the different zones of SOA, soft degradation and catastrophic failure, is given with the output characteristics of a device in Fig. 4. The stress plan for long term reliability studies can be developed from these d.c biasing and RF plots (dynamic loadline), whereas extreme or harsh conditions are defined considering very low or very high temperatures (or wide temperature range), or also radiation-rich environment (or their combination). The reliability of SiGe HBT devices is a complex issue as various physical limiting mechanisms are involved. These effects such as hotcarrier (HC), self-heating, electromigration, impact ionization or radiated ionization interact in a complex way.

In the present study, hot carrier (HC) and ionizing radiations (IR) effects on the RF Sparameters are illustrated. HC stress was developed by the application of reverse base-emitter voltage stresses [START_REF] Borgarino | Hot carrier effects in Si-SiGe HBTs[END_REF]; these conditions are relevant for reliability concerns as the base of the HBT devices is highly doped. The d.c effects of the HC stress have been investigated first by the observation of the forward Gummel plots (Fig. 5). The base current and the ideality factor both increase with the stress, whereas the collector current remains unchanged, as a typical signature of HC stress [START_REF] Vendrame | Degradation mechanisms in polysilicon emitter bipolar junction transistor for digital applications[END_REF]. LFN spectra on SIB (relative to base current fluctuations) and SIC (relative to collector current fluctuation) have been measured and modelled, both exhibiting an increase of the 1/f component (two orders of magnitude for SIB), as depicted in Fig. 6. The noise floor also features a slight decrease after the stress, that is not discussed in this paper. By a dedicated model of the different noise sources located in the base and in the emitter (bulk and space charge region), it has been evidenced how these LFN contributors evolve with stress. Defects associated to surface recombination (at the emitter perimeter) and trap-assisted tunneling recombination mechanism have been created by the hot carriers during the stress. All these conclusions are supported by electrical d.c and LFN measurements, as well as electrical and noise modellings. More detailed analysis can be found in [START_REF] Borgarino | On the effects of hot carriers on the RF characteristics of Si/SiGe heterojunction bipolar transistors[END_REF], also with numerical simulations. TCAD model has been developed and instructed to simulate the d.c and dynamic RF signatures before and after the application of the stress. Fig. 7 and Fig. 8 show the measured and simulated change in S21 gain parameter (from [S] matrix) at constant base current. Impact of stress on the transition frequency fT is given in the inset of Fig. 7. Slight differences between measurement and simulation are due to the access pads contribution on the module (losses) and phase (electrical lengths of lines), not removed from the measured S-parameters. The simulation well reproduce the experimentally observed behavior. It must be mentioned that at constant collector current, no change is noticed on S21 (simulation and measurement). Proton irradiation have been experienced on these devices, with no change in the S-parameters at constant IC. Same conclusions have been obtained on different new generations SiGe technologies, and improvements on the base process (stabilization of the doping species) and on the passivation technique have led to mature technologies, thanks to the fine determination of events occurring within the volumes or interfaces of the active zones. Similar studies have been developed on recent technologies, featuring current gain enhancement in SiGe HBTs under inverse and forward operation modes, and have been analyzed with TCAD models [START_REF] Raghunathan | Hot-carrier-damage-induced current gain enhancement (CGE) effects in SiGe HBTs[END_REF].

Other studies on SiGe and III-V HBT technologies have taken full advantage of electrical and LFN measurements, with dedicated models. Hydrogen related defects in GaInP/GaAs HBT have been evidenced, with complex kinetic affecting surface recombination and carbon passivation, which translates into enhancement and decrease periods on gain current β as reported in [START_REF] Tartarin | Hydrogen induced degradation in GaInP/GaAs HBTs revealed by low Frequency noise measurements[END_REF]. Advanced TCAD models have been developed to optimize the performance and the reliability for very high frequency application purpose [START_REF] Shröter | SiGe HBT technology: future trends and T-CAD-based roadmaps[END_REF][START_REF] Chevalier | Si/SiGe:C and InP/GaAsSb heterojunction bipolar transistors for THz applications[END_REF]. The different trade-offs induced by the device size reduction in Si/SiGe and InP/GaAsSb HBTs (extrinsic to intrinsic base link due to lateral and vertical scaling), with isolation strategies (trench isolation) and the related increase of the thermal resistance RTH are considered together. Reliability predictive models are used for the simulation of HC damage in SiGe HBT [START_REF] Moen | Predictive physics-based TCAD modeling of the mixed mode degradation mechanism in SiGeHBTs[END_REF].

Electrical and LFN measurements and models for reliability case studies on GaN HEMTs.

Wide bandgap semiconductors feature superior material properties enabling power device operation at higher temperatures, voltages, and switching speeds than current Si or GaAs technologies. While GaAs devices total achievable power had plateaued in the 2010 th after more than three decades of developments, GaN HEMT overcomes low breakdown voltages and low thermal conductivity, and makes much higher operating voltages as well as higher current densities. These combine to enable six to ten times higher output power density with GaN than with GaAs, even at high frequencies thanks to the high mobility of electrons in the GaN quantum well. Thermal effects are one of the major issue for HF devices, as the output power is proportional to the barrier thickness, which is inversely proportional to the device cut off frequency.

However, SOA limits are not at their higher potential due to still limiting reliability issues to solve: in these complex devices, it is necessary to discriminate between electrical, thermal and mechanical effects that can combine in direct and inverse relationship. A method based on variable stimuli such as temperature, time or frequency makes it possible to distinguish between these events. Fig. 3 lays the outline for the measurement and modeling techniques, mixed with destructive techniques as proposed in Fig. 2, to define a coarse and fine analysis of defects related to the application of a stress. Among drastic intrinsic limitations for RF applications, lag effects on drain current can raise from various origins below the gated zone and in the ungated drain-gate or gate-source regions, with static or dynamic behaviors that can interfere with the carrier frequency according to the load line area. Even for power devices field effect transistors, LFN measurements are of great interest to scan the spectral signatures of the defects. By defining constant temperature, d.c quiescent drain (leakage) current, or gate-drain (gate-source) voltage measurement strategies, it is possible to perform cross-experiments according to various stimuli. Frequency (LFN spectra or HF frequency dispersion [START_REF] Nsele | BroadBand Frequency Dispersion Modeling of the Output Conductance and Transconductance in AlInN/GaN HEMTs[END_REF]), time (transient d.c or pulsed lag measurements [START_REF] Tartarin | I-DLTS, electrical lag and low frequency noise measurements of trapping effects in AlGaN/GaN HEMT for reliability studies[END_REF][START_REF] Tartarin | Diagnostic Tools For Accurate Reliability Investigations of GaN Devices[END_REF]) or temperature (DLTS) measurements have been carried out on different GaN technological processes during development steps (low TRL scale 3 to 5) or qualification steps (higher TRL scale 5 to 8). Two main mechanisms can be sorted out. The first one concerns defects that affect the carriers flowing in the 2DEG channel, and not related to the gate command. The second one is attributed to the gated zone, and represent the most challenging electrical issue for the control of the electron density in the 2DEG. A listing of common phenomena is identified in [START_REF] Meneghesso | Reliability issues of GaN High Electron Mobility Transistors[END_REF]. It has been demonstrated that, even if (de)trapping centers are identified in the channel of the device, their evolution after the application of a stress is not significant (drain current spectral density SID). However, if these centers are located in the gate-source or gate-drain interfaces, or in the GaN bulk layer, they can strongly impact the stability of the RF signal (lag and other transient effects).

Fig. 9 represents the gate current LF spectra SIG on virgin devices and HTRB stressed devices. From Fig. 9 a. SIG features a strong degradation after ageing (increase by more than 3 decades). This evolution is clearly correlated to an increase of the gate leakage current during stress as reported in Fig. 9 b. The induced degradation during HTRB stress impacts the 1/f noise source (electron carrier fluctuation in the leakage conduction path), and thus the leakage current IG can be used as a dependable marker of reliability for device selection. A unified model has been developed to account for the gate Schottky barrier height inhomogeneity [START_REF] Lazar | New Approach for an Accurate Schottky Barrier Height's Extraction by I-V-T Measurements[END_REF], for the thermal evolution of the IG-VGS plots and for the trap effects on the transient characteristics IGS(time) and IDS(time) [START_REF] Lazăr | Correlation between transient evolutions of the gate and drain currents in AlGaN/GaN technologies[END_REF]. In addition to these electrical and noise characterizations, TCAD models have been developed to successfully account for the physical origin of the evolutions on drain current and gate leakage current [START_REF] Saugnon | Fully automated RF-thermal stress workbench with S-parameters tracking for GaN reliability analysis[END_REF][START_REF] Tartarin | Gate defects analysis in AlGaN/GaN devices by mean of accurate extraction of the Schottky Barrier Height, electrical modelling, T-CAD simulations and TEM imaging[END_REF]. A synoptic of the tree structure used to discriminate between the different physical parameters that are tuned to account for the degradation process is given in Fig. 10. This structure helps to find the set of parameters evolving with stress. In addition to usual IDQ or HTRB stresses, RF stresses are also of great interest to track a) SIG LFN spectra vs frequency b) SIG (@1kHz) versus IG for two HTRB stresses Similarly, the evolutions on access resistances as proposed in [START_REF] Joh | RF power degradation of GaN High Electron Mobility Transistors[END_REF] could be interpreted as the consequence of variations in the intrinsic local biasing of the device (charges accumulation) [START_REF] Tartarin | Gate defects analysis in AlGaN/GaN devices by mean of accurate extraction of the Schottky Barrier Height, electrical modelling, T-CAD simulations and TEM imaging[END_REF],

and not as a structural change in the resistive layers constitution.

5) Conclusions and perspectives

A focus on two key technologies used for high frequency applications (SiGe HBT and GaN HEMT) is achieved, demonstrating the need to develop a multi-tool strategy for fine and secured reliability analysis; it makes use of cross-experiments in spectral and transient domains, with electrical and noise models. With the SiGe HBT evolution, BiCMOS technologies and mixed-signal applications have strongly proliferated for multiple highperformance and high frequency applications. As a consequence, reliability issues has come to the forefront. The safe-operating area of SiGe HBTs is strongly defined by Hot-Carrier mechanisms, that have been presented through a collective use of experimental and simulation tools. Althrough these effects have been widely studied over the past two decades, a T-CAD based predictive approach is required to estimate the impact of HC when scaling device performance and technology. The roadmap to high frequency and high power applications takes advantage of GaN HEMT technologies, with challenges to overcome concerning self-heating dissipation by various strategies. Besides the elevated junction temperature that affects the drift velocity and mobility [START_REF] Nuttinck | Thermal analysis of AlGaN-GaN power HFETs[END_REF], and that accelerate aging of the device [START_REF] Kuball | Thermal management and device failure assessment of high-power AlGaN/GaN HFETs[END_REF], the presence of charges still need to be mastered to avoid any current collapse or even virtual gate effects, or any other mechanism inducing high frequency dispersion effects. LFN and electrical measurements are used with TCAD simulations to determine with a high degree of confidence the origin and the kinetic of the defects responsible for the degradation of the device.

Fig. 1 .

 1 Fig. 1. Overview of frequency, power and mean time to failure (MTTF) from academic studies and manufacturer production line, for SiGe Heterojuntion Bipolar Transistors and GaN High Electron Mobility Transistors (references are listed in the final relevant section).

Fig. 2

 2 Fig. 2 exposes the cooperation between experimental analysis modelling techniques developed to diagnose on failure root mechanisms and to optimize the technological process. The proposed methodology used to identify the defects

  Fisrt order thermal effects (constant P DC vs temperature or constant temperature vs P DC ) No impact of stress on the control zone (base or gate) Base or Gate defect : more measurement are needed to determine the underlying physical mechanisms N Y Transient measurements: V IN or I IN (time), V OUT or I OUT (time), gain(time) -first diagnostic Virgin and Stressed devices -find out correlation laws « in-out » Nota : then use only Long configuration (DC-like behavior). Transistor (HBT or FET) Direct Current I OUT -V OUT and I OUT -V IN measurements Signal Analyser : Short/Medium/Long (S/M/L) Plots move with Short/Medium/Long integration time ? (*)

Fig. 4 .

 4 Fig. 4. Illustration of output characteristics of SiGe HBTdevices with thermal, hot carrier and electrical limits defining the area of safe operating mode, degradation and catastrophic failure. Example of a dynamic load-line is represented, close to the SOA limit for power applications. Current and voltage are in arbitrary unit.

Fig. 5 .Fig. 6 .

 56 Fig. 5. Measured and simulated Gummel plots before and after applying the stress. The ideality factors are all experimental.

Fig. 7 .

 7 Fig. 7. Changes in the extrinsic scattering parameter S21 measured from 40 MHz to 30 GHz at constant base current IB during the stress. The inset represents the variations of the measured extrinsic transition frequency ft during the stress, measured at constant IB and constant IC.

Fig. 8 .

 8 Fig. 8. Changes in the intrinsic transmission scattering parameter S21 simulated from 40 MHz to 40 GHz at constant IB (VCE=1 V) before and after stress.

Fig. 9 .

 9 Fig. 9. a) Gate current LFN spectra for virgin (black) and HTRB stressed devices (greys) and b) linear plot from SIG vs IG for virgin and stressed devices (two HTRB stresses @ 150°C and 25°C).

MainFig. 10 .

 10 Fig. 10. Tree structure for the analysis of TCAD model before/after stress. Application to GaN HEMT technologies. Electrical parameters (green boxes) are bind to physical parameters; fixed technological parameters (grey) or stress dependent parameters (blue & cyan) are used to fit the main d.c electrical markers (orange and green) defining the output & transfer characteristics, leakage current [31]).

Fig. 11 .

 11 Fig. 11. Schematic of RF/thermal stress test setup: ① power amplifier, ② coupler, ③ four-way power divider, ④ SPDT, ⑤ attenuators, ⑥ SP4T 50Ω RF load terminations ⑦ SP4T.

Ft/Fmax GHz (gate length for FET) Output power density (@freq) MTTF in hours @Temperature Examples from academic best results

	Nota: as the
	proposed FoM are
	rarely available
	together, this
	table gathers
	some examples

Acknowledgements

Part of the work was supported by the French Ministry of Defence (DGA), the French National Research Agency (ANR), and various European Programs. The author is thankful to the partners from SiGe-microsystems, ST-microelectronics, III-V Lab, UMS, OMMIC, to our academic partners (particularly L. Bary, M. Borgarino for SiGe studies) and to our PhD students (G. Soubercaze, G. Astre, S. Karboyan, O. Lazar, D. Saugnon).


Monografías Matemáticas

García De Galdeano

Jacky Cresson

Jordy Palafox

LIE ALGEBRAS AND GEOMETRIC COMPLEXITY OF AN ISOCHRONOUS CENTER CONDITION

Keywords: Isochronous center, planar polynomial vector fields, linearisation AMS classification: 34C20, 30G20, 37G05

Using the mould formalism introduced by Jean Ecalle, we define and study the geometric complexity of an isochronous center condition. The role played by several Lie ideals is discussed coming from the interplay between the universal mould of the correction and the different Lie algebras generated by the comoulds. This strategy enters in the general program proposed by J. Ecalle and D. Schlomiuk in [4] to study the size and splitting of some Lie ideals for the linearisability problem.

§1. Introduction

In this article, we are interested in characterizing isochronous center of polynomial real vector fields which can be written using complex coordinates as (see [START_REF] Francoise | Géométrie analytique et systèmes dynamiques[END_REF]):

X = ξ(x∂ x -y∂ y ) + P(x, y)∂ x + Q(y, x)∂ y , (1)
where ξ ∈ C such that ξ 2 = -1, ∂ x = ∂ ∂x , ∂ y = ∂ ∂y , P(x, y) is a polynomial of a maximal degree d given by P(x, y) = 2≤i+ j≤d p i, j x i y j , with ȳ = x and the coefficients p i, j ∈ C satisfy q i, j = p j,i . Definition 1. A vector field is said to be isochronous if all its orbits in a neighbourhood of an equilibrium point are periodic with the same period.

We have a characterization of isochronicity (see [START_REF] Chavarriga | A survey of isochronous centers[END_REF]):

Theorem 1. Such a vector field X is isochronous if and only if X is linearisable by an analytic change of coordinates.

In this article, we use two classes of objects to study linearisation: prenormal forms and the correction introduced by J. Ecalle and B. Vallet in [START_REF] Ecalle | Correction and linearization of resonant vector fields and diffeomorphisms[END_REF]. The main characteristic of these objects is that they are computable algorithmically and thus provide explicit criterion of linearisability and moreover that they possess a rich algebraic structure coming from the use of the mould formalism. The linearisability is equivalent to the fact that any prenormal form of a given vector field is reduced to the linear part or its correction is trivial (see [START_REF] Ecalle | Correction and linearization of resonant vector fields and diffeomorphisms[END_REF]).

More precisely, let us consider a vector field in a prepared form:

X = X lin + n∈A(X) B n , (2)
where X lin = ξ(x∂ x -y∂ y ), x = ȳ, the B n are homogeneous differential operators of degree n = (n 1 , n 2 ) ∈ Z 2 , with n 1 ≥ 0, n 2 ≥ 0 and all most one of the n i is equal to

-1 i.e B n (x m ) = c n,m x n+m , c n,m ∈ C, m ∈ N 2 , x = (x 1 , x 2 ) and x m = x m 1 1 x m 2 2 .
The set of degree associated to X is denoted by A(X).

Considering a polynomial vector field of the form (1), the previous decomposition leads to three types of differential operators:

B (i-1,k-i) = x i-1 y k-i (p i,k-i x∂ x + pk-i+1,i-1 y∂ y ), B (-1,k) = p 0,k y k ∂ x , B (k,-1) = p0,k x k ∂ y , (3)
where 2 ≤ k ≤ d and 1 ≤ i ≤ k.

A prenormal form associated to X is an expression of the form:

X = X lin + n∈A * (X) M n B n ,
where n is a word obtained by concatenation of letters in A(X), A * (X) is the set of words constructed on A(X), and for For n ∈ A * (X), the quantity ω(n) is called the weight of the word n. A word n ∈ A * (X) is said to be resonant if ω(n) = 0.

n = n 1 • n 2 • ... • n r ∈ A * (X), B n = B n 1 • B n 2 • • • • • B n r
A fundamental role is played by resonant words as they are in one to one correspondence with monomials in the vector field X which forming an obstruction to linearisability. Precisely, we have: Proposition 2. The mould M • of a prenormal form satisfies M n = 0 for non-resonant word n.

In the same way, the correction of X is given by Carr(X) = • Carr • B • where the mould Carr • is defined algorithmically (see [START_REF] Ecalle | Correction and linearization of resonant vector fields and diffeomorphisms[END_REF]) and satisfies Carr ∅ = 0 and Carr n = 0 for all non-resonant word n ∈ A * (X).

Using these two objects, the linearisation is equivalent to the fact that:

• M • B • = 0, (4)
where the mould M • is either the mould Carr • of the correction or a mould Pran • of a given prenormal form. Isochronicity is then equivalent to find the conditions on the coefficients of a polynomial P such that such a formal series vanishes. The interest of the mould formalism with respect to others approaches is that it separates the part depending on the coefficients of P (the comould contribution) from the universal one depending only on the alphabet A(X) generated by the vector field (the mould part). Moreover, as already pointed out by J. Ecalle and D. Schlomiuk (see [START_REF] Ecalle | The nilpotent and distinguished form of resonant vector fields or diffeomorphisms[END_REF], §.10,p.1474), the main difficulty in order to solve the isochronicity problem is that the "ideal I generated by (finitely many) Taylor coefficients of X is unwieldy, unstructured and lacking in truly canonical bases." They propose to replace the commutative ideal I by some Lie ideals which arise naturally in the mould formalism approach (See [START_REF] Ecalle | The nilpotent and distinguished form of resonant vector fields or diffeomorphisms[END_REF], §.10,p.1475-1475 for more details).

The aim of this paper is to interpret some classical isochronous center conditions in term of Lie ideals following the general program proposed by J. Ecalle and D. Schlomiuk. We are leaded to define a notion of geometric complexity for an isochronous center condition.

The plan of the paper is as follows: by studying first the isochronous center conditions for quadratic polynomials, we observe that the uniform and holomorph isochronous centers are related to the nilpotence and triviality of the Lie algebra generated by the family of differential operators B n , n ∈ A(X) and the resonant set respectively. These conditions do not depend on the value of the underlying mould and in some way are the simplest one. We then generalize these results for polynomial vector fields with an arbitrary large degree improving on an unpublished paper of B. Schuman [START_REF] Schuman | Correction et linéarisation des champs de vecteurs polynomiaux[END_REF]. We then discuss the notion of geometric complexity for an isochronous center condition and state a conjecture. §2. Quadratic isochronous center

We consider a real quadratic planar vector field X in its complex representation on C 2 defined by: ẋ = ξx + p 2,0 x 2 + p 1,1 xy + p 0,2 y 2 , ẏ = -ξy + p0,2

x 2 + p1,1 xy + p2,0 y 2 , (5)
where y = x, p i, j ∈ C, i + j = 2, i = 0, 1, 2.

The following results can be founded in [START_REF] Francoise | Les conditions du centre pour un champ de vecteurs quadratiques[END_REF][START_REF] Zoladek | Quadratic systems with center and their perturbations[END_REF]:

Theorem 3. A real quadratic planar vector field X is an isochronous center if and only if at least one of the following conditions is satisfied:

i) p 1,1 = p 0,2 = 0 (Holomorphic isochronous center), ii) p 2,0 = p1,1 , p 0,2 = 0 (Uniform isochronous center), iii) p 2,0 = 5 2 p1,1 , |p 1,1 | 2 = 4 9 |p 0,2 | 2 , iv) p 2,0 = 7 6 p1,1 , |p 1,1 | 2 = 4|p 0,2 | 2 .
Our approach suggests to look for these conditions by studying first their consequences on the Lie algebra generated by the B n , n ∈ A(X). Lemma 4. A quadratic vector field satisfying p 2,0 = p1,1 , p 0,2 = 0, (6)

or p 1,1 = p 0,2 = 0, (7)
is such that the Lie algebra b generated by B(X) = {B n } n∈A(X) is nilpotent of order 1, i.e. for all n, m ∈ A(X) we have

[B n , B m ] = 0. ( 8
)
Using this Lemma and the structure of the series, we recover easily that these conditions correspond to an isochronous center. Indeed, we have the following general observation: Lemma 5. If the Lie algebra b is nilpotent of order 1 then any prenormal form associated to a mould Pran • is reduced to Pran

(X) = X lin + n∈A(X), ω(n)=0 Pran n B n , (9)
and the correction is given by

Carr(X) = n∈A(X), ω(n)=0
Carr n B n .

Proof. By the classical projection Theorem (see [START_REF] Serre | Lie algebras and Lie groups[END_REF] Theorem 8.1 p.28), we have

• M • B • = r≥1 1 r n∈A * (X) M n [B n ], (11)
where

[B n ] = B n for n ∈ A(X) and [B n ] = [B n r , [B n r-1 , . . . , [B n 2 , B n 1 ] . . . ]] for n = n 1 . . . n r ∈ A * (X), r ≥ 2.
Indeed the constant term of the series must be zero as this formal series corresponds to a vector field. As a consequence, under these conditions we obtain

• M • B • = r≥1 1 r n∈A(X) M n B n , (12)
The mould M • is equal to zero on non-resonant words so that the right side is reduce to a sum over resonant letters.

As already noted, the value and nature of the moulds are not important. Only the Lie algebraic structure of b is taken into account.

The proof that conditions i) and ii) lead to isochronous center is then easily deduces. Indeed, quadratic vector fields do not produce resonant letters. As a consequence, Lemma 5 implies that Pran(X) = X lin and Carr(X) = 0.

The proof of Lemma 4 goes as follows: First, as p 0,2 = 0 in the two cases, we have

B (2,-1) = B (-1,2) = 0. A simple computation gives [B (1,0) , B (0,1) ] = x 2 y p 1,1 p1,1 -p 2,0 ∂ x + p1,1 -p 1,1 + p2,0 ∂ y . ( 13
)
As a consequence, a condition for b to be nilpotent of order 1 is

p 1,1 p1,1 -p 2,0 = 0. ( 14
)
Then if p 1,1 = 0 or p1,1 = p 2,0 the Lie algebra generated by the vector fields in B(X) is nilpotent of order one.

We generalize these conditions for a homogeneous polynomial perturbation of degree d ≥ 2 in the following. §3. Lie algebras generated by polynomial vector fields

The previous Section indicates that the Lie algebra b generated by the set B(X) of comoulds associated to X plays a central role in the understanding of some center conditions. In this Section, we derive some useful results which will be used in our study of uniform and holomorph center conditions.

Descending central series and Nilpotent Lie algebra

In this section, we give some reminders about Lie algebras. We refer to [START_REF] Reutenauer | Free Lie algebras[END_REF] and [START_REF] Serre | Lie algebras and Lie groups[END_REF] for more details. Definition 3. Let g be a Lie algebra, we define its descending central series by:

C 1 (g) = g, C 2 (g) = [g, g], C n+1 (g) = [g, C n (g)] for i ≥ 2.
Due to Lemma 5, we are interested in Nilpotent Lie algebra which are defined as follows: Definition 4. A Lie algebra g is nilpotent if there exists an integer i such that C i (g).

Preliminaries

We begin with some computations on Lie brackets: Lemma 6. Let X be a polynomial vector field of degree d ≥ 2 and B(X) its associated set of homogeneous differential operators. We have for n = 2, . . . , d and i = 1, . . . , n:

[B (n,-1) , B (-1,n) ] = n|p 0,n | 2 (xy) n-1 x∂ x -y∂ y , [B (i-1,n-i) , B (n-i,i-1) ] = (n -i)p i,n-i U i -(i -1)p n-i+1,i-1 Ūi (xy) n-1 ∂ x - (n -i) pi,n-i U i + (i -1) pn-i+1,i-1 U i (xy) n-1 ∂ y , (15)
where U i = p n-i+1,i-1pi,n-i .

The previous Lemma gives a special role to the quantities U i , i = 1, . . . , n and also to the coefficients p i,n-i . Using this remark, we are able to derive a class of explicit conditions for which the associated Lie algebra b is nilpotent of order 1.

Nilpotent Lie algebras -Uniform conditions

A trivial condition in order to ensure that the Lie Bracket [B (i-1,n-i) , B (n-i,i-1) ] is zero is to set U i = 0 for i = 1, . . . , n. This condition coupled with p 0,n = 0 annihilating the bracket [B (n,-1) , B (-1,n) ] is in fact more powerful. Indeed, we have: Lemma 7. Let X be a polynomial vector field of the form (1) where P is homogeneous of degree d. If the coefficients of P satisfy p d-i+1,i-1pi,d-i = 0, i = 1, . . . , d, and p 0,d = 0, (

then the Lie algebra b is nilpotent of order 1.

Proof. The condition p 0,d = 0 implies not only that the bracket [B (n,-1) , B (-1,n) ] is zero, but also that B (n,-1) = B (-1,n) = 0. Then, the Lie algebra b is only generated by the operators B (i-1,d-i) for = 1, . . . , d. Again the condition p d-i+1,i-1pi,d-i = 0, i = 1, . . . , d implies that all the operators B (i-1,d-i) are of the form

B (i-1,d-i) = p i,d-i x i-1 y d-i X E
, where X E is the classical Euler vector field defined by X E = x∂ x + y∂ y . The Euler vector field acts trivially on each monomials. Precisely, we have X E (x i-1 y d-i ) = (d -1) x i-1 y d-i . The Lie bracket of B (i-1,d-i) with B ( j-1,d-j) is then easily computed. We obtain:

[B (i-1,d-i) , B ( j-1,d-j) ] = p i,d-i p j,d-j x i-1 y d-i X E [x j y d-j ]∂ x + X E [x j-1 y d-j+1 ]∂ y -x j-1 y d-j X E [x i y d-i ]∂ x + +X E [x i-1 y d-i+1 ]∂ y , = 0, ( 17
)
which implies that C 2 (b) = 0. This concludes the proof.

Resonant subset -Holomorphic conditions

We have pointed out the special role of the coefficients p i,n-i in the computations. We introduce the resonant subset b res of b generated by all the Lie bracket of [B n ] such that ω(n) = 0, n ∈ A(X) * . As in the previous Section, we prove the following result: Lemma 8. Let X be a polynomial vector fields of the form 1. We assume that the coefficients of P satisfies p i,n-i = 0, n = 2, . . . , d, i = 0, . . . , d -1.

(

) 18
then the subset b res of the Lie algebra b is trivial, i.e. b res = {0}.

Proof. As usual, the previous condition implies that p 0,n = 0 for all n = 2, . . . , d so that B (n,-1) = B (-1,n) = 0. We then concentrate on the operators B (i-1,n-i) for i = 1, . . . , n and n = 2, . . . , d. The assumption implies that B (i-1,n-i) = 0 for i = 2, . . . , n -1 and A resonant word must mix letters with positive weight and negative weight. Assume that a given word begins with some (0, n 1 -1), . . . , (0, n k -1) for a given k ∈ N * . One must have at least one operator of the form B (m 1 -1,0) in order to have a resonant sequence. As [B (0,n 1 -1),...,(0,n k -1) ] depends only on y∂ y and B (m 1 -1,0) of x∂ x , the corresponding Lie bracket is zero. Then all [B n ] such that ω(n) = 0 is zero. This concludes the proof.

B (0,n-1) = y n pn,0 ∂ y , and
B (n-1,0) = x n p n,0 ∂ x . (19)
As expected by J. Ecalle and D. Schlomiuk in [START_REF] Ecalle | The nilpotent and distinguished form of resonant vector fields or diffeomorphisms[END_REF], this question is more tractable than the initial one. In particular, our two examples already give conditions for which the minimal degree is attained by an isochronous center condition. As a consequence, we can look over the set of center condition to compute in each case the couple g(C). We have Conditions ensuring the nilpotent character of b or the triviality of b res are always of degree 1 as they can be read on the Lie bracket of homogeneous vector fields B n , n ∈ A(X). Other isochronous center conditions depend on the interplay between moulds and comoulds and generate polynomial identities of degree at least 2. As a consequence, we are leaded to the following conjecture: Conjecture: The number glin(d) is equal to (d, 1) corresponding to holomorphic isochronous center.

Definition 2 .

 2 denotes the usual composition of differential operators and M n ∈ C. Following J. Ecalle[START_REF] Ecalle | Singularités non abordables par la géométrie[END_REF], we denote by M • and we call mould the map from A * (X) to C defined M • [n] = M n . The corresponding map for the differential operators B n is called a comould. Let ω : (A * (X), .) → (Z, +) be the morphism defined for any letter n ∈ A(X) by ω(n) = n, λ , with λ = (ξ, -ξ). Then we have for any word n = n 1 •...•n r that ω(n) = r j=1 ω(n j ).

for n = 2 ,

 2 . . . , d. A simple computation shows that for all n, m ∈ 2, . . . , d, we have [B (0,n-1) , B (m-1,0) ] = 0. (20) Moreover, we have [B (0,n-1) , B (0,m-1) ] = y n+m-1 c n,m ∂ y , where c n,m is a constant. As a consequence, we have [B (0,n 1 -1),...,(0,n k -1) ] = y n 1 +•••+n k -1 c n 1 ,...,n k ∂ y . The same is true for [B (m 1 -1,0),...,(m r -1,0) ] replacing y by x.

g

  Holo (d) = (d, 1), g Uni (d) = (d + 1, 1) if d is even and g Uni (d) = (d + 2, 1) if d is odd. (26)

In[START_REF] Conti | Centers of planar polynomial systems. a review[END_REF], this condition is stated for the underlying real vector field.

Jacky Cresson and Jordy Palafox Laboratoire de Mathématiques et de leurs Applications de Pau, CNRS UMR 5142 Université de Pau et des Pays de l'Adour Av. de l'Université, 64013 Pau Cedex, France jacky.cresson@univ-pau.fr and jordy.palafox@univ-pau.fr

§4. Uniform isochronism

Following R. Conti [START_REF] Conti | Centers of planar polynomial systems. a review[END_REF], an isochronous center of a vector field X is said to be uniform if all the periodic orbits have the same period. This condition is equivalent to (see [START_REF] Conti | Centers of planar polynomial systems. a review[END_REF], §.19, Definition 19.1 p.28) the following relation on the polynomial P: yP(x, y) = xP(x, y).

(UI)

This conditions induces specific relations on the coefficients of the polynomial: Lemma 9. A polynomial P of degree d satisfies the uniform isochronicity condition (UI) if and only if the coefficients p i, j , i

Using this Lemma and the structure Lemma 7 for the Lie algebra generated by the comoulds in B(X), we deduce the following Lemma: Lemma 10. Let X be a polynomial vector field of the form (1) where P is homogeneous of degree d. Assume that the coefficients of P satisfy

Moreover, if d is odd with d = 2m + 1, we suppose that p m+1,m = 0. Then, the vector field X is isochronous.

Proof. By Lemma 7, the first conditions imply that the Lie algebra b is nilpotent of order 1. We deduce from Lemma 5 that any prenormal form is reduced to

If d is even, the set of resonant letters is empty and Pran(X) = X lin . If d is odd, then we have only one resonant letter in B(X) given by n m = (m + 1, m). Then a prenormal form is given by

As p m+1,m = 0 we obtain B n m = 0 and Pran(X) = X lin . This concludes the proof. §5. Holomorphic isochronous centers

Following Conti [START_REF] Conti | Centers of planar polynomial systems. a review[END_REF], a vector field is said to satisfy the Cauchy-Riemann 1 conditions if

The Cauchy-Riemann conditions impose strong constraints on the coefficients of P: Lemma 11. A polynomial P of degree d satisfies the Cauchy-Riemann conditions (CR) if and only if p i,n-i = 0 for i = 0, . . . , n -1, n = 2, . . . , d.

We deduce:

Lemma 12. Let X be a polynomial vector fields of degree d of the form [START_REF] Chavarriga | A survey of isochronous centers[END_REF]. Assume that P satisfies the Cauchy-Riemann condition. Then the vector field X is linearisable.

Proof. Formula (12) implies that any prenormal form associated to a mould Pran • can be written as Pran

By definition all the Lie brackets in (25) belong to b res . Using Lemma 8, the Cauchy-Riemann conditions imply that b res is trivial. As a consequence, all the Lie brackets [B n ] reduce to zero for n ∈ A * (X) such that ω(n) = 0. We conclude that any prenormal form reduces to the linear one and the vector field X is formally linearisable.

The conditions on the coefficients correspond to the characterizations of holomorphic isochronous centers. §6. Linearisability and complexity Following J. Ecalle and D. Schlomiuk in [START_REF] Ecalle | The nilpotent and distinguished form of resonant vector fields or diffeomorphisms[END_REF], we introduce the following problem: Minimal complexity of the linearisability problem: Let d ∈ N * , d ≥ 2 and X a polynomial vector field of degree d given by [START_REF] Chavarriga | A survey of isochronous centers[END_REF]. We denote by lin(d) the minimal number of algebraic relations depending on the coefficients of P which induce analytic linearisability. Can we determine a bound or an explicit formula for lin(d) ?

Using our approach, we understand that this number depends on the complexity of the isochronous center conditions. In particular, each condition C i where i = 1, . . . , m of an isochronous center are determined by a finite family of polynomial C i (p) of degree c i . A natural notion of complexity is given by: Definition 5. Let S ∈ C m be an algebraic set. Consider a representation (R) of S . The complexity of the representation (R) is defined by the triplet (m, P(R), C(R)), composed by the following data:

1. The dimension of the ambient space C m , 2. The number of condition in (R), 3. The maximal degree of polynomials defining the conditions in (R).

In our case, the dimension of the ambient space is fixed by d and is given by the number m(d) of coefficients of a generic polynomial of degree d. We then introduce the notion of geometric complexity for an isochronous condition: Definition 6. An isochronous condition C is said of geometric complexity (q, m) if it admits a representation made of q polynomial identities of degree at most m.

The previous notion can be used to refine the question raised by J. Ecalle and D. Schlomiuk: Minimal geometric complexity of the linearisability problem: Let d ∈ N * , d ≥ 2 and X a polynomial vector field of degree d given by [START_REF] Chavarriga | A survey of isochronous centers[END_REF]. We denote by glin(d) the minimal geometric complexity of an isochronous center condition. Can we find a bound or a formula for glin(d) ?