Mehdi Benallegue
email: mehdi.benallegue@aist.go.jp

Abdelaziz Benallegue
email: abdelaziz.benallegue@uvsq.fr

Rafael Cisneros
email: rafael.cisneros@aist.go.jp

Yacine Chitour
email: yacine.chitour@lss.supelec.fr.

VELOCITY-AIDED IMU-BASED ATTITUDE ESTIMATION

This paper addresses the problem of estimating the attitude of a rigid body, which is subject to high accelerations and equipped with inertial measurement unit (IMU) and sensors providing the body velocity (expressed in the reference frame attached to the body). That issue can be treated differently depending on the level of confidence in the measurements of the magnetometer of the IMU, particularly with regard to the observation of the inclination component with respect to the vertical direction, rendering possible to describe the interaction with gravity. Two cases are then studied: either (i) the magnetometer is absent and only the inclination can be estimated, (ii) the magnetometer is present, giving redundancy and full attitude observability. In the latter case, the presented observer allows to tune how much the inclination estimation is influenced by the magnetometer. All state estimators are proposed with proof of almost global asymptotic stability and local exponential convergence. Finally, these estimators are compared with state-of-the-art solutions in clean and noisy simulations, allowing recommended solutions to be drawn for each case.

Introduction

The orientation, or attitude, of a mechanical system in the world is often an important part of its dynamical state and constitutes sometimes the most important variable determining the dynamics such as in the case of drones [START_REF] Gene F Franklin | Feedback control of dynamic systems[END_REF] or legged robots [START_REF] Wieber | Modeling and Control of Legged Robots[END_REF]. While it is possible for fixed-based robots to reconstruct the orientation of any link using the joint position, this is not possible for mobile robots, and specifically floating-base ones. However, a set of sensors is usually dedicated to the estimation of the orientation. These sensors, usually grouped in a set called inertial measurement units, measure the linear acceleration, including the gravitational one, the angular velocity and sometimes the magnetic field measurement, all expressed in the frame of the sensor. Thanks to this set of measurements an efficient estimation for the attitude can be built, but only when the system has negligible linear accelerations compared to gravity [START_REF] Mahony | Nonlinear complementary filters on the special orthogonal group[END_REF][START_REF] Martin | An invariant observer for earth-velocity-aided attitude heading reference systems[END_REF].

The constraint that the inertial accelerations must be negligible compared to gravity can be limiting to low dynamics motion or simply impossible to hold, especially when the system is subject to impacts such as during the case of bipedal walking. A dynamical model of the system can be used to predict the accelerations and compensate for them. This prediction can be based on the forces models, either in the case of unmanned aerial vehicles [START_REF] Martin | The true role of accelerometer feedback in quadrotor control[END_REF][START_REF] Martin | A semi-global modelbased state observer for the quadrotor using only inertial measurements[END_REF] or legged robots [START_REF] Bloesch | State Estimation for Legged Robots -Consistent Fusion of Leg Kinematics and {IMU}[END_REF][START_REF] Mifsud | Estimation of Contact Forces and Floating Base Kinematics of a Humanoid Robot Using Only Inertial Measurement Units[END_REF]. However this solution is specific to every dynamical system and requires to identify many dynamical parameters.

Another solution is to "aid" the inertial measurement unit (IMU) with independent measurements able to clear the acceleration ambiguity, such as the position in the world frame provided by GPS [START_REF] Hua | Attitude estimation for accelerated vehicles using GPS/INS measurements[END_REF][START_REF] Roberts | On the attitude estimation of accelerating rigid-bodies using GPS and IMU measurements[END_REF] or linear velocity, either expressed in the world frame [START_REF] Martin | An Invariant Observer for Earth-Velocity-Aided Attitude Heading Reference Systems[END_REF], the local frame of the sensor [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF][START_REF] Allibert | Velocity aided attitude estimation for aerial robotic vehicles using latent rotation scaling[END_REF][START_REF] Benallegue | Tilt estimator for 3D non-rigid pendulum based on a tri-axial accelerometer and gyrometer[END_REF], or a mixture of both [START_REF] Hua | Riccati nonlinear observer for velocity-aided attitude estimation of accelerated vehicles using coupled velocity measurements[END_REF], sometimes while reconstructing the velocity itself [START_REF] Allibert | Estimating body-fixed frame velocity and attitude from inertial measurements for a quadrotor vehicle[END_REF]. For instance, the velocity-aided solution allows to reconstruct the attitude efficiently with proven Lyapunov convergence. In this work, we similarly assume that the velocity in the sensor local frame is available. This velocity can be provided by a sensor such as Doppler effect radars. It can be provided also by the measurements of the gyrometer in the presence of a known anchor in the environment. This is for example the case of humanoid robots in contact with the environment, because the contact point position and velocity in the sensor frame are known [START_REF] Benallegue | Tilt estimator for 3D non-rigid pendulum based on a tri-axial accelerometer and gyrometer[END_REF]. There are a few estimators considering the same case. In [START_REF] Allibert | Velocity aided attitude estimation for aerial robotic vehicles using latent rotation scaling[END_REF] a velocity-aided estimator with proof of convergence has been presented, but there were possible cases of singularities if the scaling factor reaches zero. In [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] two estimators with proof of convergence have been presented, including one without gain condition, but we presented in [START_REF] Benallegue | Tilt estimator for 3D non-rigid pendulum based on a tri-axial accelerometer and gyrometer[END_REF] a slight improvement of the estimator proposed in [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] where better performance and simple convergence analysis are obtained. In [START_REF] Martin | A global exponential observer for velocity-aided attitude estimation[END_REF], a global estimator has been presented, but the globality has been reached at the cost of breaking the normality constraint of the gravity and magnetic field direction vectors. It presented also a projected rotation matrix on SO [START_REF] Benallegue | Tilt estimator for 3D non-rigid pendulum based on a tri-axial accelerometer and gyrometer[END_REF], which can be discontinuous in case of singularities.

We propose here to extend the work presented in [START_REF] Benallegue | Tilt estimator for 3D non-rigid pendulum based on a tri-axial accelerometer and gyrometer[END_REF] and at the same time the work of [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] and [START_REF] Martin | A global exponential observer for velocity-aided attitude estimation[END_REF] with a set of improved estimators. The choice of the estimator should be made mainly according to the availability and quality of the magnetometer measurements. The estimators include (i) an estimator for the tilt (roll and pitch angles) that does not use the magnetometer and (ii) a set of complete attitude estimators based on our confidence in the magnetometer. We first introduce a new tilt observer called a "two-step state observer" which operates in two steps: the first one provides an intermediate estimate x 2 ∈ R 3 of x 2 while the second step furnishes the recommended estimate x2 ∈ S 2 of x 2 based on x2 . The expected efficiency of this estimator is that its two constitutive steps are decoupled: the first one insures global exponential convergence of x 2 towards x 2 while the second one is an S 2 -constrained complementary-filter estimation also providing better robustness to disturbances. Finally, this estimator is used to obtain the complete rotation in the same way as the authors quoted. All full rotation estimators share a common simple structure, proof of asymptotic convergence and good overall performance. The quality of the estimate is evaluated by the comparative simulation with a state-of-the-art solution.

2. Problem statement 2.1. Frames and measurements. The problem we address is the estimation of the attitude of accelerated rigid body vehicles moving in 3D-space. We denote W the world frame and L the local frame of the sensor. This attitude estimation has to rely on an IMU consisting in a three-axial accelerometer, gyrometer and magnetometer, and using a measurement of the velocity of the sensor. The accelerometer provides y a the sum of the gravitational field and the linear acceleration of the sensor, the gyrometer provides y g measuring the angular velocity ω of the IMU and the magnetometer provides y m the measurement of the unit vector m along the Earth's magnetic field, all these signals are expressed in the sensor frame L. The velocity sensor provides y v the linear velocity v of the local frame L with respect to the world W, but expressed in L y v =v, (2.1)

y g =ω, (2.2)
y a =S(ω)v + v + g 0 R T e z , (2.3)
y m =R T m, (2.4)
R, g 0 , e z and m are respectively the orientation of the IMU with respect to the world, the standard gravity constant, a unit vector collinear with the gravitational field, expressed in W and directed upward, and a unit vector aligned with the earth's magnetic field expressed in W. Finally ω is the angular velocity of the sensor expressed in L such that

(2.5) Ṙ = RS(ω),
where the orientation R is the attitude we wish to estimate using these available measurements and the function S is the skew-symmetric matrix operator allowing to perform cross-product.

State definition and dynamics.

Let us define the following state variables:

x 1 ∆ = v, (2.6) x 2 ∆ = R T e z , (2.7) x 3 ∆ = R T m, (2.8)
where x 1 ∈ R 3 , x 2 ∈ S 2 and x 3 ∈ S 2 , with the set S 2 ⊂ R 3 being the unit sphere centered at the origin, and defined as (2.9)

S 2 ∆ = x ∈ R 3 / x = 1 .
The variables x 1 and x 3 are measured using y v and y m , even if they are noisy. On the contrary, x 2 is the tilt which cannot be obtained algebraically from the measurements.

From equations (2.3) and (2.6) we get

ẋ1 = -S(ω)x 1 + y a -g 0 R T e z . (2.10)
This, together with the time-differentiation of x 2 and x 3 using equation (2.5), provide us with the following state dynamic equations (2.11)

     ẋ1 = -S(ω)x 1 + y a -g 0 x 2 , ẋ2 = -S(ω)x 2 , ẋ3 = -S(ω)x 3 .
The system (2.11) is suitable for the observer synthesis.

In the next sections we show how to reconstruct the attitude in two cases:

• In Section 3: There is no magnetometer and only the tilt defined by x 2 = R T e z can be estimated.

• In Section 4: The magnetometer is available and the full rotation matrix R can be reconstructed. The presence of the magnetometer provides us with some redundancy but the magnetometer's measurements may be not reliable. This requires us to be able to tune how much the magnetometer interferes with tilt estimation.

2.3.

Basic facts and notation. We next introduce notations and recall basic properties used in the developments below where v, w and u are vectors and R ∈ SO(3) a rotation matrix where I denotes the 3 × 3 identity matrix.

S (v) S (w) = wv T -v T w I, (2.12) S (v) S (w) S(v) = -v T w S(v), (2.13) R n i=1 S (v i ) R T = n i=1 S (Rv i ) , (2.14) S (S (v) w) = S (v) S (w) -S (w) S (v) , (2.15) = wv T -vw T , (2.16) S 3 (v) = -v 2 S (v) , (2.
For

n ≥ 1, define Υ n ∆ = R 3n ×S ez and Υ * n ∆ = R 3n ×S * ez with S ez ∆ = z ∈ R 3 | (e z -z) ∈ S 2 and S * ez ∆ = S ez \ {2e z }.
Finally we will say that a dynamical system (D) ẋ = f (x) defined on a differential manifold X is almost globally asymptotically stable with respect to an equilibrium point x 0 of f if (D) is (Lyapunov) locally stable with respect to x 0 and there exists an open an dense subset X 1 of X such that every trajectory starting in X 1 converges asymptotically to x 0 .

Tilt estimation

In this section we ignore the signals of the magnetometer, either because it is unavailable or because the magnetic field is not steady. In this case the whole orientation cannot be observed. Nevertheless, we show hereinafter how we can have an efficient estimation of tilt x 2 = R T e z .

We present a novel observer called "Two-steps state observer" which is designed in two steps: the first step provides x 2 ∈ R 3 an intermediate estimate of x 2 ; and the second provides x2 ∈ S 2 the recommended estimate of x 2 based on the intermediate one x 2 . The expected efficiency of this estimator is that it relies on two stages, the first independent one is given by x 2 which is globally exponentially converging to x 2 in an efficient way, and the second is a S 2 -constrained complementary-filter estimation providing continuity and better robustness to disturbances. Indeed, the global exponential convergence of the error of the first stage actually leads to the violation of the normality constraint of R T e z . Furthermore, the simple normalization of this estimation may lead undefined output and unbounded time-derivatives when the norm is close to zero. This can cause a problem when continuity of the estimation is required. A simple solution is to add the second stage to maintain the constraint of the tilt estimation in S 2 while keeping bounded velocities.

3.1. Two-steps first order state observer designed in R 3 × S 2 . The simplest two-steps estimator can be described as follows, (3.1)

     ẋ1 = -S(y g )x 1 + y a -g 0 x 2 , x 2 = -α1 g0 (y v -x1 ) , ẋ2 = -S (y g -γS(x 2 )x 2 ) x2 ,
where α 1 and γ are positive scalar gains.

If the initial value of x2 is in S 2 , then the dynamics of the last equation ensures that the norm of this vector remains in time constantly equal one. The initial value for x1 on the other side could be anywhere in R 3 .

Using the estimation errors defined as

x 2 ∆ = x 2 - x 2 = x 2 + α1 g0 (x 1 -x1 ) = p 1 and x2 ∆ = x 2 -x2
, and equation (2.11) we get the error dynamics as

(3.2) ṗ1 = -S(ω)p 1 -α 1 p 1 , ẋ2 = -S(ω)x 2 + γS 2 (x 2 ) (x 2 -p 1 ) .
To run the analysis of errors, we set z p1 = Rp 1 and z 2 = Rx 2 . Noticing Rx 2 = e z -z 2 , one gets

(3.3) żp1 = -α 1 z p1 , ż2 = γS 2 (e z -z 2 ) (z 2 -z p1 ) .
This new error dynamics is autonomous and defines a time-invariant ordinary differential equation (ODE). If one defines the state ξ 1 ∆ = (z p1 , z 2 ) ∈ Υ 1 one can write (3.3) as ξ1 = F 1 (ξ 1 ) where F 1 gathers the right-hand side of (3.3) and defines a smooth vector field on Υ 1 .

We now turn to the convergence analysis of (3.3) and we get the following.

Theorem 1. The time-invariant ODE defined by (3.3) verifies the following (1) The state space is equal to Υ 1 , it admits two equilibrium points namely the origin (0, 0) and (0, 2e z ) and all trajectories of (3.3) converge to one of the two equilibrium points.

(a) The system (3.3) is almost globally asymptotically stable with respect to the origin, which is locally exponentially stable. (b) For every compact set K of Υ * 1 and positive number > 0, there exists (α 1 , γ) such that trajectories of (3.3) starting in K converge exponentially to the origin with an exponential rate larger than or equal to .

The proof of the theorem is given in Section 7.1.

Remark 2. The estimator for the tilt x 2 operates in two decoupled steps: the first one shows that the artificial state x 2 estimates x 2 (the dynamics of the error term z p1 is independent of the rest of the system dynamics) and then, in the second step, one brings back x 2 on S 2 through x2 .

We must now compare our tilt estimator with previous works. We start by considering the seminal work [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] and the tilt estimator introduced there, that we actually recall below in (3.8). The price to pay in the present paper with respect to (3.8) is the extra state x 2 but it has the advantage of not being constrained to S 2 anymore. This is clearly put forward when one compares the error dynamics given by (3.3) and (3.10). The decoupling in (3.3) between the errors z p1 and z 2 not only allows one to have better convergence results with respect to (3.8) (much simpler convergence analysis, no conditions on the gains and arbitrary rate of exponential convergence) but also to improve the robustness of the estimation to noise.

The second major reference for tilt estimation is that of [START_REF] Martin | A global exponential observer for velocity-aided attitude estimation[END_REF] where the authors provide an estimator in R 3 × R 3 instead of R 3 × S 2 (cf. the variables v and γ). The error system turns out to be (essentially) linear and time invariant with, therefore, the best convergence properties. However, the tilt estimator γ does not belong to S 2 and that may create singularity issues (i.e., γ may be equal to zero during the estimation or becomes collinear to another important unit vector) when one is interested by tilt estimation only or uses such an estimator to get an estimator of the total rotation R by, for instance, the TRIAD method [START_REF] Malcolm | Three-axis attitude determination from vector observations[END_REF].

In conclusion, our tilt estimator x2 combines the good convergence properties of the tilt estimator of [START_REF] Martin | A global exponential observer for velocity-aided attitude estimation[END_REF] with the fact that it remains on S 2 , like the tilt estimator of [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF].

3.2. Two-steps n th order state observer designed in R 3n ×S 2 . An interesting way to comprehend the estimator (3.1) is by noting that the dynamics of x 2 has a first order exponential convergence to x 2 and that the dynamics of x2 is a complementary filter of x 2 . Therefore, we can extend this feature to higher order of exponential convergence while keeping the same twosteps structure tilt estimator. Let n ≥ 2 be an integer. The n-th order observer is designed on R 3n × S 2 , where n is the order of the filter for the first step of the estimator. Increasing the order of linear filtering allows one to reduce the effect of the noises on the signals of the accelerometer y a and the linear velocity y v .

Define

p n ∆ = y v -x1 = x 1 -x1 .
Then the two-steps n-th order observer is given by

(3.4)          ẋ 2 = -S (y g ) x 2 -α1 g0 p 2 , ṗi = -S (y g ) p i + p i+1 , (i = 2, • • • , n -1) ẋ1 = -S (y g ) x1 + y a + n i=2 α i p i -g 0 x 2 , ẋ2 = -S (y g -γS (x 2 ) x 2 ) x2 .
Here, the gains α i , (i = 1, . . . , n) are positive and chosen so that the polynomial

s n + α n s n-1 + α n-1 s n-2 + α n-2 s n-3 + ... + α 2 s + α 1 is Hurwitz. Moreover,
x1 = x 1 -x1 = p n , x 2 = x 2 -x 2 = p 1 and x2 ∆ = x 2 -x2 , we get the error dynamics as (3.5)          ṗ1 = -S (ω) p 1 + α1 g0 p 2 , ṗi = -S (ω) p i + p i+1 , (i = 2, • • • , n -1) ṗn = -S (ω) p n -g 0 p 1 - n i=2 α i p i ẋ2 = -S (ω) x2 + γS 2 (x 2 ) x2 -γS 2 (x 2 ) p 1 .
To run the error analysis, we set z p1

∆ = Rp 1 , z pi ∆ = α1 g0 Rp i (i = 2, • • • , n) and z 2 ∆ = Rx 2 . Then one gets (3.6)      żpi = z pi+1 , (i = 1, • • • , n -1) żpn = - n i=1 α i z pi , ż2 = γS 2 (e z -z 2 ) (z 2 -z p1 ) .
This new dynamics is autonomous as for the first order case and defines a time-invariant ordinary differential equation (ODE) on Υ n . Similarly to the previous state estimator, if one defines the state

ξ n ∆ = (z p1 , • • • , z pn , z 2 ) ∈ Υ n , one can write (3.6) as ξn = F n (ξ n )
where F n gathers the right-hand side of (3.6) and defines a smooth vector field on Υ n .

Note that the first n lines of (3.6) constitute a separate tilt estimator defined in R 3n , which is similar in the case n = 2 to the one provided in [START_REF] Martin | A global exponential observer for velocity-aided attitude estimation[END_REF]. We show hereinafter the convergence and the performances of this estimation which are similar to the two-steps firstorder tilt estimator.

Theorem 3. The time-invariant ODE defined by (3.6) verifies the same statements as in Theorem1 up to changing Υ 1 and Υ * 1 by Υ n and Υ * n , the first zero in the equilibrium points now belonging to R 3n and α 1 changed by (α 1 , . . . , α n ).

The proof of the theorem is given in Section 7.2. The higher order of the convergence of x 2 allows one to improve the robustness of the estimation to noise by increasing the order n, as shown in simulations.

Remark 4. In the above construction, it is worth noticing that the gain-coefficients α i , 1 ≤ i ≤ n can be chosen time-varying and the above construction remains unchanged till (3.6). One can, therefore, use the observer given in [START_REF] Chitour | Time-varying high-gain observers for numerical differentiation[END_REF] which improves the performances of the complementary filter regarding a possible peaking phenomenon and noise.

3.3. One-step tilt observer designed in R 3 × S 2 . In this section, we provide a tilt observer first given in [START_REF] Benallegue | Tilt estimator for 3D non-rigid pendulum based on a tri-axial accelerometer and gyrometer[END_REF], which is a slight improvement of the estimator proposed in [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF]. It is designed in R 3 × S 2 in one step by using the available measures y g , y a and y v and is given by

(3.7) ẋ1 = -S(y g )x 1 -g 0 x2 + y a + αx 1 , ẋ2 = -S(y g + γS(x 2 )x 1 )x 2 ,
where α and γ are positive scalar gains which verify the condition γg 0 ≤ α 2 and x1 and x2 are the estimations of x 1 and x 2 .

The observer proposed in [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] is recalled as

(3.8) ẋ1 = -S(y g )x 1 -g 0 x2 + y a + k v 1 x1 -k v 2 S 2 (x 2 ) x1 , ẋ2 = -S(y g + k r 1 S(x 2 )x 1 )x 2 , where k v 1 , k v 2 and k r 1 are positive scalar gains which verify the condition k r 1 g 0 ≤ k v 1 k v 2 .
We can notice that our proposed observer can be obtained from the one of [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] by taking

k v 2 = 0, k v 1 = α, k r 1 = γ and the condition on the gains becomes k r 1 g 0 ≤ (k v 1 ) 2 instead of k r 1 g 0 ≤ k v 1 k v 2 . Using the errors x1 ∆ = x 1 -x 1 and x2 ∆ = x 2 -x 2
as well as ω = y g , a time-differentiation of these expressions provides us with the following error dynamics:

(3.9) ẋ1 = -S(ω)x 1 -αx 1 -g 0 x2 , ẋ2 = -S(ω)x 2 -γS 2 (x 2 )x 1 ,
The error dynamics of the observer of [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] is given by (3.10)

ẋ1 = -S(ω)x 1 -αx 1 -g 0 x2 + k v 2 S 2 (x 2 ) x1 , ẋ2 = -S(ω)x 2 -γS 2 (x 2 )x 1 ,
To run the analysis of errors, we define z i ∆ = Rx i . We notice also that R(x 2 + x2 ) = e z which leads to Rx 2 = e z -z 2 , we obtain this new error dynamics of our proposed observer

(3.11) ż1 = -αz 1 -g 0 z 2 , ż2 = -γS 2 (e z -z 2 )z 1 ,
We do the same for the observer of [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF], we get the error dynamics as

(3.12) ż1 = -αz 1 -g 0 z 2 + k v 2 S 2 (e z -z 2 ) z 1 , ż2 = -γS 2 (e z -z 2 ) z 1 .
These new error dynamics are autonomous. In fact, if we define the following state vector ξ ∆ = (z 1 , z 2 ) and the state space Υ 1 , we can write these errors dynamics as ξ = F (ξ) where F defines smooth vector fields on Υ 1 .

Almost global asymptotic stability with respect to the origin (0, 0) is obtained in the case of the proposed observer with the gain condition γg 0 ≤ α 2 and the proof is conducted in the same way as in [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] by considering only one Lyapunov function candidate given by

V ∆ = αz 1 + g 0 z 2 2 2 + g 2 0 z 2 2 2 . (3.13)
The time derivative of (3.13) in view of (3.11) yields

V = -α (1 -G 0 ) αz 1 + g 0 z 2 2 + αg 2 0 G 0 z T 2 S 2 (e z )z 2 -αG 0 (αz 1 + g 0 z 2 ) T (e z -z 2 ) 2 , (3.14)
where G 0 = γg0 α 2 ≤ 1. This Lyapunov function makes the convergence analysis much easier than that given in [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF].

Attitude estimation observer

The measurements of the magnetometer provide the direction of the magnetic field expressed in the local frame of the sensor. Usually most of the measurement is constituted with the earth natural magnetic field, which provides bi-dimensional data on the attitude of the sensor, providing then enough inputs to reconstruct the full attitude and having then some redundancy with the accelerometer for tilt estimation. However, sometimes due to the proximity of sources of interference, the magnetometer's measurements could lack the necessary reliability to let it influence the critical tilt estimation, but remains the best measurement available to reconstruct the orientation around the vertical direction. In this case the solution is to use an estimation allowing to tune the influence of the magnetometer on the tilt. 4.1. Design of the attitude observer. Let R ∈ SO(3) denote the estimate of R. The proposed non-linear observer takes advantage of the estimator of x 2 designed into R 3 given by (3.4) and the attitude estimator proposed by Mahony et. al [START_REF] Mahony | Nonlinear complementary filters on the special orthogonal group[END_REF], and it is given by (4.1)

       Ṙ = RS(y g -σ), σ = ρ 1 S( RT e z )x 2 + ρ 2 S( RT m)y m +µ RT e z RT e z T S RT m y m .
where ρ 1 , ρ 2 and µ are positive scalar gains and x 2 is given by the the first stage of any order of the two-step tilt estimator from Section 3, for example with (3.1).

In the case where ρ 2 = 0, we recover an estimate of the total rotation with decoupled tilt in an essentially similar way as that of [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] where the magnetometer has no influence on the tilt. On the contrary, if µ = 0, the corresponding estimator is closer to that of [START_REF] Mahony | Nonlinear complementary filters on the special orthogonal group[END_REF] and the estimator fully uses the redundancy.

Let R = R RT be the attitude estimation error. A time-differentiation of the expression of (4.1) and the use of equation ( 3.3) provides us with the following error dynamics:

Ṙ = RS (σ) , (4.2)
where σ is given by

σ = I + µ ρ2 e z e T z ρ 1 S (e z ) RT e z + ρ 2 S (m) RT m -ρ 1 S(e z ) RT z p1 . (4.3)
Using unit-quaternions instead of elements of SO(3), one associates Q and Q with the rotations R and R respectively, and similarly the unit-quaternion error Q = (q 0 , q) = Q Q-1 with the attitude estimation error R. Here, q0 ∈ R and q ∈ R 3 are the scalar and the vector components of Q respectively. We can, therefore, write Lemma 2). Set ∆ = (q 0 I -S (q)) W ρ q. The error dynamics written as a quaternion error dynamics is now given by (4.5)

R = I + 2q 0 S (q) + 2S 2 (q) , ρ 1 S (e z ) RT e z + ρ 2 S (m) RT m = -2 (q 0 I -S (q)) W ρ q, (4.4) with W ρ ∆ = -ρ 1 S 2 (e z ) -ρ 2 S 2 (m) being a positive- definite symmetric matrix ([20],
           q0 = qT I + µ ρ2 e z e T z , + 1
2 ρ 1 qT S(e z ) I -2q 0 S (q) + 2S 2 (q) z p1 , q = -(q 0 I + S (q)) I + µ ρ2 e z e T z -1 2 ρ 1 (q 0 I + S (q)) S(e z ) I -2q 0 S (q) + 2S 2 (q) z p1 . The above equation together with the first equation of (3.3) define a time-invariant ordinary differential equation (ODE) and, by considering the state ξ ∆ = z p1 , Q and the state space Υ ∆ = R 3 × S 3 , one can write (3.3) and (4.5) as ξ = F (ξ) where F gathers the right-hand side of (4.5) and defines a smooth vector field on Υ . We analyze this dynamics in the next section.

Stability analysis. Let us consider the following positive-definite differentiable function

V ∆ = ρ 2 1 α 1 z p1 2 + 2q T W ρ q, (4.6)
which is clearly radially unbounded.

Theorem 5. The time-invariant ODE defined by (4.5) verifies the following.

(1) Its equilibrium points are

Ω 1 = {(0, (±1, 0))} , Ω 2 = {(0, (0, ±v jρ )), j = 1, 2, 3} ,
where v jρ are unit eigenvectors of W ρ for 1 ≤ j ≤ 3.

(2) All trajectories of (4.5) converge to one of the equilibrium points defined in item 1. (3) The set equilibrium Ω 1 which corresponds to the equilibrium point z p1 = 0, R = I is asymptotically stable with a domain of attraction containing the domain (4.7)

V c ∆ = ξ = z p1 , Q ∈ Υ | V (ξ) < 2λ min (W ρ ) .
(4) The equilibria of the set Ω 2 are unstable and the system is almost globally asymptotically stable with regard to Ω 1 . The proof is given in Section 7.3. Remark 6. The magnetic field measurements y m can also be filtered using an additional unconstrained state on the unit sphere in the same way as done for the tilt in order to improve robustness to noise. Remark 7. In the above, we have chosen, for the simplicity of the analysis, to estimate the intermediate state x 2 with the two-steps first order state observer given by (3.1). One can also rely on the two-steps n th order state observer given in (3.6). For the corresponding stability analysis, one replaces the z p1 2 term in the Lyapunov function V given in (4.6) by ψ T n P α ψ n given in (7.2).

Simulations

We show hereinafter results of the estimators in a simulated environment.

Signal generation and initialization.

In this section, we present simulation results showing the effectiveness of the proposed estimators. We generated the signal ω and v with trigonometric functions and generated the trajectory of R by integration (see Figure 5.1), then we simulated the signals of the accelerometer y a , the gyrometer y g and the magnetometer y m such that

m = 1 √ 2 (1,
= 20, α 1 = 2 √ γg 0 2nd order γ = 20, α 1 = 2 √ γg 0 2 , α 2 = 2 2 √ γg 0 3rd order γ = 20, α 1 = 2 √ γg 0 3 , α 2 = 3 2 √ γg 0 2 , α 3 = 3 2 √ γg 0 Table 2.
Parameters of tested tilt estimators were added to the four measurements, the accelerometer y a the gyrometer y g , the velocity sensor y v , and the normalized magnetometer y m to which a stronger noise and a bias have been added to make it unreliable and unsuitable to influence tilt estimation. The detail of the noise properties is summarized in Table 1.

For each tested estimator the initial state was set to R3 = 2 m×ez m×ez m×ez m×ez T -I , which corresponds to an undesired equilibrium. The velocity estimation was initialized to the current sensor value (for instance x1 (0) = x 1 (0)).

Comparison between two stage tilt estimators.

The first test is to compare the tilt estimators presented in Section 3. Specifically, the first order, the second order and the third order tilt estimators were compared for the perfect and the noisy measurements. The estimators were designed to have the same (multiple) pole. The parameters are detailed in Table 2.

The result of the simulation with perfect measurements is shown in Figure 5.2 where we compare the errors produced by the estimations x2 but also the intermediate estimations x 2 . We see that the intermediate estimation errors converge exponentially to zero while the estimation itself remains in the undesired equilibrium. We can see that the first order estimator is obviously the fastest followed by the other orders.

However, the more interesting case of the noisy one displayed in Figure (5.3). We see then that with higher orders of the estimator better filtering is provided. We see also that the sphere constraint of the final estimate x2 allows to reduce the noise by removing the components which are orthogonal to the constraints. Nevertheless the difference between the second and the third order is small enough to consider that the second order estimator is a good trade-off between complexity and speed on one side and filtering quality on the other. Therefore, in the following simulations we will use to feed the attitude estimator in (4.1) with x 2 and then compare it with state-of-the-art approaches. (1) The attitude estimator described in Section 4, with ρ 2 = 0 which means it is decoupled to avoid any impact of the magnetometer on the tilt estimation. We refer to it by hierarchic. (2) The estimator described in Section 4, with µ = 0 using redundancy, referred to as Invariant. (3) The estimator in R 3 × S 2 in Section 3.3, providing only tilt estimation, and referred to as Benallegue 2017. (4) The estimator designed by Hua et al, named "Observer 2" in [START_REF] Hua | Stability analysis of velocity-aided attitude observers for accelerated vehicles[END_REF] reported in (3.8), that we refer to as Hua 2016. (5) The estimator described in the preprint [START_REF] Martin | A global exponential observer for velocity-aided attitude estimation[END_REF] by Martin et al, named Martin 2016, which is Plot showing the estimation error, in the case of noisy measurements, of the tilt vector x 2 for the three orders of tilt estimators. For each order we show both the intermediate estimation x 2 and the final one x2 . The three top images show the three components of the vector difference error z 2 = R (x 2 -x2 ), the 4th plot shows the evolution of the angle between the tilt x 2 and its estimation x2 , and the bottom part shows an enlarged plot of the second [START_REF] Allibert | Estimating body-fixed frame velocity and attitude from inertial measurements for a quadrotor vehicle[END_REF][START_REF] Allibert | Velocity aided attitude estimation for aerial robotic vehicles using latent rotation scaling[END_REF] of the angle error.

based on an estimator equivalent to 2nd order x 2 of Sec. 3.2 and another exponential estimator of the tilt and using TRIAD [START_REF] Malcolm | Three-axis attitude determination from vector observations[END_REF] to reconstruct the attitude. The estimation is designed for the tilt to depend only on the accelerometer and the yaw angle only on the magnetometer. shows the evolution of the angle between the tilt x 2 and its estimation x2 . Note that we use the names of the second column of Table 3.

Each estimator provides a specific tilt estimator. Note that the tilt estimation of the invariant observer is the only one that requires magnetometer's measurements. The corresponding gains were designed to have the most equivalent behavior possible, regarding ture and the considered errors. These estimators, as well as their tilt component and the gains used, are summarized in Table 3 1 . 5.3.1. Perfect measurements. Figure 5.4 shows the evolution of the tilt error for the five tilt estimators. The first estimator to converge is the one of R 3 ×R 3 , namely 2nd order x 2 of Sec. 3.2 which is not constrained to the unit sphere, this is because the starting position is not an equilibrium point for this vector. However, the normalization of this vector gives a discontinuous trajectory visible at the bottom plot showing the angle 1 The tilt estimation column relates the different estimators to their tilt estimation component and the gains column gives the gain values adopting the notation used in each corresponding cited document. error. This is the estimation used in [START_REF] Martin | A global exponential observer for velocity-aided attitude estimation[END_REF]. The next estimator to converge is the invariant one, this is due to the fact that this estimator uses also the measurement of the magnetometer to speedup the convergence. After that the estimator in R 3 × R 3 × S 2 of Section 3.2 is the next to quickly converge, while staying continuous and constrained on the unit sphere. The other estimators converge later, especially the estimation of Hua 2016.

Figure 5.5 shows the evolution over time of the estimation error of the vector m p ∆ = e z × x 3 × e z which is orthogonal to e z but pointing at the same horizontal direction as x 3 . The error is shown as an angle which can be interpreted as a "yaw angle error" when the tilt error is small. In this figure we see that the estimation of Martin 2016 [START_REF] Martin | A global exponential observer for velocity-aided attitude estimation[END_REF], is discontinuous at another instant than the discontinuity of the tilt, which means that the attitude had two discontinuities while converging. The invariant estimator converges fast, taking full profit from the redundancy. The Hierarchic, was the next estimator to converge. We see finally that the estimation error of Hua 2016 moved at second 2 to zero. However, this does not correspond to the convergence of the estimator since it took the tilt estimation 4 more seconds to converge (see Figure 5.4). This means that it only went from an undesired equilibrium to another one. Note that some angles increase and then decrease, and this happens because of the tilt estimation converging at the same time and the orthogonality constraint being respected. From this plot we see that the noise allowed the estimators to instantly leave the repulsive undesired equilibrium. Then most estimators except for Hua 2016 converge in less than half a second, the unconstrained R 3 × R 3 being the fastest. After the convergence of all the estimators, we see in the enlarged plot that the estimator in R 3 × R 3 × S 2 has the lowest tilt estimation error angle. The dynamics of the estimators in R 3 × R 3 and R 3 × S 2 have identical local behavior near the desired equilibrium and are almost superimposed in the steady behavior. Interestingly the invariant observer gives worse estimations, that is because it involved the unreliable magnetometer measurements which downgrade the performances. We see in table 4 the mean value of the tilt error angle over 8 seconds after the second 2 of the simulation. The constrained R × R 3 × S 2 the best estimations and the invariant and Hua 2016 both give the worst ones.

k v 1 = k v 2 = α, k r 1 = k r 2 = γ Martin 2016 [16] R 3 × R 3 (2nd order x 2 of Sec. 3.2) L = K = α 2 , M = µ
Figure 5.7 shows the evolution of the estimation of R T m p together with a zoom on the 6-th second of the simulation. We see that the estimations converge in the first second except for Hua 2016. The high level of noise in the magnetometer produces a poor estimation quality, but in the steady behavior a difference can be shown between observers. This can be quantitatively assessed by looking at Table 5 showing the average error angle values in the interval [2s,10s]. Martin 2016 has low quality estimations because the estimation of yaw is performed independently from the measurements of the accelerometer. The other estimators take profit from the better reliability of the tilt estimation and provide a relatively similar level of performance with a slight advantage to the hierarchic estimator.

Note that a behavior between the hierarchic and the invariant estimator can be obtained by choosing values of µ and ρ 2 appropriately, especially that small values of ρ 2 provide better theoretical convergence guarantees without downgrading excessively the quality of the estimation.

Conclusion

We have presented a set of attitude estimators using the measurements of an accelerometer, a gyrometer, a magnetometer, and a linear velocity expressed in the local frame. These estimators are intended to be used in specific cases, mostly related to the reliability of the magnetometer for tilt estimation. Indeed, the magnetometer can be either reliable, unreliable, or totally unavailable. For instance, an invariant complementary filter has good performances when the magnetometer is reliable but is disturbed when it is not. Among the estimators, we developed a second order complementary filter for the tilt, and we augmented it with an attitude estimator allowing to tune how much we use the magnetometer measurement in the tilt estimation. We have assessed the performances of these estimators through simulations of perfect and noisy measurements. 3) leaves invariant Υ 1 since along its trajectories, e z -z 2 keeps a constant norm equal to one. Moreover, it admits two equilibrium points (0, 0) and (0, 2e z ). Let us now consider the following positive-definite differentiable function

V 1 : Υ 1 → R + V 1 ∆ = 1 2α 1 z T p1 z p1 + 1 2γ z T 2 z 2 (7.1) then the time derivative of V 1 is given by V1 = -z p1 2 + z T 2 S 2 (e z -z 2 ) z 2 -z T 2 S 2 (e z -z 2 ) z p1
If we use w = S(e z )z 2 , we can write

V1 ≤ - z p1 w 1 -1 2 -1 2 1 z p1 w ≤ 0,
and V1 < 0 if (z p1 , w) = 0. This is equivalent to (z p1 , z 2 ) not being an equilibrium point, i.e., all trajectories of (3.3) converge to one of the two equilibrium points defined previously.

2) The linearized system at (0, 0) is Hurwitz yielding that (0, 0) is locally exponentially stable. At (0, 2e z ), the linearized system clearly admits two positive real eigenvalues. We can conclude that the system (3.3) is almost globally asymptotically stable with respect to the origin (0, 0). Moreover, the set of points of Υ 1 for which V 1 has values less than V 1 (0, 2e z ) = 2/γ is clearly included in the basin of attraction of (0, 0).

3) Let K be a compact set in Υ * 1 and > 0. It is easy Plot showing the estimation error of the tilt vector x 2 for the five tilt estimators while the measurements were noisy. The three top images show the three components of the vector difference error z 2 = R (x 2 -x2 ), the 4th plot shows the evolution of the angle between the tilt x 2 and its estimation x2 , and the bottom part shows an enlarged plot of the 6th second of the angle error. Note that we use the names of the second column of Table 3. to see that for every (z p1 , z 2 ) ∈ K one has that

V 1 (z p1 , z 2 ) ≤ C 1 (K) 2α 1 + 2 γ C 2 (K), with C 2 (K) < 1.
By the remark at the end of the argument of item 2) and by choosing α 1 large enough, one gets K is in the basin of attraction of (0, 0).

To obtain the last statement, we first prove that there exists C 3 (K) > 0 such that w ≤ C 3 (K) z 2 for (z p1 , z 2 ) ∈ K. It is enough to check that for z 2 small enough. If one writes z 2 = (e T z z 2 )e z + z ⊥ 2 one gets that z 2 ≤ 2 z ⊥ 2 (by using that z 2 ∈ S ez ) and w = z ⊥ 2 for z 2 small, hence the claim. Next, one deduces that there exists

C(K) > 0 such that V1 ≤ -C(K)( z p1 2 + z 2 2 ),
for trajectories starting in K (and staying in a compact neighborhood of K in the basin of attraction of (0, 0)).

Then V1 ≤ -2C(K) min(α 1 , γ)V 1 . By taking α 1 and γ large enough, one gets the conclusion.

7.2. Proof of Theorem 3. The argument is similar to that of Theorem 1. For that purpose consider the Hurwitz n × n matrix in companion form A α = J n -ae T n , J n stands for the n-th Jordan block, a = (α 1 , . . . , α n ) T and e n = (0, • • • , 0, 1) T . Then set M α = A α ⊗I 3×3 and ψ n = (z p1 , . . . , z pn ) ∈ R 3n . Note that the n first equations in (3.6) can be written ψn = M α ψ n . Let P α be the positive definite real symmetric matrix, unique solution of the Lyapunov equation M T α P α + P α M α = -I 3n×3n . Recall that P α ≤ C reα , where C n is a universal positive constant and re α > 0 is the minimum of -Re(λ), Re stands for the real part and λ is any eigenvalue of the A α , cf. [START_REF] Horn | Topics in Matrix Analysis[END_REF]. One now considers the Lyapunov function (7.2) V n = ψ T n P α ψ n + 1 2γ z T 2 z 2 .

We now follow exactly the argument of Theorem 1 and replace z p1 by ψ n to get the conclusion.

7.3. Proof of Theorem 5. Let us prove the four items of the theorem.

1) The equilibria are calculated by solving the equation ξ = 0. The solutions of this equation system are given by (z p1 = 0, = 0). We know from ( [START_REF] Tayebi | Inertial vector measurements based velocity-free attitude stabilization[END_REF], Lemma 3) that = 0 is equivalent to (q 0 = ±1, q = 0) or (q 0 = 0, q = ±v ρ ) where v ρ is one of the unit eigenvectors of W ρ . This completes the proof of item 1.

2) Using the error dynamics given by (4.5), the time derivative of V is then given by (7.3) V = -2ρ 2 1 z p1 2 + 4q T W ρ q, which can be developed into V = -2ρ 2 1 z p1 2 (7.4) -4q T W ρ (q 0 I + S (q)) I + µ ρ 2 e z e T z + 2ρ 1 qT W ρ (q 0 I + S (q)) S(e z ) RT z p1 . (7.5) Using the definition of the vector , we obtain The right-hand side of the above inequality is a quadratic form in ( z p1 ,

V = -2ρ
) which is clearly negative definite. One easily verifies that V < 0 if ξ = z p1 , Q is not an equilibrium. Since (4.5) is autonomous and V is radially unbounded, one can use Lasalle's invariance theorem. Therefore, every trajectory converges asymptotically to a trajectory along which V ≡ 0. Since V is non-increasing, V (ξ) < 2λ min (W ρ ) at t = 0, implies that q(t) < 1 for every t ≥ 0. Since the trajectory converges to one of the equilibrium points, it must be one with (z p1 = 0, q = 0) which corresponds to Ω 1 because this is the only one contained in V c .

4) The undesired equilibria characterized by q0 = 0 are given by X = (0, (0, v ρ )). Let us show that X = (z p1 = 0, (q 0 = 0, q = v ρ )) is unstable. The linearized error dynamics around the unstable equilibrium X = (0, (0, v ρ )) is given by (7.7) ξ = Aξ, with A given in Equation (7.8) in the box on the following page. It is clear that there is at least one positive eigenvalue of the matrix A. Thus, there exists an unstable manifold of dimension at least one in neighborhoods of the Ω 2 = {(0, (0, ±v jρ )), j = 1, 2, 3}, and since all trajectories converge to an equilibrium point, then (4.5) is almost globally asymptotically stable with respect to the two equilibrium points Ω 1 = {(0, (±1, 0))} which correspond to z p1 = 0, R = I . This completes the proof.

Let us consider v ⊥ ρ = S (v ρ ) e z , we can write the matrix A as (7.8)

A =         -α 1 I 0 0 -1 2 ρ 1 v ⊥T ρ λ ρ 1 + µ ρ2 v T ρ e z 2 -µ ρ2 v T ρ e z v ⊥ ρ T (λ ρ I -W ρ ) 1 2 ρ 1 S (v ρ ) S (e z ) I -2v ρ v T ρ -λ ρ µ ρ2 e T z v ρ v ⊥ ρ I + µ ρ2 v ⊥ ρ v ⊥T ρ (λ ρ I -W ρ )        

17

 17

  )

Figure 5 . 1 .

 51 Figure 5.1. Plot showing the real state of the system. On the top the orientation is shown in Euler angles and on the bottom, the velocity x 1 is shown in its three components.

Figure 5 . 2 . 5 . 3 .

 5253 Figure 5.2. Plot showing the estimation error of the tilt vector x 2 for the three orders of tilt estimators. For each order we show both the intermediate estimation x 2 and the final one x2 . The three top images show the three components of the vector difference error z 2 = R (x 2 -x2 ) and the bottom plot shows the evolution of the angle between the tilt x 2 and its estimation x2 .

Figure 5 . 3 .

 53 Figure 5.3. Plot showing the estimation error, in the case of noisy measurements, of the tilt vector x 2 for the three orders of tilt estimators. For each order we show both the intermediate estimation x 2 and the final one x2 . The three top images show the three components of the vector difference error z 2 = R (x 2 -x2 ), the 4th plot shows the evolution of the angle between the tilt x 2 and its estimation x2 , and the bottom part shows an enlarged plot of the second[START_REF] Allibert | Estimating body-fixed frame velocity and attitude from inertial measurements for a quadrotor vehicle[END_REF][START_REF] Allibert | Velocity aided attitude estimation for aerial robotic vehicles using latent rotation scaling[END_REF] of the angle error.

Figure 5 . 4 .

 54 Figure 5.4. Plot showing the estimation error of the tilt vector x 2 for the five tilt estimators. The three top images show the three components of the vector difference error z 2 = R (x 2 -x2 ) and the bottom plot shows the evolution of the angle between the tilt x 2 and its estimation x2 . Note that we use the names of the second column of Table3.

Figure 5 . 5 .

 55 Figure 5.5. Plot showing the evolution of the angle between R T m p and its estimation RT m p for the five full attitude estimators.

5. 3 . 2 .

 32 Noisy Measurements. Figure5.6 shows the tilt estimation error with the difference and the angle, similarly to Figure5.4, with an additional enlarged sample plot of the behavior after the convergence.

  (4.1) with ρ 2 = 0) R 3 × R 3 × S 2 (2nd order x2 of Sec. 3.2) γ = 20, α = 2 √ γg 0 = 28.0143, Invariant (Observer (4.1) with µ = 0) Invariant (uses magnetometer) µ = 20 (hierarchic) or ρ 2 = 20 (invaria Benallegue 2017 [3] of Sec. 3.3 R 3 × S 2 , (provides tilt only) ρ 1 = 20, α 1 = γg 0 , α 2 = 2 √ γg 0 Hua 2016 [10] Hua 2016 (Observer (3.8))

7 . Appendix 7 . 1 .

 771 Proof of Theorem 1. 1) One easily checks that the time-invariant ODE defined by(3.

Figure 5 . 6 .

 56 Figure 5.6.Plot showing the estimation error of the tilt vector x 2 for the five tilt estimators while the measurements were noisy. The three top images show the three components of the vector difference error z 2 = R (x 2 -x2 ), the 4th plot shows the evolution of the angle between the tilt x 2 and its estimation x2 , and the bottom part shows an enlarged plot of the 6th second of the angle error. Note that we use the names of the second column of Table3.

Figure 5 . 7 .

 57 Figure 5.7. Plot showing the evolution of the angle between R T m p and its estimation RT m p for the five full attitude estimators under noisy measurements.

  x1 and x2 are estimations of x 1 and x 2 respectively and x 2 is an intermediate estimation of x 2 . Using the estimation errors defined as

Table 1 .

 1 0, 1) T . These signals were used in two cases, ideal signals and noisy ones. For the noisy signals, Gaussian noises Description of the noise parameters.

	Measurement Noise std.	Bias
	Accelero. y a 0.31 m/s	2	(0 0 0) T
	Gyro. y g	0.1 rad/s	(0 0 0) T
	Magneto. y m	0.71		(0.2 0.2 0.2)	T
	Velocity y v	0.31m/s	(0 0 0)
	Order	Parameters
	1st order	γ

T

Table 3 .

 3 Summary of compared estimators in the simulations.

	Tilt estimation Mean tilt error angle
	R 3 × R 3 × S 2	0.0442 rad
	Invariant	0.1543 rad
	R 3 × S 2	0.0748 rad
	Hua 2016	0.1960 rad
	R 3 × R 3	0.0749 rad

Table 4 .

 4 Average tilt error angles during 8 seconds after the convergence of the estimators.

	Attitude est. Mean m p error angle
	Hierarchic	0.2374 rad
	Invariant	0.2511 rad
	Hua 2016	0.2671 rad
	Martin 2016	0.3036 rad

Table 5 .

 5 Average R T m p estimation error angles during 8 seconds after the convergence of the estimators.


Marie Guionie

Aurélien Thorette

Marco Romanelli

Emmanuel Pinsard

Laurent Lablonde

Benoît Cadier

Marc Vallet

Marc Brunel

Dual-polarization DFB fiber laser stabilized by frequency-shifted feedback

teaching and research institutions in France or abroad, or from public or private research centers.

Dual-polarization DFB fiber laser stabilized by frequency-shifted feedback

Marie Guionie 1 , Aurélien Thorette Distributed-feedback (DFB) fiber lasers can sustain the oscillation of two orthogonal polarizations at different frequencies, due to the fiber Bragg grating UV-induced birefringence [START_REF] Loh | 1.55 µm phase-shifted distributed feedback fibre laser[END_REF]. Such dual-frequency fiber lasers (DFFL) are promising as heterodyne sources but, in the context of microwave photonics, the stabilization of the beat frequency becomes necessary. Locking the beat frequency of a dual-frequency laser to a local oscillator by optical frequency-shifted feedback (FSF) has been proved to be efficient [START_REF] Kervevan | Beat-note jitter suppression in a dual-frequency laser using optical feed-back[END_REF], and can be applied in principle to any dual-polarization laser. Here we investigate this FSF method to an Er 3+ -doped DFB fiber laser, resulting in an all-fibered 1.5μm system emitting a stabilized beat note at 1 GHz.

The set-up is depicted in Fig. 1(a). The DFFL is pumped at 976 nm through the WDM and emits at 1547 nm two orthogonal polarizations (with associated frequencies ν x and ν y ). The free-running beat note is around 1 GHz; its line-width is equal to 3 kHz and presents a drift of a few MHz over hours [3]. The polarizationmaintaining FSF loop is based on a polarization beam splitter whose output ports are closed on one another after passing through the frequency shifter, here an intensity modulator EOM driven at f LO , and an optical isolator. The polarization controller PC permits to select one of the polarization states, say ν x (ν x < ν y ), which is reinjected into the laser. The EOM thus generates two sidebands at frequencies ν x ± f LO . Hence, for a small detuning Δν = f LO -(ν y -ν x ), the re-injected field contains an optical sideband resonant with ν y . Intracavity coupling can then lead to frequency locking between ν y and the ν x +f LO component of the reinjected field. This results in the transfer of the spectral purity from f LO to the beat frequency (ν y -ν x ). The observed dynamics depend on both detuning Δν and injection strength Γ. First, we find a frequency locking between the beat note and the synthesizer, delimited by the experimental blue points (see Fig. 1(b)), over a 4 MHz range at maximum injection strength. Then, at larger excitation rates, the laser may display a complicated, possibly chaotic dynamics on the Δν > 0 side (delimited by the red points). In the locking zones, the stabilized beat frequency exhibits a narrow linewidth, smaller than 1 Hz (see Fig. 1(c)). The associated phase noise is measured to be -100 dBc/Hz at 1 kHz from carrier. The laser remains phase-locked for days in laboratory conditions.

Such stabilization scheme based on an EOM is advantageous in terms of carrier tuning, compactness and integration, and could find application in, e.g., radio-over-fiber. To this aim, similar studies on DFFL at higher frequencies are under investigation. Finally, it is interesting to note that the linewidth enhancement factor of such lasers could be measured using this stabilization scheme [START_REF] Thorette | Linewidth enhancement factor measurement based on FM-modulated optical injection: application to rare-earth-doped active medium[END_REF].

Fig. 1 (

 1 Fig. 1 (a) Experimental scheme. The optional EDFA permits here to control the re-injected power. (b) Normalized injection strength Γ with respect to the detuning Δν at two different excitation ratios η. Colored areas delimit the different regimes. White: unlocked, blue: locked, red: intensity modulated. (c) RF Spectrum of the locked DFFL beat.

  1 , Marco Romanelli 1 , Emmanuel Pinsard 2 , Laurent Lablonde 2 , Benoît Cadier 2 , Marc Vallet 1 , Marc Brunel 1

	1. Univ Rennes, CNRS, Institut FOTON -UMR 6082, F-35000 Rennes, France
	2. iXBlue Photonics, Rue Paul Sabatier, 22300 Lannion, France