Sébastien Briot
email: sebastien.briot@irccyn.ec-nantes.fr

Vigen Arakelian
email: vigen.arakelyan@insa-rennes.frjean-paul.le-baron@insa-rennes.fr

Jean-Paul Le Baron

Shaking Force Minimization of High-Speed Robots via Centre of Mass Acceleration Control

Keywords: Shaking force, balancing, high-speed robots, optimal motion planning

This paper deals with the problem of shaking force balancing of high-speed manipulators. The known solutions of this problem are carried out by an optimal redistribution of moving masses which allows the cancellation or the reduction of the variable loads on the manipulator frame. In this paper an innovative solution is developed which is based on the optimal control of the robot links centre of masses. Such a solution allows the reduction of the acceleration of the total mass centre of moving links and, consequently, the considerable reduction in the shaking forces. The efficiency of the suggested method is illustrated by the numerical simulations carried out for different trajectories: for examined planar two and three links serial manipulators the shaking force reduction reaches up to 77%. This approach is also a more appealing alternative to conventional balancing methods because it allows the reduction of the shaking force without counterweights. As a result, the input torques are also decreased, which is shown using dynamic simulation software.

Introduction

A mechanical system with unbalance shaking force/moment transmits substantial vibration to the frame. Thus, a primary objective of the balancing is to cancel or reduce the variable dynamic loads transmitted to the frame and surrounding structures. Different approaches and solutions devoted to this problem have been developed and documented for one degree of freedom mechanisms [START_REF] Lowen | Balancing of linkagesan Update[END_REF], [START_REF] Arakelian | Shaking force and shaking moment balancing of mechanisms: a historical review with new examples[END_REF]. A new field for their applications is the design of mechanical systems for fast manipulation, which is a typical problem in advanced robotics.

The balancing of a mechanism is generally carried out by two steps: (i) the cancellation (or reduction) of the shaking force and (ii) the cancellation (or reduction) of the shaking moment. Traditionally, the cancellation of the shaking force transmitted to the manipulator frame can be achieved via adding counterweights in order to keep the total centre of mass of moving links stationary [START_REF] Lowen | Balancing of linkagesan Update[END_REF], via additional structures [START_REF] Lowen | Balancing of linkagesan Update[END_REF], [START_REF] Van Der Wijk | Synthesis method for linkages with centre of mass at invariant link point -Pantograph based mechanisms[END_REF] or by elastic components [START_REF] Alici | Optimum force balancing with mass distribution and a single elastic element for a five-bar parallel manipulator[END_REF].

With regard to the shaking moment balancing of manipulators, the following approaches were developed: (i) balancing by counter-rotations [START_REF] Berkof | Complete force and moment balancing of inline four-bar linkages[END_REF]- [START_REF] Van Der Wijk | Comparative analysis for low-mass and low-inertia dynamic balancing of mechanisms[END_REF], (ii) balancing by adding four-bar linkages [START_REF] Ricard | On the design of reactionless parallel manipulators[END_REF]- [START_REF] Jiang | Dynamic optimization of reactionless four-bar linkages[END_REF], (iii) balancing by creating redundant mechanism which generates optimal trajectories of moving links [START_REF] Papadopoulos | Design and motion planning for a zeroreaction manipulator[END_REF]- [START_REF] Nenchev | Reaction nullspace control of flexible structure mounted manipulator systems[END_REF] (vi) balancing by prescribed rotation of the end-effector [START_REF] Fattah | On the design of reactionless 3-DOF planar parallel mechanisms[END_REF]- [START_REF] Briot | Complete shaking force and shaking moment balancing of the position-orientation decoupled PAMINSA manipulator[END_REF] and (vii) balancing by adding an inertia flywheel rotating with a prescribed angular velocity [START_REF] Arakelian | Design of planar 3-DOF 3-RRR reactionless parallel manipulators[END_REF], [START_REF] Van Der Wijk | Active dynamic balancing unit for controlled shaking force and shaking moment balancing[END_REF].

In the present paper we consider a simple and effective balancing method, which allows the considerable reduction of the shaking force of non-redundant manipulators without adding counterweights. It is based on the optimal control of the acceleration of the total mass centre of moving links. To the best of the authors' knowledge, this problem is addressed for the first time.

Minimization of the Shaking Forces via an Optimal Motion Planning of the Total Mass

Centre of Moving Links

Definition of the optimal trajectory

The shaking forces f sh of a manipulator can be written in the form:

sh iS m   fx (1)
where  i m is the total mass of the moving links of the manipulator and S x is the acceleration of the total mass centre. The classical balancing approach consists in adding counterweights in order to keep the total mass centre of moving links stationary. In this case, S x = 0 for any configuration of the mechanical system. But, as a consequence, the total mass of the manipulator is considerably increased. Thus, in order to avoid this drawback, in the present study, a new approach is proposed, which consists of the optimal control of the total mass centre of moving links. Such an optimal motion planning allows the reduction of the total mass centre acceleration and, consequently, the reduction of the shaking force.

Classically, manipulator displacements are defined considering either articular coordinates q or Cartesian variables x. Knowing the initial and final manipulator configurations at time t 0 and t f , denoted as q 0 = q(t 0 ) and q f = q(t f ), or x 0 = x(t 0 ) and x f = x(t f ) , in the case of the control of the Cartesian variables, the classical displacement law may be written in the form:

  ( ) ( ) q t s t    f 0 0 q q q q (2a) or   ( ) ( ) x t s t    f 0 0 x x x x ( 2b
)
where s q (t) and s x (t) may be polynomial (of orders 3, 5 and higher), sinusoidal, bang-bang, etc. motion profiles [START_REF] Khalil | Modeling, identification and control of robots[END_REF].

From expression (1), we can see that the shaking force, in terms of norm, is minimized if the norm S

x of the masses centre acceleration is minimized along the trajectory. This means that if the displacement x S of the manipulator centre of masses is optimally controlled, the shaking force will be minimized. As a result, the first problem is to define the optimal trajectory for the displacement x S of the manipulator centre of masses.

For this purpose, let us consider the displacement x S of a point S in the Cartesian space. First, in order to minimize the masses centre acceleration, the length of the path followed by S should be minimized, i.e. point S should move along a straight line passing through its initial and final positions, denoted as x S0 and x Sf , respectively.

Then, the motion profile used on this path should be optimized. It is assumed that, at any moment during the displacement, the norm of the maximal admissible acceleration the point S can reach is constant and denoted as max

S

x . Taking this maximal value for the acceleration into consideration, it is known that the motion profile that minimize the time interval (t 0 , t f ) for going from position x S0 = x S (t 0 ) to position x Sf = x S (t f ) is the "bang-bang" profile [START_REF] Khalil | Modeling, identification and control of robots[END_REF], given by (Fig. 1a) Consequently, if the time interval (t 0 , t f ) for the displacement between positions x S0 and x Sf is fixed, the "bang-bang" profile is the trajectory that minimizes the value of the maximal acceleration max S x . Thus, in order to minimize S x for a displacement during the fixed time interval (t 0 , t f ), the "bang-bang" profile has to be applied on the displacement x S on the manipulator total mass centre.

( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) S S
           f0 xx (4)

Observations about the modification of the optimal trajectory for taking into account the actuators properties

It should be mentioned that the given "bang-bang" profile (Fig. 1a) is based on theoretical considerations. In reality, the actuators are unable to achieve discontinuous efforts. Therefore, this motion profile should be modified by a trapezoidal profile (Fig. 1b) in order to take into account the actuators properties in terms of maximal admissible effort variations.

For a given time interval (t 0 , t f ), the trapezoidal profile, as we define it, is characterized by two parameters: t 1 , t 2 (Fig. 1b). In order to find the optimal values for t 1 and t 2 , the following problem should be considered:

under the constraints

  max max i i d dt    (6)
where  i is the input effort of the actuator i and   max i  the maximal admissible input effort variation for the actuator i. This problem is highly non linear, therefore it can be solved by numerical optimization methods. It should be mentioned that in the illustrative examples given in section 3, the trapezoidal profile taking into account the actuators properties has been found using the optimisation function "fgoalattain" of Matlab.

Expression of the manipulator coordinates as a function of the mass centre parameters

Once the displacement of the manipulator centre of masses is defined, the second problem is to find the articular (or Cartesian) coordinates corresponding to this displacement. For this purpose, let us consider a manipulator composed of n links. The mass of the link i is denoted as m i (i = 1, …, n) and the position of its centre of masses as x Si . Once the articular coordinates q or Cartesian variables x are known, the values of x Si may easily be obtained using the manipulator kinematics relationships. As a result, the position of the manipulator centre of masses, defined as

1 1 n S i Si i tot m m    xx , where 1 n tot i i mm    (7)
may be expressed as a function of x or q. But, in order to control the manipulator, the inverse problem should be solved, i.e. it is necessary to express variables q or x as a function of x S . Here, two cases should be distinguished:

(i) dim(x S ) = dim(q), i.e. the manipulator has got as many actuators as controlled variables for the displacements x S of the centre of masses (two variables for planar cases, three variables for spatial problems). In such case, the variables q or x can be directly expressed as a function of x S using (7), i.e. q = f(x S ).

(ii) dim(x S ) < dim(q), i.e. the manipulator has got more actuators than controlled variables. In such case, the problem is under-determined as there are more parameters in variables q or x than in x S . In order to solve it, let us consider that p 0 parameters of vector q 0 (or x 0 ) and p f parameters of vector q f (or x f ) are fixed. In a first task, it is necessary to define the m-p 0 and m-p f other parameters of the initial and final manipulator configurations (m = dim(q)). The way to fix it is to find the manipulator initial and final configurations, taking into account the p 0 initial and p f final fixed parameters, that will allow minimizing the norm of the vector SS  f0 xx , i.e. the length of the displacement of the manipulator centre of masses. Then, the second task is to choose m-k articular variables among the m possible of vector q (k = dim(x S )). These m-k variables, denoted as q m-k will be controlled using some classical displacement law given at [START_REF] Arakelian | Shaking force and shaking moment balancing of mechanisms: a historical review with new examples[END_REF] or can be used in order to minimize some other performance criteria, such as the shaking moments or some other interesting performance criterion (see section 3.2). The k other variables, denoted as q k , should be expressed as a function of x S and q m-k using (7), i.e. q k = f(x S , q m-k ).

In order to demonstrate the proposed balancing method, two illustrative examples are given in the following section.

Illustrative examples

The planar 2R serial manipulator

Let us consider the shaking force minimization of a 2R serial manipulator (Fig. 2). This manipulator is controlled using two rotary actuators having two input parameters which are denoted as q 1 and q 2 . For simulations the following parameters have been used:

-l OA = 0.5 m, l AB = 0.3 m, where l OA and l AB are the lengths of segments OA and AB, respectively; -r 1 = 0.289, where l OS1 = r 1 l OA and r 2 = 0.098, where l AS2 = r 2 l AB , l OS1 and l AS2 being the lengths of segments OS 1 and AS 2 , respectively. The mass and inertia parameters are:

-m 1 = 24.4 kg and m 2 = 8.3 kg, where m i is the mass of element i (i = 1, 2); -m tool = 5 kg, where m tool is the payload; -I 1 = 1.246 kg.m² and I 2 = 0.057 kg.m², where I i is the axial moment of inertia of element i.

Let us now express the articulated joint positions q = [q 1 , q 2 ] T as a function of the position x S of the manipulator centre of masses. From [START_REF] Arakelian | Complete shaking force and shaking moment balancing of linkages[END_REF], we obtain: x q q q q m r l m l r l y q q q q mm q q q m ll q q q m

                                               
x [START_REF] Herder | A counter-rotary counterweight for lightweight dynamic balancing[END_REF] This relationship leads to:

2 2 2 1 1 1 1 2 ( cos ) ( sin ) 0 S eq S eq eq x l q y l q l      (9)
where 1

1 1 2 ( ) / eq tool OA tot l m r m m l m    and 2 2 2 ( ) / eq tool AB tot l m r m l m  .
Replacing cos q 1 and sin q 1 by 22 11

(1 ) /(1 ) tt  and 2 [START_REF] Gosselin | Synthesis and design of reactionless treedegree-of-freedom parallel mechanisms[END_REF] 2 /(1 ) tt  (t 1 = tan(q 1 /2)), respectively, and developing (9), we obtain:

2 2 2 1 1 2 tan b b c a q ca            (10)
In [START_REF] Ricard | On the design of reactionless parallel manipulators[END_REF], the sign ± stands for the two possible working modes of the manipulator (for simulations, the working mode with the "+" sign is used). Once q 1 is known, q 2 may easily be found from (8): x l q

        . ( 12
)
Let us now test the proposed approach. In order to show the efficiency of this optimal planning, several trajectories are tested. These trajectories are defined as follows. First, the maximal inscribed square inside of the workspace is found (Fig. 3). For this manipulator, it is a square of length 0.55 m, of which centre E is located at x = 0 m and y = 0.475 m. Then, in order to avoid problems due to the proximity of singular configuration, the tested zone is restricted to a square centred in E of edge length equal to 0.45 m (in grey on Fig. 3). Finally, we discretize each edge into four segments delimited by the points P i (i = 1 to 16). The tested trajectories will be the segments P 1 P 13 , P 2 P 12 , P 3 P 11 , P 4 P 10 , P 5 P 9 , P 15 P 7 , P 14 P 8 and P 13 P 9 . Each trajectory will have duration of 0.5 s and, for each trajectory, three different kinds of motion profiles are applied:

1. a fifth order polynomial profile is applied on the displacement of the manipulator endeffector; 2. a "bang-bang" profile is applied on the displacement of the manipulator centre of masses; 3. a trapezoidal acceleration variation is applied on the displacement of the manipulator centre of masses, taking into account that, for each actuator, the input effort variation is limited by 3.10 4 Nm/s. The displacements of the end-effector and manipulator links centre of masses for the trajectory P 5 P 9 are shown in Fig. 4. These trajectory parameters are implemented into ADAMS software and it is computed the variations of shaking forces. Fig. 5 presents the shaking force transmitted by the manipulator for trajectory P 5 P 9 . The obtained results for the whole paths are summarized in table 1. It is shown that the optimal trajectory planning ("bang-bang profile") allows the reduction of the shaking force from 36% up to 76.7 %. Moreover, it appears that for given actuator parameters, the minimizations obtained in the cases of the "bang-bang" and trapezoidal profiles are very close (less than 1 %). It is due to the fact that the actuators can apply high input effort variations during a displacement. However, such a result depends on the actuator power capacity and it will be variable for each type of actuator.

Obviously, the rate of reduction depends on the design parameters of the robot. For each system, it will be different Let us now consider the second example.

The planar 3R serial manipulator

This manipulator is controlled using three rotary actuators (Fig. 6), with three input parameters which are denoted as q 1 , q 2 and q 3 . The link parameters are the following:

-l OA = 0.5 m, l AB = 0.3 m, l BC = 0.1 m, where l OA , l AB and l BC are the lengths of segments OA , AB and BC, respectively; -r 1 = 0.289, where l OS1 = r 1 l OA , r 2 = 0.098, where l AS2 = r 2 l AB , and r 3 = 0.5, where l BS3 = r 3 , l BC , l OS1 , l AS2 and l BS3 being the lengths of segments OS 1 , AS 2 and BS 3 , respectively. Its mass and inertia parameters are: -m 1 = 24.4 kg, m 2 = 8.3 kg and m 3 = 2 kg, where m i is the mass of element i (i = 1, 2,

3); -m tool = 5 kg, where m tool is the payload; -I 1 = 1.246 kg.m², I 2 = 0.057 kg.m² and I 3 = 0.025 kg.m², where I i is the axial moment of inertia of element i. In order to have the possibility to control the manipulator, let us express the relation between the articulated joint positions q = [q 1 , q 2 , q 3 ] T and the position x S of the manipulator centre of masses. From (7), we obtain: x q q q q m r l m l r l y q q q q mm q q q m l l r l q q q m m l m q q q     . In [START_REF] Papadopoulos | Design and motion planning for a zeroreaction manipulator[END_REF], there are three unknowns q 1 , q 2 , q 3 for two fixed parameters x S and y S . Therefore, as mentioned in section 2, a way to solve this problem is to consider that one parameter, for example , is used to minimize some objective function. Then the expressions of q 1 , q 2 and q 3 can be found as a function of x S , y S and . In the remainder of the paper, angle  is used in order to minimize the shaking moment m sh of the robot. Obviously, if necessary, it can be replaced by another criterion, such the energy, the torques, etc.

                                                        x 1 1 2 1 1 2 cos cos( ) cos sin sin( ) sin OA AB BC q q q ll q q q                        (13)
(13) leads to the following loop closure equation: (( cos ) cos ) (( sin ) sin ) 0 S eq eq S eq eq eq x l l q y l l q l

        (14)
where Replacing cos q 1 and sin q 1 by 22 11

1 1 1 2 3 ( ) / eq tool OA tot l m r m m m l m     ,
(1 ) /(1 ) tt  and 2 112 /(1 ) tt  (t 1 = tan(q 1 /2)), respectively, and developing [START_REF] He | Optimal motion planning of parallel redundant Mechanisms woth Shaking Force Reduction[END_REF], we obtain:

2 2 2 1 1 2 tan b b c a q ca            (15)
where [START_REF] Papadopoulos | Design and motion planning for a zeroreaction manipulator[END_REF]

( cos )

eq S eq a l x l

    , ( 16a
) 13 2 ( sin ) eq S eq b l y l     , ( 16b
) 2 2 2 2 3 3 1 2 ( cos ) ( sin )
S eq S eq eq eq c x l y l l l

       . ( 16c
)
In expression [START_REF] Nenchev | Reaction nullspace control of flexible structure mounted manipulator systems[END_REF], the sign ± stands for the two possible working modes of the manipulator (for simulations, the working mode with the "+" sign is used). Once q 1 is known, q 2 and q 3 may easily be found from ( 13): y l l q qq x l l q

          . ( 17
) 3 1 2 q q q     (18)
Let us now test the proposed approach with this manipulator. The tested trajectories are defined as follows. First, the maximal inscribed square inside of the workspace, for any endeffector orientation, is found (Fig. 7). For this manipulator, it is a square of length 0.375 m, of which centre E is located at x = 0 m and y = 0.487 m. Then, in order to avoid problems due to the proximity of singular configuration, the tested zone is restricted to a square centred in E of edge length equal to 0.3 m (in grey on Fig. 7). Finally, we discretize each edge into four segments delimited by the points P i (i = 1 to 16). The tested trajectories will be the segments P 1 P 13 , P 2 P 12 , P 3 P 11 , P 4 P 10 , P 5 P 9 , P 15 P 7 , P 14 P 8 and P 13 P 9 . It should be noted that in this case there is an independent parameter which can be defined from complementary condition describing the orientation of the end-effector. For numerical simulations, it is chosen to begin the tested trajectories with an end-effector orientation  0 = 0 deg and to finish it at  f = 120 deg. The simultaneous minimization of the shaking force and the shaking moment cannot be done without using an optimization algorithm in order to solve the following problem:

max( ) min sh m   (19)
under the constraints 00 ( ) , ( )

ff tt      (20a) 0 ( ) ( ) 0 f tt   (20b) 0 ( ) ( ) 0 f tt   (20c) 0 () SS t  0 xx , () S f S t  f xx (20d)
Several motion profiles for  can be tested. Here it is proposed to use polynomials. Our observations showed that the polynomial function that makes it possible to obtain optimal results is of degree 8. Each trajectory will have duration of 0.5 s and, for each trajectory three different kinds of motion profile are applied:

1. a fifth order polynomial profile is applied on the displacement (translation and rotation) of the manipulator end-effector; 2. a "bang-bang" profile is applied on the displacement of the manipulator centre of masses and the angle  is optimized in order to minimize the shaking moment; 3. a trapeze acceleration profile is applied on the displacement of the manipulator centre of masses, taking into account that, for each actuator, the input effort variation is limited by 3.10 4 Nm/s; the trajectory for angle  optimized in the previous case is used in order to compute the actuator displacements.

The displacements of the end-effector and manipulator links centre of masses for the trajectory P 15 P 7 are shown in Fig. 8. Fig. 9 presents the shaking force and shaking moment for the path P 15 P 7 . The obtained results for the whole paths are summarized in table 2 and 3. It is shown that the optimal trajectory planning ("bang-bang" profile) allows the reduction of the shaking forces from 48 % up to 62.2 %. Moreover, with a simultaneous optimal control of angle , the shaking moment can be reduced from 37.2 % up to 61 %.

As previously mentioned, these results depend on the design parameters of the used robot. For another manipulator, they will be different. But, in any case the shaking force and moment shall be decreased.

Table 2. Maximal value of the shaking force norm for the tested trajectories on the 3R serial manipulator.

Followed path

P

Observations about input torques

The main drawback of the shaking force balancing by counterweights is the increase of the inertia of moving links caused by adding masses, and consequently, the increase of input torques. The advantage of the suggested balancing method is in the fact that the shaking forces are only reduced by optimal control of moving links, without adding counterweights. It results in the fact that the input torques are considerably lower than in the case of balancing by counterweights. To illustrate this advantage for examined 2R serial manipulator, three kinds of simulations have been carried out using dynamic simulation software: (a) unbalanced manipulator carrying out a straight line trajectory along P 5 P 9 (Fig. 4) using a fifth order polynomial motion profile; (b) manipulator balanced by counterweights along the same trajectory1 ; (c) manipulator controlled via optimal centre of mass displacement between points P 5 and P 9 .

The obtained results are given in Figure 11. The software simulations showed that in comparison with mass balanced manipulator a 92% reduction in input torque is achieved, and, in comparison with unbalanced manipulator a 73% reduction in input torque is achieved.

Finally, we would like to mention that the method proposed in this paper focused exclusively on the force balancing because it is carried out by optimisation of the trajectory of the manipulator's center of mass. However, as was shown above, it also allows the reduction of the shaking moment and the input torques. Such a result has been observed for many simulated manipulators. But it is not possible to pretend in any way that this will be true for any manipulator.

Conclusions

In this paper, we have presented a new approach, based on an optimal trajectory planning, which allows the considerable reduction of the shaking force. This simple and effective balancing method is based on the optimal control of the acceleration of the manipulator centre of masses. For this purpose, the "bang-bang" profile has been used. The aim of the suggested method consists in the fact that the manipulator is controlled not by applying end-effector trajectories but by planning the displacements of the total mass centre of moving links. The trajectories of the total mass centre of moving links are defined as straight lines and are parameterized with "bang-bang" profile. Such a control approach allows the reduction of the maximum value of the centre of mass acceleration and, consequently, the reduction in the shaking force. It should be mentioned that such a solution is also very favourable for reduction of input torques because it is carried out without adding counterweights. The proposed balancing method has been illustrated via two examples. The numerical simulations showed that considerable reduction in shaking force and input torques were achieved.

Tables captions

Table 1. Maximal value of the shaking force norm for the tested trajectories of the 2R serial manipulator. Table 2. Maximal value of the shaking force norm for the tested trajectories on the 3R serial manipulator. Table 3. Maximal value of the shaking moment for the tested trajectories on the 3R serial manipulator.

Figure 1 .

 1 Figure 1. Motion profiles used for the shaking force minimization.

Figure 2 .

 2 Figure 2. Schematics of the 2R serial manipulator.

Figure 3 .Figure 4 .

 34 Figure 3. The tested trajectories of the 2R serial manipulator.

Figure 5 .

 5 Figure 5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line).

Figure 6 .

 6 Figure 6. Schematics of the 3R serial manipulator.

  where

Figure 7 .Figure 8 .

 78 Figure 7. The tested trajectories of the 3R serial manipulator.

Figure 9 .

 9 Figure 9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line).

Figure 10 .

 10 Figure 10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line).

(a) input torque 1 (b) input torque 2 Figure 11 .

 1211 Figure 11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Figures captionsFigure 1 .

 1 Figures captions

Figure 2 .Figure 3 .

 23 Figure 1. Motion profiles used for the shaking force minimization.Figure 2. Schematics of the 2R serial manipulator. Figure 3. The tested trajectories of the 2R serial manipulator.Figure 4. Manipulator end-effector displacements along the trajectory P 5 P 9 : (a) for case 1 and (b) for cases 2 and 3.Figure5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure6. Schematics of the 3R serial manipulator. Figure7. The tested trajectories of the 3R serial manipulator. Figure8. Manipulator end-effector displacements along the trajectory P 15 P 7 : (a) for case 1 and (b) for optimal cases 2 and 3. Figure9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Figure 4 .

 4 Figure 1. Motion profiles used for the shaking force minimization.Figure 2. Schematics of the 2R serial manipulator. Figure 3. The tested trajectories of the 2R serial manipulator.Figure 4. Manipulator end-effector displacements along the trajectory P 5 P 9 : (a) for case 1 and (b) for cases 2 and 3.Figure5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure6. Schematics of the 3R serial manipulator. Figure7. The tested trajectories of the 3R serial manipulator. Figure8. Manipulator end-effector displacements along the trajectory P 15 P 7 : (a) for case 1 and (b) for optimal cases 2 and 3. Figure9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Figure 5 .

 5 Figure 1. Motion profiles used for the shaking force minimization.Figure 2. Schematics of the 2R serial manipulator. Figure 3. The tested trajectories of the 2R serial manipulator.Figure 4. Manipulator end-effector displacements along the trajectory P 5 P 9 : (a) for case 1 and (b) for cases 2 and 3.Figure5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure6. Schematics of the 3R serial manipulator. Figure7. The tested trajectories of the 3R serial manipulator. Figure8. Manipulator end-effector displacements along the trajectory P 15 P 7 : (a) for case 1 and (b) for optimal cases 2 and 3. Figure9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Figure 6 .

 6 Figure 1. Motion profiles used for the shaking force minimization.Figure 2. Schematics of the 2R serial manipulator. Figure 3. The tested trajectories of the 2R serial manipulator.Figure 4. Manipulator end-effector displacements along the trajectory P 5 P 9 : (a) for case 1 and (b) for cases 2 and 3.Figure5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure6. Schematics of the 3R serial manipulator. Figure7. The tested trajectories of the 3R serial manipulator. Figure8. Manipulator end-effector displacements along the trajectory P 15 P 7 : (a) for case 1 and (b) for optimal cases 2 and 3. Figure9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Figure 7 .

 7 Figure 1. Motion profiles used for the shaking force minimization.Figure 2. Schematics of the 2R serial manipulator. Figure 3. The tested trajectories of the 2R serial manipulator.Figure 4. Manipulator end-effector displacements along the trajectory P 5 P 9 : (a) for case 1 and (b) for cases 2 and 3.Figure5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure6. Schematics of the 3R serial manipulator. Figure7. The tested trajectories of the 3R serial manipulator. Figure8. Manipulator end-effector displacements along the trajectory P 15 P 7 : (a) for case 1 and (b) for optimal cases 2 and 3. Figure9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Figure 8 .Figure 9 .

 89 Figure 1. Motion profiles used for the shaking force minimization.Figure 2. Schematics of the 2R serial manipulator. Figure 3. The tested trajectories of the 2R serial manipulator.Figure 4. Manipulator end-effector displacements along the trajectory P 5 P 9 : (a) for case 1 and (b) for cases 2 and 3.Figure5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure6. Schematics of the 3R serial manipulator. Figure7. The tested trajectories of the 3R serial manipulator. Figure8. Manipulator end-effector displacements along the trajectory P 15 P 7 : (a) for case 1 and (b) for optimal cases 2 and 3. Figure9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Figure 10 .

 10 Figure 1. Motion profiles used for the shaking force minimization.Figure 2. Schematics of the 2R serial manipulator. Figure 3. The tested trajectories of the 2R serial manipulator.Figure 4. Manipulator end-effector displacements along the trajectory P 5 P 9 : (a) for case 1 and (b) for cases 2 and 3.Figure5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure6. Schematics of the 3R serial manipulator. Figure7. The tested trajectories of the 3R serial manipulator. Figure8. Manipulator end-effector displacements along the trajectory P 15 P 7 : (a) for case 1 and (b) for optimal cases 2 and 3. Figure9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Figure 11 .

 11 Figure 1. Motion profiles used for the shaking force minimization.Figure 2. Schematics of the 2R serial manipulator. Figure 3. The tested trajectories of the 2R serial manipulator.Figure 4. Manipulator end-effector displacements along the trajectory P 5 P 9 : (a) for case 1 and (b) for cases 2 and 3.Figure5. Variations of the shaking forces in the case of the trajectory P 5 P 9 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure6. Schematics of the 3R serial manipulator. Figure7. The tested trajectories of the 3R serial manipulator. Figure8. Manipulator end-effector displacements along the trajectory P 15 P 7 : (a) for case 1 and (b) for optimal cases 2 and 3. Figure9. Variations of the shaking forces in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure10. Variations of the shaking moment in the case of the trajectory P 15 P 7 : case 1 (black full line), case 2 (black dashed line) and case 3 (grey full line). Figure11. Manipulator input torques for trajectory P 5 P 9 corresponding to the three simulated models: (i) unbalanced manipulator carrying out a straight line trajectory of the end-effector using a fifth order polynomial motion profile (black full line); (ii) manipulator balanced by counterweights along the same trajectory (grey full line); (iii) manipulator controlled via optimal centre of mass displacement (black dashed line).

Table 1 .

 1 Maximal value of the shaking force norm for the tested trajectories of the 2R serial manipulator. Followed path P 1 P 13 P 2 P 12 P 3 P 11 P 4 P 10 P 5 P 9 P 15 P 7 P 14 P 8 P 13 P 9 max(|| ||)

	f (N)	sh	Case 1 Case 2 Case 3	194.7 165.3 178.8 178.2 155.3 218.7 201.3 195.4 121.0 85.0 47.5 41.5 44.3 136.3 121.7 111.8 124.0 86.0 48.0 41.8 44.7 137.7 123.2 113.1
	% of	Cases 2/1	37.8	48.6	73.4	76.7	71.5	37.8	39.5	42.8
	reduction	Cases 3/1	36.3	48.0	73.2	76.5	71.2	37.0	38.9	42.1

Table 3 .

 3 1 P 13 P 2 P 12 P 3 P 11 P 4 P 10 P 5 P 9 P 15 P 7 P 14 P 8 P 13 P 9 max(|| ||) Maximal value of the shaking moment for the tested trajectories on the 3R serial manipulator. Followed path P 1 P 13 P 2 P 12 P 3 P 11 P 4 P 10 P 5 P 9 P 15 P 7 P 14 P 8 P 13 P 9

	f (N)	sh	Case 1 Case 2 Case 3	158.7 144.9 133.6 119.6 103.1 321.3 266.3 238.7 81.8 77.2 68.8 58.1 47.9 121.1 105.5 92.2 82.3 77.6 69.0 58.3 48.0 122.3 106.5 93.0
	% of	Cases 2/1	48.4	46.7	48.5	51.4	53.5	62.3	60.3	61.4
	reduction	Cases 3/1	48.1	46.4	48.3	51.2	53.4	61.9	60.0	61.0

The counterweights are located at 0.2 m and 0.35 m from the joints centres A and O, respectively.


Chapitre d'équation 1 Section 1

Abstract-Off-line robot dynamic identification methods are mostly based on the use of the Inverse Dynamic Identification Model (IDIM), which calculates the joint force/torque that is linear in relation to the dynamic parameters, and on the use of linear least squares technique to calculate the parameters (IDIM-LS technique). The joint forces/torques are calculated as the product of the known control signal (the current reference) by the joint drive gains. Then it is essential to get accurate values of joint drive gains to get accurate identification of inertial parameters. In this paper it is proposed a new method for the identification of the total joint drive gains in one step using available joint sampled data given by the standard controller of the moving robot and using only the weighted mass of a payload, without any CAD values of its inertial parameters. A new inverse dynamic model calculates the current reference signal of each joint j that is linear in relation to the dynamic parameters of the robot, to the inertial parameters of a known mass fixed to the end-effector, and to the inverse of the joint j drive gain. This model is calculated with current reference and position sampled data while the robot is tracking one reference trajectory without load on the robot and one trajectory with the known mass fixed on the robot. Each joint j drive gain is calculated independently by the weighted LS solution of an over-determined linear systems obtained with the equations of the joint j. The method is experimentally validated on an industrial Stäubli RX-90 robot.

I. INTRODUCTION

EVERAL schemes have been proposed in the literature to identify the dynamic parameters of robots [START_REF] De | Parameters identification of robots manipulators via sequential hybrid estimation algorithms[END_REF]- [START_REF] Lu | Experimental determination of dynamic parameters of robotic arms[END_REF]. Most of the dynamic identification methods have the following common features: -the use of an Inverse Dynamic Identification Model (IDIM) which calculates the joint force/torque linear in relation to the dynamic parameters, -the construction of an over-determined linear system of equations obtained by sampling IDIM while the robot is tracking some trajectories in closed-loop control, -the estimation of the parameter values using least squares techniques (LS). This procedure is called the IDIM-LS technique. The experimental works have been carried out either on prototypes in laboratories or on industrial robots and have shown the benefits in terms of accuracy in many cases. Good results can be obtained provided two main conditions are satisfied: -a well-tuned derivative band-pass filtering of joint position is used to calculate the joint velocities and accelerations, -the accurate values of joint drive gains g  are known to calculate the joint force/torque as the product of the known control signal calculated by the numerical controller of the robot (the current references) by the joint drive gains [START_REF] Restrepo | Calibration of drive chain of robot joints[END_REF]. This needs to calibrate the drive train constituted by a current controlled voltage source amplifier with gain i G which supplies a permanent magnet DC or a brushless motor with torque constant t K coupled to the link through direct or gear train with gear ratio N .

Because of large values of the gear ratio for industrial robots, ( N >50), joint drive gain,

it g NG K  
, is very sensitive to errors in i G and t K which must be accurately measured from special, time consuming , heavy tests, on the drive chain [START_REF] Restrepo | Calibration of drive chain of robot joints[END_REF].

In this paper it is proposed a new method for the identification of the total joint drive gains in one step, using current reference and position sampled data while the robot is tracking one reference trajectory without load fixed on the robot and one trajectory with a known mass fixed on the robot whose inertial parameters are measured l. Each joint j drive gain is calculated independently by the weighted LS solution of an over-determined linear system obtained with the equations of the joint j .

The method is experimentally validated on an industrial Stäubli RX-90 robot.

The paper is organized as follows: section 2 recalls the dynamic modelling and identification procedures. Section 3 deals with the new modelling and identification method for the robot drive gains parameters. Section 4 presents the experimental validations. Finally, section 5 gives the conclusion.

II. USUAL INVERSE DYNAMIC MODELS AND IDENTIFICATION

A. Dynamic Modelling

In the following, all mechanical variables are given in SI units on the joint side. All forces/torques, positions, velocities and accelerations have a conventional positive sign

New Method for Global Identification of the Joint Drive Gains of Robots using a Known Payload Mass

Maxime Gautier and Sébastien Briot S in the same direction. That defines a motor convention for the mechanical behaviour. The inverse dynamic model (IDM) of a rigid robot composed of n moving links calculates the motor force/torque vector, as a function of the generalized coordinates and their derivatives. It can be obtained from the Newton-Euler or the Lagrangian equations [START_REF] Khalil | Modeling, identification and control of robots[END_REF]:

dyn in f     (1)
where: τ dyn is the ( 1) n vector of dynamic forces/torques due to the inertial, centrifugal, Coriolis, and gravitational effects:

( ) ( ) ( ) dyn = M q q +C q,q q +Q q  (2) q, q, q, are respectively the ( 1) n vectors of generalized joint positions, velocities and accelerations, ( M q ) is the () nn  robot inertia matrix, C( q,q ) is the () nn  matrix of centrifugal and Coriolis effects, Q( q ) is the ( 1) n vector of gravitational forces/torques. τ in is the ( 1) n input electromagnetic force/torque vector of the drive chain, in SI units on the joint side.

()

in idm 0 = g v v      (3)
where:

idm v  is the ( 1) n vector of current references of the current amplifiers, v τ0 is the ( 1)
n vector of amplifier offsets, g τ is the () nn  matrix of the drive gains,

it g NG K   (4)
N is the () nn  gear ratios matrix of the joint drive chains such as, m q Nq  , with m q the ( 1) K are data that are generally given by the manufacturers, but with a quite important inaccuracy. They can also be independently identified using the methods presented in [START_REF] Restrepo | Calibration of drive chain of robot joints[END_REF]. However, these methods involve special tests on the joint drive train of the robot which are heavy and time consuming procedures.

n
f  is the ( 1)
n vector of the loss force/torque due to frictions and motor iron losses, eddy currents and hysteresis effect. Usually, it is approximated with a viscous friction coefficient and a Coulomb friction force/torque: sign

f v c coff F q F ( q ) F      (5)
where:

v F is the () nn  diagonal matrix of viscous parameters, c F is the () nn  diagonal matrix of dry friction parameters, and sign(.) denotes the sign function, coff F is a ( 1) n vector of asymmetrical Coulomb friction force/torque between positive and negative velocities.

Thus (1) becomes: sign( )

dyn idm v c off F q F q        (6)
where:

off coff 0 F g v    is the ( 1) n vector of offset
force/torque that regroups the amplifier offset and the asymmetrical Coulomb friction coefficient.

idm idm gv   
defines the motor force/torque which includes the amplifier offset component.

Then the inverse dynamic model (IDM) is given by: = ( ) + ( , ) idm τ M q q N q q (7) ( , ) N q q regroups the ( 1)

n vector of centrifugal,
Coriolis, gravitational, friction and offset forces/torques:

( )= ( ) ( ) sign( ) v c off N q, q C q,q q Q q F q F q      (8)

B. Identification Model

The modified Denavit and Hartenberg notation allows obtaining a dynamic model [START_REF] Khosla | Parameter identification of robot dynamics[END_REF] that is linear in relation to a set of standard dynamic parameters,

st χ : ( ) ( ) idm st st st q,q,q, q,q,q      (9)
where: given by 12 ...

st ( q,q,q )  is the   st nn  jacobian matrix of
T T T n T st st st st        : j st
 is the (14×1) vector of standard dynamic parameters of joint and link j :

j j T st j j j j j j j j j j j j j off XX XY XZ YY YZ ZZ MX MY MZ M Ia Fv Fc        (10)
where:

j j j j j j
XX , XY , XZ , YY , YZ , ZZ are the 6 components of the inertia matrix of link j at the origin of frame j . The columns (: )

st ,i  of the matrix () st q,q,q  can be obtained using the recursive algorithm of Newton-Euler, which calculates idm  (7) in terms of the standard dynamic parameters, such that:

(: ) ( with 1 0 for )

ik st idm st st ,i q,q,q, , , k i         (10)
To increase the efficiency of this algorithm, we use the customized symbolic technique [START_REF] Khalil | Modeling, identification and control of robots[END_REF], [START_REF] Khalil | SYMORO+: A system for the symbolic modelling of robots[END_REF].

The identifiable parameters are the base parameters which are the minimum number of dynamic parameters from which the dynamic model can be calculated. They are obtained from the standard inertial parameters by regrouping some of them by means of linear relations [START_REF] Mayeda | Base parameters of manipulator dynamic models[END_REF], which can be determined for the serial robots using simple closed-form rules [3], [START_REF] Khalil | Modeling, identification and control of robots[END_REF], or by numerical method based on the QR decomposition [START_REF] Gautier | Numerical calculation of the base inertial parameters[END_REF].

The minimal dynamic model can be written using the b n base dynamic parameters  as follows:

()

idm q,q,q     ( 11
)
where  is obtained from st  by eliminating the columns corresponding to the non identifiable parameters.

Because of perturbations due to noise measurement and modelling errors, the actual force/torque  differs from idm τ by an error, e , such that:

()

idm τ e q,q,q e        ( 12
)
where:

gv    ( 13
)
v  is the ( 1) n vector of the actual current references of the current amplifiers. Equation ( 12) represents the Inverse Dynamic Identification Model (IDIM).

C. Identification of the Dynamic Parameters

The off-line identification of the base dynamic parameters  is considered, given measured or estimated off-line data for τ and   q, q, q , collected while the robot is tracking some planned trajectories.

(

) q, q, q in (12) are estimated with ( ) ˆq, q, q , respectively, obtained by band-pass filtering the measure of q [START_REF] Gautier | Dynamic identification of robots with power model[END_REF]. The principle is to sample the identification model [START_REF] Featherstone | Dynamics[END_REF] at a frequency m f in order to get an over-determined linear system of ()

m nn  equations and b n unknowns such that:     fm fm fm ˆŶ τ W q,q,q χ ρ  (14)
In order to cancel the high frequency torque ripple in   This parallel filtering procedure can be carried out with the Matlab decimate function [START_REF] Gautier | Dynamic identification of robots with power model[END_REF]. It is obtained:

    ˆŶ τ W q,q,q χ ρ  ( 15
)
ρ is the ( x1) r vector of errors, with md r n n / n  ,   ˆŴ q, q, q is the ( x ) b rn observation matrix. In Y and W , the equations of each joint j are grouped together such that: Using the base parameters and tracking "exciting" reference trajectories, a well conditioned matrix W is obtained. The LS solution χ of ( 15) is given by:

11 ( ) ( ) ( ) ( ) TT T n T T n T Y Y ... Y ,W W ... W          (16)
    1 TT χ W W W Y W Y    (17)
Standard deviations i   , are estimated assuming that W is a deterministic matrix and  , is a zero-mean additive independent Gaussian noise, with a covariance matrix C  , such that: T2 ()

r CE ρρ I     (18)
E is the expectation operator and I r , the () rr  identity matrix. An unbiased estimation of the standard deviation

  is: 2 2 () Ŷ -W r b    (19)
The covariance matrix of the estimation error is given by:

T 2 T 1 [( )( ) ] ( ) ˆˆĈE χ χ χ χ W W         . () i 2 ˆˆĈ i,i     is the i th diagonal coefficient of ˆĈ  (20)
The relative standard deviation ri %   is given by:

100 ri i ˆî %      , for i  ≠ 0 (21)
The OLS can be improved by taking into account different standard deviations on joint j equations errors [START_REF] Gautier | Dynamic identification of robots with power model[END_REF]. Each equation of joint j in ( 15) is weighted with the inverse of the standard deviation of the error calculated from OLS solution of the equations of joint j , given by:

  j j j j Y τ W χ ρ  (22)
This weighting operation normalises the errors in [START_REF] Mayeda | Base parameters of manipulator dynamic models[END_REF] and gives the weighted LS (WLS) estimation of the parameters.

III. GLOBAL IDENTIFICATION OF THE JOINT DRIVE GAINS

A. Inverse Dynamic identification Model of the robot and the payload

The payload is considered as a link 1 n  fixed to the link n of the robot. The model [START_REF] Khalil | Symbolic calculation of the base inertial parameters of closed-loop robots[END_REF] becomes:

idm uL kL uL kL                (23)
where:

kL  is the ( 1)
L n  vector of the known inertial parameters of the payload; The model ( 12) becomes:

uL kL uL kL e                (24)
Taking into account that the joint j force/torque depends only on parameters j ,n  of links j to n , the IDIM joint j model is given by:

Equation ( 26) is the minimal model IDIM used to identify the drive gain of joint j .

B. Identification of the drive gains

Considering now that the robot has carried out two trajectories: (a) without the payload and (b) with the payload fixed to the end-effector, the sampling and filtering of the model IDIM (26) can be written as: where: j V is the vector of j v  samples, j a W is the observation matrix of joint j in the unloaded case, j b W is the observation matrix of joint j in the loaded case, j uL W is the observation matrix for joint j corresponding to the unknown payload inertial parameters, j kL W is the observation matrix for joint j corresponding to the known payload inertial parameters, *  is the weighted LS solution of the over-determined system (27). The last coefficient of vector *  gives the estimation of the inverse of the joint j drive gain. Using this global approach, it is possible to only use the knowledge we can have on the payload mass value L M in order to identify the robot drive gains. This mass can generally be accurately measured using a weighting machine. It is assumed that parameter L M appears explicitly in the expression of several joint torques (at least one; if not, the approach cannot be achieved), i.e. 0 j kL L WM  for all these joints.

For all of them, the drive gains can be identified using the previous approach. Then, using global identification models for the same joints, other exciting payload inertial parameters (such as the static moments

L MX , L MY or L MZ
) can be estimated and will be used in order to identify the drive gains for the other joints.

The advantage of this method is its simplicity. But its main drawback is that the identification is carried out sequentially and the estimation errors will be cumulated. As a result, it could be preferable to have a complete knowledge of all payload inertial parameters (e.g. using data obtained from CAD software) in order to decouple the drive gain identification for each joint. But, as it will be shown in the next section, the experimental results confirm the efficiency of the proposed approach.

IV. CASE STUDY

A. Description of the RX 90 kinematics

The Stäubli RX-90 robot (Fig. 1) has a serial structure with six rotational joints. Its kinematics is defined using the modified Denavit and Hartenberg notation (MDH) [START_REF] Khalil | A new geometric notation for open and closed loop robots[END_REF]. In this notation, the link j fixed frame is defined such that the j z axis is taken along joint j axis and the j x axis is along the common normal between j z and j1 z  (Fig. 1). The geometric parameters defining the robot frames are given in Table 1. The payload is denoted as the link 7. The parameter 0 j   , means that joint j is rotational,



.The main advantage of using the MDH notation is that the identifiable inertial parameters can be determined symbolically using simple closed-form rules [3], [START_REF] Khalil | Modeling, identification and control of robots[END_REF], [START_REF] Khalil | Symbolic calculation of the base inertial parameters of closed-loop robots[END_REF] and it is possible to generalize the proposed methods for tree structure robots and closed loop robots [START_REF] Hollerbach | Model Identification[END_REF], [START_REF] Khalil | Modeling, identification and control of robots[END_REF], [START_REF] Featherstone | Dynamics[END_REF], [START_REF] Gautier | Numerical calculation of the base inertial parameters[END_REF].

B. Identification of the drive gains

The method is validated using a calibrated payload (Fig. 2). Its mass has been measured with a weighting machine ( L M = 7.025 Kg± 0.050 Kg). The other parameters have been estimated using CAD software. They are given in table 2. Their values are accurate due to the simplicity of the payload shape (Fig. 2).

For the RX 90, it can be shown that: -the torques of joints 1 to 3 depends on all payload parameters; -the torques of joints 4 and 5 do not depend on parameter ZZ are too small to allow the correct calibration of the drive gains of joint 5 and 6.

Therefore, the drive gains of joints 1 to 4 only are identified. In order to identify the other drive gains, more exciting payload should be designed, especially with large values for parameters L MX and L MY . Staübli's control software provides directly the value of the joint torque. Therefore, it is proceeded to a calibration of the drive gains compared with a priori values used in the controller which are taken equal to 1. A first identification of both the RX 90 dynamic parameters and the load inertial parameters is achieved using the IDIM-LS method with the model (25), as detailed in [START_REF] Khalil | Identification of the payload inertial parameters of industrial manipulators[END_REF] (case 1). The results are presented in table 2. It can be observed that the mass L M and the first moment L MZ are poorly identified with the manufacturer's gains.

Then, the new drive gains for joints 1 to 4 are identified using our new approach (case 2). The identified drive gains are presented in table 3. It can be shown that there is an average of 8% error with respect to the a drive gains.

A new identification of the payload dynamic parameters is achieved with the identified gains. The results are presented in table 2. It can be observed that the mass L M and the first moment L MZ are better identified using the new gains. In table 4, the quality of identification for the two studied cases in terms of reduction of the norm of the error vector ρ is shown. With the identified gains, the quality of identification is largely improved.

Figure 3 shows that the reconstructed torques with the identified gains for joints 1 to 4 are close to the actual ones.

In order to show the sensitivity of the identification results to the variation of the drive gain values, a 20% error on these gains (with respect to those identified in case 2) is simulated (case 3). The quality of identification is clearly worst in this case (table 4). The identified payload parameters are given in table 2. Their values are far from reality. V. CONCLUSION This paper has presented a new method for the identification of the total drive gains for robot joints in one step. This is a method easy to implement which does not need any special test or measurement on elements inside the joint drive train. It is based on a IDIM-LS technique using current reference and position sampled data while the robot is tracking one reference trajectory without load fixed on the robot and one trajectory with a known payload fixed on the robot, whose inertial parameters are measured or calculated by a CAD model. The method has been experimentally validated on a Stäubli RX-90 robot. With the new identified drive gains, the identification of the total dynamic model of the robot has been improved and another payload has been accurately identified. This shows the effectiveness of the method.

  Manuscript received March 15, 2011. M. Gautier is with the Institut de Recherche en Communication et Cybernétique de Nantes (IRCCyN) and with the University of Nantes, 44321 Nantes France (phone: +33(0)240376960; fax: +33(0)240376930; email: Maxime.Gautier@irccyn.ec-nantes.fr).S. Briot is with the Institut de Recherche en Communication et Cybernétique de Nantes (IRCCyN), 44321 Nantes France (e-mail: Sebastien.Briot@irccyn.ec-nantes.fr).

  , MZ are the 3 components of the first moment of link j . j M is the mass of link j , j Ia is a total inertia moment for rotor and gears of actuator j .

fm Y τ

  and to window the identification frequency range into the model dynamics, a parallel decimation procedure low-pass filters in parallel fm Y and each column of fm W and resamples them at a lower rate, keeping one sample over d n .

  equations of joint j . The ordinary LS (OLS) solution χ minimizes ρ  .

uL

  is the ((10) 1) L n  vector of the unknown inertial parameters of the payload, kL  is the()

  are the respective non zero elements of line j in (24), j v  is the actual joint j current reference given by the numerical controller of the robot.Equation (25) can be rewritten as:

  means that the link 7 is fixed on the link 6 .

 6 Since all the joints are rotational then j  is the position variable of joint j . It should be mentioned that the variables of joints 2 and 3 in the present notation are obtained from the measured values of joints 2 and 3 using the relations, 2

Fig. 2 . 3 x 4 , x 5 4 Fig. 1 .

 234541 Fig. 2. The 7.025 Kg payload

TABLE 4

 4

			QUALITY OF IDENTIFICATION.
		Error norm 	Relative Error norm Ŷ 	ˆ 
	Case 1	77.28	0.075643	1.16293
	Case 2	76.11	0.0744996	1.14535
	Case 3	85.85	0.0840385	1.29288
	ˆŶW  	is the minimal norm of error, ˆ  is given by (19).

TABLE 3 IDENTIFIED

 3 DRIVE GAINS.

	j g 		0.9401 0.9619 1.1148 0.9146
	2   8.17e-2 3.57e-2 3.72e-2 5.03e-2 i 
		ri % 	4.3%	1.9%	1.7%	2.8%

 Joint 1 Joint 2 Joint 3 Joint 4

TABLE 2 IDENTIFICATION

 2 OF THE PAYLOAD DYNAMIC PARAMETERS. Parameter A priori value Estimated value 2

	i   		ri % 


Off-line robot dynamic identification methods are mostly based on the use of the Inverse Dynamic Identification Model (IDIM), which calculates the joint force/torque that is linear in relation to the dynamic parameters, and on the use of linear least squares technique to calculate the parameters (IDIM-LS technique). The joint forces/torques are calculated as the product of the known control signal (the current reference) by the joint drive gains. Then it is essential to get accurate values of joint drive gains to get accurate identification of inertial parameters. In this paper it is proposed a new method for the identification of the total joint drive gains in one step, using available joint sampled data given by the standard controller of the moving robot and using CAD or measured values of the inertial parameters of a known payload. A new inverse dynamic model calculates the current reference signal of each joint j that is linear in relation to the dynamic parameters of the robot, to the inertial parameters of a known payload fixed to the end-effector, and to the inverse of the joint j drive gain. This model is calculated with current reference and position sampled data while the robot is tracking one reference trajectory without load on the robot and one trajectory with the known payload fixed on the robot. Each joint j drive gain is calculated independently by the weighted LS solution of an over-determined linear systems obtained with the equations of the joint j. The method is experimentally validated on an industrial Stäubli RX-90 robot.

I. INTRODUCTION

EVERAL schemes have been proposed in the literature to identify the dynamic parameters of robots [START_REF] De | Parameters identification of robots manipulators via sequential hybrid estimation algorithms[END_REF]- [START_REF] Lu | Experimental determination of dynamic parameters of robotic arms[END_REF]. Most of the dynamic identification methods have the following common features: -the use of an Inverse Dynamic Identification Model (IDIM) which calculates the joint force/torque linear in relation to the dynamic parameters, -the construction of an over-determined linear system of equations obtained by sampling IDIM while the robot is tracking some trajectories in closed-loop control, -the estimation of the parameter values using least squares techniques (LS). This procedure is called the IDIM-LS technique. The experimental works have been carried out either on Manuscript received March 15, 2011. M. Gautier is with the University of Nantes and the Institut de Recherche en Communications et Cybernétique de Nantes (IRCCyN) 44321 Nantes France (phone: +33240376960; fax: +33240376930; e-mail: Maxime.Gautier@irccyn.ec-nantes.fr).

S. Briot is with the CNRS and IRCCyN, 44321 Nantes France (e-mail: Sebastien.Briot@irccyn.ec-nantes.fr).

This work is supported by the French grant ANR-ARPEGE 2010-SEGI-003-COROUSSO.

prototypes in laboratories or on industrial robots and have shown the benefits in terms of accuracy in many cases. Good results can be obtained provided two main conditions are satisfied: -a well-tuned derivative band-pass filtering of joint position is used to calculate the joint velocities and accelerations, -the accurate values of joint drive gains g  are known to calculate the joint force/torque as the product of the known control signal calculated by the numerical controller of the robot (the current references) by the joint drive gains [START_REF] Restrepo | Calibration of drive chain of robot joints[END_REF]. This needs to calibrate the drive train constituted by a current controlled voltage source amplifier with gain i G which supplies a permanent magnet DC or a brushless motor with torque constant t K coupled to the link through direct or gear train with gear ratio N .

Because of large values of the gear ratio for industrial robots, ( N >50), joint drive gain,

i t g NG K  
, is very sensitive to errors in i G and t K which must be accurately measured from special, time consuming , heavy tests, on the drive chain [START_REF] Restrepo | Calibration of drive chain of robot joints[END_REF], [START_REF] Corke | In situ measurement of robot motor electrical constants[END_REF].

In this paper it is proposed a new method for the identification of the total joint drive gains in one step, using current reference and position sampled data while the robot is tracking one reference trajectory without load fixed on the robot and one trajectory with a known payload fixed on the robot, whose inertial parameters are measured or calculated by a CAD model. Each joint j drive gain is calculated independently by the weighted LS solution of an overdetermined linear system obtained with the equations of the joint j .

The method is experimentally validated on an industrial Stäubli RX-90 robot.

The paper is organized as follows: section 2 recalls the dynamic modelling and identification procedures. Section 3 deals with the new modelling and identification method for the robot drive gains parameters. Section 4 presents the experimental validations. Finally, section 5 gives the conclusion.

II. USUAL INVERSE DYNAMIC MODELS AND IDENTIFICATION

A. Dynamic Modelling

In the following, all mechanical variables are given in SI

New Method for Global Identification of the Joint Drive Gains of Robots using a Known Inertial Payload

Maxime Gautier and Sébastien Briot S units on the joint side. All forces/torques, positions, velocities and accelerations have a conventional positive sign in the same direction. That defines a motor convention for the mechanical behaviour. The inverse dynamic model (IDM) of a rigid robot composed of n moving links calculates the motor force/torque vector, as a function of the generalized coordinates and their derivatives. It can be obtained from the Newton-Euler or the Lagrangian equations [START_REF] Khalil | Modeling, identification and control of robots[END_REF]:

dyn in f      (1)
where: τ dyn is the ( 1) n  vector of dynamic forces/torques due to the inertial, centrifugal, Coriolis, and gravitational effects:

( ) ( ) ( ) dyn = M q q+C q,q q+Q q

    (2)
q, q, q,   are respectively the ( 1) n  vectors of generalized joint positions, velocities and accelerations, ( M q ) is the ( ) n n  robot inertia matrix,

C( q,q )
 is the ( ) n n  matrix of centrifugal and Coriolis effects, Q( q ) is the ( 1) n  vector of gravitational forces/torques.

τ in is the ( 1) n  input electromagnetic force/torque vector of the drive chain, in SI units on the joint side.

( ) It should be mentioned that the components of matrices N , i G and t K are data that are generally given by the manufacturers, but with a quite important inaccuracy. They can also be independently identified using the methods presented in [START_REF] Restrepo | Calibration of drive chain of robot joints[END_REF]. However, these methods involve special tests on the joint drive train of the robot which are heavy and time consuming procedures. f  is the ( 1) n  vector of the loss force/torque due to frictions and motor iron losses, eddy currents and hysteresis effect. Usually, it is approximated with a viscous friction coefficient and a Coulomb friction force/torque: sign

in idm 0 = g v v      (3) where
f v c c o f f F q F (q) F        (5)
where

: v F is the ( ) n n  diagonal matrix of viscous parameters, c F is the ( ) n n
 diagonal matrix of dry friction parameters, and sign(.) denotes the sign function, coff F is a [START_REF] De | Parameters identification of robots manipulators via sequential hybrid estimation algorithms[END_REF] n  vector of asymmetrical Coulomb friction force/torque between positive and negative velocities. Thus (1) becomes: sign( )

dyn idm v c off F q F q          (6)
where:

off coff 0 F g v      is the ( 1)
n  vector of offset force/torque that regroups the amplifier offset and the asymmetrical Coulomb friction coefficient.

idm idm g v    
defines the motor force/torque which includes the amplifier offset component.

Then the inverse dynamic model (IDM) is given by: = ( ) + ( , )

idm τ M q q N q q   (7)
( , ) N q q  regroups the ( 1) n  vector of centrifugal, Coriolis, gravitational, friction and offset forces/torques:

( )= ( ) ( ) sign( ) v c o f f N q, q C q,q q Q q F q F q           (8)

B. Identification Model

The modified Denavit and Hartenberg notation allows obtaining a dynamic model ( 7) that is linear in relation to a set of standard dynamic parameters, st χ :

( ) ( ) idm st st st q,q,q, q,q,q          (9)
where:

st ( q,q,q )    is the   ...

T T T n T st st st st          : j st
 is the (14×1) vector of standard dynamic parameters of joint and link j :

j j T s t j j j j j j j j j j j j j o f f XX XY XZ YY YZ ZZ MX MY MZ M Ia Fv Fc        (10)
where:

j j j j j j
XX , XY , XZ , YY , YZ , ZZ are the 6 components of the inertia matrix of link j at the origin of frame j . The columns (: ) st ,i  of the matrix ( ) st q,q,q    can be obtained using the recursive algorithm of Newton-Euler, which calculates idm  (7) in terms of the standard dynamic parameters, such that:

(: ) ( with 1 0 for ) i k st idm st st ,i q,q,q, , , k i           (10)
To increase the efficiency of this algorithm, we use the customized symbolic technique [START_REF] Khalil | Modeling, identification and control of robots[END_REF], [START_REF] Khalil | SYMORO+: A system for the symbolic modelling of robots[END_REF]. The identifiable parameters are the base parameters which are the minimum number of dynamic parameters from which the dynamic model can be calculated. They are obtained from the standard inertial parameters by regrouping some of them by means of linear relations [START_REF] Mayeda | Base parameters of manipulator dynamic models[END_REF], which can be determined for the serial robots using simple closed-form rules [START_REF] Gautier | Direct calculation of minimum set of inertial parameters of serial robots[END_REF], [START_REF] Khalil | Modeling, identification and control of robots[END_REF], or by numerical method based on the QR decomposition [START_REF] Gautier | Numerical calculation of the base inertial parameters[END_REF].

The minimal dynamic model can be written using the b n base dynamic parameters  as follows:

( ) idm q,q,q       (11)
where  is obtained from st  by eliminating the columns corresponding to the non identifiable parameters. Because of perturbations due to noise measurement and modelling errors, the actual force/torque  differs from idm τ by an error, e , such that:

( )

idm τ e q,q,q e          (12)
where:

g v     ( 13
)
v  is the ( 1) n  vector of the actual current references of the current amplifiers. Equation ( 12) represents the Inverse Dynamic Identification Model (IDIM).

C. Identification of the Dynamic Parameters

The off-line identification of the base dynamic parameters  is considered, given measured or estimated off-line data for τ and   q, q, q   , collected while the robot is tracking some planned trajectories. ( ) q, q, q   in [START_REF] Khalil | Symbolic calculation of the base inertial parameters of closed-loop robots[END_REF] are estimated with ( ) ˆq, q, q   , respectively, obtained by band-pass filtering the measure of q [START_REF] Gautier | Dynamic identification of robots with power model[END_REF]. The principle is to sample the identification model [START_REF] Khalil | Symbolic calculation of the base inertial parameters of closed-loop robots[END_REF] at a frequency m f in order to get an over-determined linear system of ( )

m n n  equations and b n unknowns such that:     fm f m f m ˆŶ τ W q,q,q χ ρ     (14)
In order to cancel the high frequency torque ripple in   fm Y τ and to window the identification frequency range into the model dynamics, a parallel decimation procedure low-pass filters in parallel fm Y and each column of fm W and resamples them at a lower rate, keeping one sample over d n . This parallel filtering procedure can be carried out with the Matlab decimate function [START_REF] Gautier | Dynamic identification of robots with power model[END_REF]. It is obtained:

    ˆŶ τ W q,q,q χ ρ     (15)
ρ is the ( x1) r vector of errors, with

m d r n n / n   ,   ˆŴ q, q, q   is the ( x ) b r n observation matrix.
In Y and W , the equations of each joint j are grouped together such that: The ordinary LS (OLS) solution χ minimizes ρ  .

        1 1 T T T T T T n n Y Y ... Y ,W W ... W               (16)
Using the base parameters and tracking "exciting" reference trajectories, a well conditioned matrix W is obtained. The LS solution χ of ( 15) is given by:

    1 T T χ W W W Y W Y     (17)
Standard deviations i   , are estimated assuming that W is a deterministic matrix and  , is a zero-mean additive independent Gaussian noise, with a covariance matrix C  , such that:

T 2 ( ) r C E ρρ I      (18)
E is the expectation operator and I r , the ( ) r r  identity matrix. An unbiased estimation of the standard deviation

  is: 2 2 ( ) Ŷ -W r b      (19)
The covariance matrix of the estimation error is given by:

T 2 T 1 [( )( ) ] ( ) ˆˆĈ E χ χ χ χ W W         . ( ) i 2 ˆĈ i,i      is the i th diagonal coefficient of ˆĈ  (20)
The relative standard deviation ri %   is given by:

100 ri i ˆˆi %       , for i  ≠ 0 (21)
The OLS can be improved by taking into account different standard deviations on joint j equations errors [START_REF] Gautier | Dynamic identification of robots with power model[END_REF]. Each equation of joint j in ( 15) is weighted with the inverse of the standard deviation of the error calculated from OLS solution of the equations of joint j , given by:

  j j j j Y τ W χ ρ   (22)
This weighting operation normalises the errors in ( 15) and gives the weighted LS (WLS) estimation of the parameters.

III. GLOBAL IDENTIFICATION OF THE JOINT DRIVE GAINS

A. Inverse Dynamic identification Model of the robot and the payload

The payload is considered as a link 1 n  fixed to the link n of the robot. The model [START_REF] Khalil | A new geometric notation for open and closed loop robots[END_REF] becomes:

  idm L L             (23)
Where:

L χ is the (10 1)  vector of the inertial parameters of the payload which can be calculated with a CAD model of the payload.

L  is the   10 n  jacobian matrix of idm τ , with respect to the vector L χ . The model ( 12) becomes:

  L L e              (24)
Taking into account that the joint j force/torque depends only on parameters j ,n  of links j to n , the IDIM joint j model is given by:

j ,n j j j j,n j j L L g v e                     ( 25
)
where j j ,n j L , ,    , are the respective non zero elements of line j in (24), j v  is the actual joint j current reference given by the numerical controller of the robot.

Equation ( 25) can be rewritten as:

1 j ,n j j j ,n j j j L L j g v e g g                     (26)
Equation ( 26) is the minimal model IDIM used to identify the drive gain of joint j .

B. Identification of the drive gains

Considering now that the robot has carried out two trajectories: (a) without the payload and (b) with the payload fixed to the end-effector, the sampling and filtering of the model IDIM (26) can be written as: 0

* j j , n j j j j j * * a j j j j b L L W / g V / gW W W 1/ g                       (27)
where: ( )

j j V v  is the vector of j v  samples, j a
W is the observation matrix of joint j in the unloaded case, j b W is the observation matrix of joint j in the loaded case, j L W is the observation matrix for joint j corresponding to the known payload inertial parameters, *  is the weighted LS solution of the over-determined system (27). The last coefficient of vector *  gives the estimation of the inverse of the joint j drive gain. Using this global approach and assuming that all the exciting payload dynamic parameters for each joints are known (using experimental measurements or computed with a CAD software), the drive gain identification can be carried out independently for each joint.

The advantage of this method is its simplicity compared with the usual measurement of each gain of the drive train elements. As the drive gains identification is decoupled, i.e. the identification is carried out for each joint independently, there is no cumulative estimation errors for the drive gains.

Its main drawback is that the identification requires the knowledge of all the payload inertial parameters that can be difficult to get with a sufficient accuracy, even using CAD software. But the experimental results obtained in section IV confirm the efficiency of the proposed approach.

IV. CASE STUDY

A. Description of the RX 90 kinematics

The Stäubli RX-90 robot (Fig. 1) has a serial structure with six rotational joints. Its kinematics is defined using the modified Denavit and Hartenberg notation (MDH) [START_REF] Khalil | A new geometric notation for open and closed loop robots[END_REF]. In this notation, the link j fixed frame is defined such that the j z axis is taken along joint j axis and the j x axis is along the common normal between j z and j 1 z  (Fig. 1). The geometric parameters defining the robot frames are given in Table 1. The payload is denoted as the link 7. The parameter 0 j   , means that joint j is rotational, j  and j d parameterize the angle and distance between j 1 z  and j z along j 1 x  , respectively, whereas j  and j r parameterize the angle and distance between j 1 x  and j x along j z , respectively. For link 7, 2 j   means that the link 7 is fixed on the link 6. Since all the joints are rotational then j  is the position variable of joint j . We note that the variables of joints 2 and 3 in our notation are obtained from the measured values of joints 2 and 3 using the relations, 2

2 2staubli /       2 3 3staubli /      .
The main advantage of using the MDH notation is that the identifiable inertial parameters can be determined symbolically using simple closed-form rules [START_REF] Gautier | Direct calculation of minimum set of inertial parameters of serial robots[END_REF], [START_REF] Khalil | Modeling, identification and control of robots[END_REF], [START_REF] Khalil | Symbolic calculation of the base inertial parameters of closed-loop robots[END_REF] and it is possible to generalize the proposed methods for tree structure robots and closed loop robots [START_REF] Hollerbach | Model Identification[END_REF], [START_REF] Khalil | Modeling, identification and control of robots[END_REF], [START_REF] Featherstone | Dynamics[END_REF], [START_REF] Gautier | Numerical calculation of the base inertial parameters[END_REF].

B. Identification of the drive gains

The method is validated using a calibrated payload (Fig. 2). Its mass has been measured with a weighting machine ( L M = 7.025 Kg ± 0.050 Kg). The other parameters have been estimated using CAD software. They are given in table 2. Their values are accurate due to the simplicity of the payload shape (Fig. 2).

For the RX 90, it can be shown that: -the torques of joints 1 to 3 depends on all payload parameters; -the torques of joints Then, the new drive gains for joints 1 to 4 are identified using our new approach (case 2). The identified drive gains are presented in table 4. It can be shown that there is an In table 5, the quality of identification for the two studied cases in terms of reduction of the norm of the error vector ρ is shown. With the identified gains, the quality of identification is largely improved.

Figure 3 shows that the reconstructed torques with the identified gains for joints 1 to 4 are close to the actual ones.

In order to show the sensitivity of the identification results to the variation of the drive gain values, a 20% error on these gains (with respect to those identified in case 2) is simulated (case 3). The quality of identification is clearly worst in this case (table 5). The identified payload parameters are given in table 3. Their values are far from reality.

Finally, in order to validate the new drive gain values, a new payload of 4.275 Kg ± 0.050 Kg, is identified. The identified mass is L M = 4.34 Kg, with 2 i   = 0.149, which is very close to the exact value. It is possible to conclude that the drive gains have been well identified.

V. CONCLUSION

This paper has presented a new method for the identification of the total drive gains for robot joints in one step. This is a method easy to implement which does not need any special test or measurement on elements inside the joint drive train. It is based on a IDIM-LS technique using current reference and position sampled data while the robot is tracking one reference trajectory without load fixed on the robot and one trajectory with a known payload fixed on the robot, whose inertial parameters are measured or calculated by a CAD model. The method has been experimentally validated on a Stäubli RX-90 robot. With the new identified drive gains, the identification of the total dynamic model of the robot has been improved and another payload has been accurately identified. This shows the effectiveness of the method.

  , MZ are the 3 components of the first moment of link j . j M is the mass of link j , j Ia is a total inertia moment for rotor and gears of actuator j . j Fv , j Fc , j off  are the coefficients of the v F , c F , off  , matrices, respectively.

W

  represent the m d n / n equations of joint j .

4 Fig. 1 .

 41 Fig. 1. Link frames of the RX-90 robotTABLE I GEOMETRIC PARAMETERS OF THE RX-90 ROBOT WITH THE PAYLOAD j  j  j

  4 and 5 do not depend on parameter L M ; -the torques values of joint 6 do not depend on parameters L M and L MZ . Moreover, our experiments have shown that parameters are too small to allow the correct calibration of the drive gains of joint 5 and 6. Therefore, the drive gains of joints 1 to 4 only are identified. In order to identify the other drive gains, more exciting payload should be designed, especially with large values for parameters L MX and L MY .Staübli's control software provides directly the value of the joint torque. Therefore, it is proceeded to a calibration of the drive gains compared with a priori values used in the controller which are taken equal to 1. A first identification of both the RX 90 dynamic parameters and the load inertial parameters is achieved using the IDIM-LS method with the model (25), as detailed in[START_REF] Khalil | Identification of the payload inertial parameters of industrial manipulators[END_REF] (case 1). The results are presented in table 3. It can be observed that the mass L M and the first moment L MZ are poorly identified with the manufacturer's gains.

Fig. 3 .

 3 Fig.3. Torque validation with identified drive gains, Case 2.

Fig. 2 .

 2 Fig. 2. The 7.025 Kg payload

  :

	idm v  is the ( 1) n 	vector of current references of the
	current amplifiers,
	v τ0 is the ( 1) n  vector of amplifier offsets,
	g τ is the ( n n  matrix of the drive gains, )
	i t g NG K  				(4)
	N is the ( n n 	)	gear ratios matrix of the joint drive
	chains such as, m q 		Nq  , with m q  the ( 1) n 	vector of
	motor velocities on the motor side,
	i G is the ( n n 	)	static gains diagonal matrix of the
	current amplifiers, assuming a large current loop bandwidth,
	usually close to 1KHz.
	t K is the ( n n 	)		diagonal matrix of the electromagnetic
	motor torque constants.

TABLE 5

 5

							TABLE 3
		IDENTIFICATION OF THE PAYLOAD DYNAMIC PARAMETERS.
		Parameter A priori value Estimated value 2   i 		ri % 
	Case 1	MZ L			1.003	0.957	4.36e-2	2.3%
		M L			7.025	6.61	1.55e-1	1.2%
	Case 2	MZ L			1.003	0.936	4.18e-2	2.2%
		M L			7.025	7.09	1.55e-1	1.1%
	Case 3	MZ L			1.003	0.837	4.41e-2	2.6%
		M L			7.025	8.61	1.96e-1	1.1%
	i  is the standard deviation (20) and 	ri %   its relative value (21)
							TABLE 4
						IDENTIFIED DRIVE GAINS.
					Joint 1 Joint 2 Joint 3 Joint 4
		j g 		0.9238 0.9533 1.1032 0.8737
		2   8.03e-2 3.55e-2 3.68e-2 4.81e-2 i 
				ri % 	4.3%	1.9%	1.7%	2.8%
					QUALITY OF IDENTIFICATION.
		Error norm 	Relative Error norm Ŷ 	ˆ 
	Case 1	77.28			0.075643	1.16293
	Case 2	76.26			0.0746481	1.14648
	Case 3	86.35			0.0845219	1.3003
	ˆŶ W   		is the minimal norm of error, ˆ  is given by (19).


• Produced database and cross-map.

•CEReS•••

  Adaption of Comet's pyrolysis process to PCB -Initial orientative bench-scale tests -Pilot scale tests -Phoenix reactor • Produced and characterised char and hydrocarbon outputs -Hydrocarbon -18% original mass; requires removal of Si prior to use in co-generation engine -Charfour size categories: +8 mm; 2-8 mm; 75 µm-2 mm; -75 µm -Bromine -74% recovered in quench water • >90% Ag, Au and Cu concentrated in two fractions -~47% input; concentration factor ~2 Adaption of Comet's pyrolysis process to PCB -Initial orientative bench-scale tests -Pilot scale tests -Phoenix reactor • Produced and characterised char and hydrocarbon outputs -Hydrocarbon -18% original mass; requires removal of Si prior to use in co-generation engine -Charfour size categories: +8 mm; 2-8 mm; 75 µm-2 mm; Selection of bioleaching consortia -Two enriched from Janina waste, two existing bioleaching consortia -Selected 48°C "TW48" consortium from Janina waste • Biolixiviant capable of leaching char • Huge problem with corrosion -caused by Cl content of waste… • Decide on appropriate strategy washing or using different waste? Selection of bioleaching consortia -Two enriched from Janina waste, two existing bioleaching consortia -Selected 48°C "TW48" consortium from Janina waste • Initial tests show biolixiviant capable of leaching char • Huge problem with corrosion -caused by Cl content of waste… • Decide on appropriate strategy washing or using different waste?

Integrated Process Concept PCB AMD generating Coal Production Wastes Preprocessing & catalytic cracking Lixiviant production (H + , Fe 3+ ) Char leaching PLS for metal recovery Solid residue processing Raffinate Non-AMD generating solids Existing E- waste handling programme Coal Mine site Liquid hydrocarbon fuels, halogen brine, C- rich concentrates Metal (Cu) and separate concentrates (precious metals, Sn-Pb, REEs, Ga and Ta) Benign residues, Fe-based products i ii iii

		Introduction Introduction CEReS :: Progress
	Raw Materials	Coal Mine site
	• Received funding from the Research Fund for Coal and -Selected Janina Spiral Tails (~12% Py) • Use Poland as case study region AMD generating Coal • Analysed four coal waste streams Production Wastes
	Steel (RFCS) -Detailed mineralogical characterisation	PCB
	Existing E-• Prove technical feasibility of individual unit processes waste handling
	• Three-year R&D project from June 2016; project budget • Bioprospecting -At (mini) pilot-scale programme
	~3.2M€ • Eight partners from five countries -Two bioleaching consortia (30°C & 48°C) • Integrate through modelling and simulation -Demonstrate viability of concept -Microbial ecology of Janina waste dump (NGS)
	-University of Exeter (Coordinator) -BRGM; Caspeo -Include economic assessment • Characterised AMD-generating potential -Université de Liège; Comet Traitements • Evaluate environmental benefit vs "do nothing" scenario UK FR BE -"Missing Acid" question…
	-GIG; TAURON Wydobycie -LCA	PL
	-University of Cape Town	RSA
		ceres.biohydromet.net

i ceres.biohydromet.net Raw Materials • Analysed WEEE processing in Poland -Three PCB categories; selected low-grade • Not as required; ~35% PCB use Comet own stock • Detailed characterisation & analysis -Metal/value breakdown -Au and PM highest values; Cu significant.

rd October, 2019 Fukuoka, Japan

Process Integration

• Process simulator compiled in USIM PAC