Dynamic parameters identiﬁcation of a haptic interface for a helicopter ﬂight simulator

. The haptic interface force feedback is one of the key factors for a reliable ﬂight simulation. This paper addresses the design and control implementation of a simple joystick-like haptic interface to be used for a helicopter ﬂight simulator. The expression of the haptic interface force is obtained by dynamic analysis of the haptic interface operation. This paper proposes a new strategy aiming at avoiding the use of an expansive and complex force/torque sensor. Accordingly, speciﬁc dynamic model is implemented by including Stribeck friction to describe the friction moment. Experimental data are processed as based on a genetic algorithm for identifying the dynamic parameters in the Stribeck friction model. This allows to obtain the friction moment parameters of the haptic interface, as well as the torque distribution due to gravity and the rotational inertia parameters of the haptic interface for the calculation of the haptic interface force. Experimental tests are carried out and results are used to validate the proposed dynamic model and dynamic parameter identiﬁcation method and demonstrate the effectiveness of the proposed force feedback while using a cheap photoelectric sensor instead of an expansive force/torque sensor.


Introduction
Due to the high risks and costs of training on a real helicopter, pilots often need to complete a preliminary training on a helicopter flight simulator. Real flight training is carried out only after a successful training on a simulator. In helicopter flight simulators, the joystick is a key element, which is used to achieving a reliable helicopter driving while providing a force feedback feeling, which needs to properly mimic the force feedback feelings of a real helicopter flying (Prachyabrued and Robert, 2018). In the real flight control process, the pilot's hand will feel the feedback force and make judgments, then adjust the flight status through the feedback force. In the whole process, the current force of the haptic interface needs to be used as the simulation calculation of the feedback force. Thus, its measurement accuracy is particularly important. Most of the operating forces are measured by using torque sensors. The papers Cavallo et al. (2004), Rinaldi and Beckham (1983), Hutton et al. (1994), and Fergani et al. (2016) installed a pull-up/pressure sensor at the end of the rod to read the applied force at each moment. However, torque sensor occupies large space in both mechanical installation and wiring, which affects the normal operation of the pilots. This paper proposes a different strategy aiming at avoiding the use of an expansive and complex force/torque sensor. Namely, the proposed approach is based on establishing a reliable dynamic model of the haptic interface to replace the role of the force/torque sensor. Accordingly, a proper dynamic model and a parameter identification strategy have been proposed. In a servo system, there are complex rolling friction and sliding friction phenomena between gear train, bearing, input/output shaft and seal. Therefore, the modelling of friction torque is very important, which directly affects the dynamic accuracy of the haptic interface. A variety of empirical models have been proposed for friction calculation of servo system, which can be divided into two phases: static friction model and dynamic friction model (Bona and Indri, 2005).
The most widely used static friction model is the Stribeck model. Several papers analyzed the friction of a servo system by establishing its Stribeck friction model (Lichun and Pavelescu, 1982;Xu et al., 2011;Chen et al., 2016;Márton and Lantos, 2009;Broel-Plater et al., 2018). Considering the advantages and drawbacks of the existing models in literature, it has been decided to implement a Lugre model to analyze the friction of the servo system, as proposed in Canudas de Wit et al. (1995), Wang et al. (2019), Azizi andYazdizadeh (2019), Freidovich et al. (2010), and Ishikawa et al. (2010). Main advantage of Lugre model is its relatively high accuracy, while its implementation is relatively complex, since it introduces an unknown quantity which can be difficult to be measured directly and requires a high precision servo system. The Stribeck model is currently used in many fields to model the friction moment. Since there are unknown parameters in the haptic interface dynamics model, such as gravity moment, Stribeck friction moment model parameters and rotational inertia, these parameters need to be identified. The main available parameter identification methods are tracing methods, least square methods, and genetic algorithms. The genetic algorithms are very effective evolutionary random search methods, which are proven in literature to be very effective in managing non-linear problems with a high level of computational robustness and calculation parallelization. Accordingly, genetic algorithms can be seen as the most convenient approach for addressing this identification problem. In this paper, a dynamic analysis is carried out on the a 2-dof flight haptic interface for a helicopter simulator. A nonlinear Stribeck model is established for the friction torque of the servo system. A curve fitting and genetic algorithm method are used for processing the data to identify the gravity torque, four parameters of the Stribeck model and the rotational inertia. Finally, the available parameters identified are brought into the dynamic model to calculate the haptic interface force in real time and an experimental validation is carried out to proof the effectiveness of the proposed approach.

The proposed haptic interface and its operating principles
This paper adopts the helicopter simulator at laboratory in Yanshan University which is shown in Fig. 1. The helicopter flight simulator is mainly composed of four modules: a helicopter simulation cockpit system, a virtual scene system, a dynamic motion system and an audio system. The simulator's haptic interface is shown in Fig. 2. The haptic interface is driven by two servomotors, and two disc speed reducers are used for speed regulation. The motor adopts Kollmorge CKM04 AC permanent magnet synchronous servo motors, and its corresponding AKD series servo driver supports Modubus TCP communication. The parameters of the servo mo- tors and the reducers are shown in Tables 1 and 2. The hardware system composition of the haptic interface is shown in Fig. 3. First, the computer is connected to the switch, and the switch is connected with two motor drivers respectively through Enternet port. Then, C + + programmed code has been written to establish communication between servo motor and computer through Modbus TCP protocol. When the driver operates the haptic interface, the computer can obtain the current electricity value from the servo motor at each time. Then, the driving torque of the servo motor is obtained. At the same time, the photoelectric encoder of the servo motor can measure the stick's position and speed at the current time, which can be used for the dynamic calculation of the haptic interface.

Dynamic analysis of the haptic interface
First, coordinate system for the haptic interface is established for the haptic interface, as shown in Fig. 4. The intersection point of the two motor axes is defined as origin O. The two motor axes are set as x axis and y axis respectively. z axis is perpendicular to the xy plane.
The free-body diagram is shown in Fig. 5. The angle θ , θ x , θ y formed between the haptic interface and the z axis is the output angle of reducer 1 and reducer 2 respectively. the relationship of the three parameters can be expressed as follows: When the haptic interface moves, the control handle, the haptic interface and the base of the haptic interface move at an angle θ , the motor 1, the reducer 1 and the rotating mechanism move at an angle θ y .
In Fig. 5, m is the mass of the haptic interface system, l G is the distance between the center of mass of the haptic interface and the origin, and l is the distance between the top of the stick and the origin. τ a1 and τ a2 are the output torques   by motor 1 and 2 through the disc reducer respectively; τ f1 and τ f2 are the friction torques on the x axis and y axis of the servo system respectively. The manual force of the pilot in the direction of x axis, y axis and z axis defined as F x , F y , F z . The Lagrange method is adopted to conduct dynamic analysis of the system. The Lagrange function can be formulated as In Eq.
(2), E k is the kinetic energy of the haptic interface system; E p is the potential energy of the haptic interface system;q 1 ,q 2 in the haptic interface system represent the output angles θ x , θ y of the two reducers respectively.  The total kinetic energy E k for the haptic interface system can be expressed as: In Eq. (3), J 1 is the sum of the moment of inertia of motor 1, reducer 1 and rotating mechanism; J 2 is the sum of the moment of inertia of the control handle, haptic interface and haptic interface base. The potential energy of the haptic interface E p can be expressed as: In combination with Fig. 5 and Eq. (1), the following equation can be obtained: According to the comprehensive Eqs.
(2)-(5), can be obtained: In Eq. (6), τ 1 and τ 2 are the torques exerted by the pilot on the x axis and y axis during the control process respectively. Especially, the torque output τ a by the motor through the disc reducer is related to the reduction ratio i of the reducer and the output torque T a of the motor, as: T a can be deduced from Eq. (7): In Eq. (8), P n is polar logarithm; ψ r is the flux chain of rotor magnetic pole; L d , L q is the armature inductance of axis d and axis q; I d , I q is the current of d axis and q axis.
Since the armature inductance of d axis and q axis of the motor used same value, namely: From Eqs. (7), (8) and (9) can be obtained as follows: At this time, the torque and the q axis current have a linear relationship, which can be expressed as: Combined with Eqs. (6) and (11), and using mgl G replaced by the gravity moment T G , the following equation can be obtained: . (12) Thus, the manual force of the pilot can be expressed as:

Friction torque model of the haptic interface
In the haptic interface structure, there are complex frictions within the servo motor and between the reducer. The friction torque generated by these frictions has a great influence on the accurate modelling of the haptic interface. Therefore, the friction torque at each moment must be accurately obtained to ensure the accuracy of the calculated force and torque at each moment.
The classic Coulomb friction model friction indicated that friction force is related to the positive pressure acting on the object. However, due to the lubrication in the servo system, the influence of viscous friction must be considered. Stribeck proposed to model the servo friction as a combination of Coulomb friction and viscous friction. He has given a quantitative formulation for friction with respect to velocity (Balogh and Krstic, 2004). Similarly, the Stribeck formula can be converted into the formula of friction torque regarding motor speed as proposed in (Iwasaki et al., 1999): In Eq. (14), ω is the motor speed; τ c is Coulomb friction torque; τ s is the maximum static friction torque; ω s is Stribeck speed; B ω is the coefficient of viscous friction. The Stribeck curve is related to the speed of the servo system. In the low-speed area, it is affected by Coulomb friction and viscous friction, and the friction torque decreases. This area is called boundary lubrication friction area. With the increase of speed, the influence of viscous friction is much higher than Coulomb friction. At this time, the friction torque is approximately proportional to the rotational speed. This operation condition is defined as liquid lubrication condition.
There are two problems in applying the Stribeck friction model to the control of haptic interface force control: on the one hand, it is difficult to describe the change of haptic interface friction torque when the velocity is zero; on the other hand, the change of velocity near the zero point will cause a sudden change of friction torque, which will lead to a jitter in the reversing process. Therefore, approximate treatment of the low speed phase of the Stribeck friction model is required. This paper proposes to use a sigmoid function to deal with Stribeck friction model (Ciliza and Tomizuka, 2007).
The Stribeck friction torque model improved by sigmoid function can be expressed as follows: According to the requirements of the model and accuracy, the value of γ is defined as 500 in the proposed model.

Identification of dynamic parameters
3.1 Identification of the gravity torque model parameter According to Eq. (13), the total force applied to the haptic interface is directly related to the magnitude of the gravity torque T G . In order to obtain an accurate model, it is necessary to identify the gravity torque T G . This is achieved with an experiment conducted on the motor 1 when its movement is on xz plane. Set the haptic interface to rotate at a constant speed, and thenθ x is zero. The haptic interface was not moved in the y direction, so θ y is zero. The driver has no additional force to the haptic interface, accordingly also τ 2 is zero. Since the velocity is constant, τ f1 is constant. In addition, sin θ x is relatively small in practical application, accordingly θ x can be used instead. The experimental process is shown in Fig. 6 and Eq. (6) can be converted into: At this time, θ x maintain a liner relationship with I q1 have a linear relationship. Discretizing Eq. (16) one can obtain:  A following experiment consists in setting the motor speed at 40 rpm (Rotation per minute), make the haptic interface rotate at a constant speed, take the values of θ x and I q1 in multiple groups. A curve fitting of the results in such experiments is shown in Fig. 7. The identified gravity torque is 11.14 N m. The toolbox of Curve Fitting Tool in Matlab is used to get the result.

Identification of friction torque model parameters via genetic algorithm
When we choose Stribeck friction model, it contains four unknowns: τ c , τ s , ω s and B ω , these four parameters need to be identified. Similarly, the motor 1 is set as the speed mode and a specific motor speed is given to make the haptic interface rotate at a constant speed. At the same time the driver has no additional force to the haptic interface, in this case, θ y , θ x and τ 1 are equal to zero. Therefore, according to Eq. (6), the following equation can be obtained: According to Sect. 2.1, the current value, rotation angle and rotation speed of the motor can be read at each moment. It is can be concluded from Eq. (18), the friction torque τ f1 at the corresponding rotational speedθ x can be obtained. The value of friction torque τ f1 can be obtained by the different output angleθ x .   Genetic algorithm is an algorithm designed according to Darwinian evolution theory: through continuous natural selection, survival of the fittest, so as to obtain the most adaptable results (Liu, 2006). The process of parameter identification mainly includes generation of initial population, selection of race parameters, setting of fitness function, selection, crossover, variation and termination.
The search scope of parameters τ c , τ s , ω s and B ω can obtained by tracing method. The iteration interval is [4-5], [4][5], [0-0.1] and [20-40] respectively. The population size was set as 100, and the crossover probability and mutation probability were set as 0.4 and 0.1 respectively. The selection of GA parameters has been carried out by considering the peculiarities of the specific identification problem as well as several preliminary identification tests. In particular, after preliminary identification tests we selected 0.4 as crossover parameter and 0.1 as variation parameters. The selected crossover parameter gives a preference to parental individuals on filial individuals. The selected variation value is motivated by a limited influence of variation on this specific parameter identification. Numerical simulations have been carried out also with other values of the crossover parameter (e.g. 0.4 and 0.8). Results led to similar outcomes of the identification procedure.
The multi-group friction torques obtained through Eq. (18) can be expressed as: In Eq. (8),τ c ,τ s ,ω s andB ω are the four parameters to be identified respectively, andτ f1 is the theoretical friction torque calculated by the motor speed. Take the target function: The actual identification task of genetic algorithm is to solve the minimum value of objective function. The minimum target value is defined as 10 −6 here, and the maximum evolutionary algebra is 10 000 generations. In each generation, the population is sorted in ascending fitness order, and the optimal solution of this generation is recorded as the elite. Parental individuals were selected in the fitness sequence for crossover, and the offspring were mutated. Finally, when the fitness of a certain generation of elite meets the requirements (J ≤ 10 −6 ) or the evolutionary algebra is more than 10 000 generations, the whole selection process of genetic algorithm is terminated.
The process of genetic algorithm is shown in Fig. 8, and its identification process is shown in Fig. 9. Figure 10 shows the comparison between the theoretical friction moment calculated by the identified parameters and the actual friction moment. The measurement result of the friction torque at the moment ofθ x < 0 is the same as the moment ofθ x > 0, but the symbol is opposite, so it was omitted. Here we only figure out the absolute value of the friction torque. The identified parameters at the end of program operation are shown in Table 1. The identification of friction model parameters in the y axis direction was the same as that in the x axis. The identification results are listed in Table 3. Several cases have been considered in the identification process for the friction torque model parameters to obtain a reliable mean value of all the parameters. Considering the characteristics of the proposed identification process we can estimate a confidence interval of 0.0001 % (the minimum target value is in Eq. 20), which can be considered suitable for the specific application, also in comparison with commonly available accuracy of commercial sensors.

Identification of rotational inertia model parameters
According to Eq. (14), the Stribeck friction torque model can be simplified as: Firstly, the motor 1 is analyzed. When the pilot does not operate the haptic interface, Eq. (6) can be converted into: Thus, the structure block diagram of Eq. (21) can be described and is shown as in Fig. 11. After sorting out Fig. 11, the simplified structure block diagram of motor driven haptic interface is shown in Fig. 12.
According to Fig. 12, the continuous function can be expressed as: The zero order holder method is adopted for Z transformation in Eq. (22), and the following equation can be obtained: In which, By discretizing Eq. (24), the following equation can be obtained: In which, In Eq. (27), τ c+ represents positive Coulomb friction and τ c− represents negative Coulomb friction. Set the torque mode through the motor driver, set the specific current size, and drive the haptic interface to rotate. When the motor speed reaches a certain value, take multiple groups of current value I q1 and rotation angle θ x .
For the obtained data, the genetic algorithm in Sect. 3.2 can be used to identify the moment of inertia. The estimated J 2 value among search range [0.1-1.1]. The population size was set as 100, and the crossover probability and mutation probability were set as 0.4 and 0.1 respectively. The selection of GA parameters has been carried out by considering the peculiarities of the specific identification problem as well as several preliminary identification tests. In particular, the selected crossover parameter has been defined to give a preference to parental individuals on filial individuals. The selected low variation value is motivated by a limited influence of variation on this specific parameter identification.
The multiple groups of rotation anglesθ x (k) calculated by Eq. (26) can be expressed as: The minimum target value is given as 10 −4 , and the maximum evolution algebra is 10 000 generations. The evolution of its optimal solution is shown in Fig. 13. The J 2 of the recognized moment of inertia is 0.2887 kg m 2 . Performming the same operation on the motor 2 as above, and the identified moment of inertia J 1 + J 2 is 0.2989 kg m 2 . Therefore, J 1 is 0.0102 kg m 2 .

Experimental validation
Based on the experimental platform in Fig. 2, the driver is connected to the switch through the X11 terminal, and the switch is connected to the computer. Use Visual Studio C++ for programming. The computer communicates with the driver via Modbus TCP. Specifically, the program adopts the dynamic mapping function of Modbus. Modbus allows the driver to map any fixed register address to a new register address, and read-write access to the remap parameters by reordering the sequence block. The computer gives instructions to read or change the motor current state by reading or modifying the value in the register.
In this paper, the sampling time of the motor is given as 15 ms, and the quadrature axis current I q1 and I q2 size are set as 0. After noise removal, the rotation speed ω 1 , ω 2 and position θ m1 , θ m2 of the two motors at each moment were measured. ω 1 and θ m1 can be converted to the reducer 1. the reducer 1 output angular speedθ x and position θ x are corresponding to them respectively; ω 2 and θ m2 can be converted to the reducer 2. The reducer 2 output angular speedθ y and position θ y are corresponding to them respectively; by substituting ω 1 and ω 2 into Eq. (14), the friction torque τ f1 and τ f2 of two degrees of freedom on the haptic interface at the current moment can be obtained. The angular accelerationθ x , θ y of two degrees of freedom of the haptic interface can be calculated by difference method, and they can be expressed as: In Eq. (29), T is the sampling time of the motor. At this time, the driver's force on the haptic interface can be obtained through Eq. (13).
In the experiment, the pilot manipulates the haptic interface diagonally and nearly sinusoidal, as shown in Fig. 14a, and records the control force after the kinetic solution. The calculation result of haptic interface force in the whole process is shown in Fig. 14b.
Then take the take-off stage and the helicopter hovering to the left for example.
For the control process in the take-off stage, after lifting the total pitch bar, the pilot slowly and uniformly pulls the haptic interface backward. The haptic interface needs to be operated negatively along the x axis to make the helicopter climb steadily, and holds the haptic interface steadily after pulling back to a certain degree. When the helicopter climbs to the target altitude, lay down the haptic interface. The haptic interface needs to be operated forward along the x axis to make the helicopter attitude. The specific control process is shown in Fig. 15. From Fig. 15b, it can be seen that the pilot in t 1 moment is pulling the haptic interface in a negative direction towards the x axis and F x , F z both are negative. In the t 2 stage, the pilot holds the haptic interface for a certain seconds. At this time, because the motor speed is zero, F x and F z only equal to gravity. F x , F z are all constant positive values. In the t 3 stage, the pilot pulls the haptic interface forward along the x axis, therefore, both F x , F z are positive. At the t 4 stage, the haptic interface returns to its original position, and the values of F x , F z back zero. In addition, because the driver did not apply force in the direction of the y axis, F y was constant zero in the whole process.
For the control process of helicopter hovering to the left, the pilot slowly and uniformly pulls the haptic interface. The haptic interface needs to be operated forward along the y axis so that the helicopter fuselage tilts to the left. After pulling the haptic interface to a certain extent, the pilot stabilizes the haptic interface. At this time, with the operation of the helicopter pedal, the helicopter can make circular circling. The specific control process is shown in Fig. 16. As can be seen from Fig. 16b, in the t 1 stage, the pilot pulls the haptic interface forward along the y axis, which F y is positive and F z negative. In the t 2 stage, the pilot holds the haptic interface firmly, which F y is negative value and F z positive value. F x was constant zero in the whole process.
It can be seen from the experiments that the control process curve of take-off stage and left hover stage is consistent   with the actual driving situation. It can be indicated that proposed identification method for the haptic interface model parameters are feasible and can be applied in the helicopter flight simulator.

Conclusion
In this paper we have addressed the dynamic modelling and parameter identification of a joystick-like haptic interface to be used for helicopter flight simulators. Careful attention has been addressed at the Stribeck friction model and its parameters for an accurate modelling of servomotors. A specific genetic algorithm has been implemented for an accurate identification of the Stribeck friction parameters. A full dynamic model has been developed and implemented as based on Lagrange formulation. Experimental tests show that the proposed identification procedure and haptic interface allow a control of both take-off and left hover stage with results being consistent with the actual helicopter driving. Accordingly, the proposed haptic interface can provide a reliable force feedback even without using expansive force sensors. As future work, we aim at improving the proposed model for a more accurate modelling of friction especially when the velocity is close to zero. Future activities will also include further experimental testing. Data availability. All data used in this paper can be obtained on request from the corresponding author.
Author contributions. DZ and JZ wrote the whole paper, HY and TN designed the experiment and dealt with data, GC and SY revised the paper.
Competing interests. The authors declare that they have no conflict of interest.