Maciej Bednarczyk

Keywords: Robotique collaborative, interactions physiques Homme-robot, commande en impédance, procédures médicales assistées par robot Collaborative robotics, physical human-robot interaction, impedance control, robot-assisted medical procedures

Advanced control of interactions in robotics.

Application to robot-assisted medical procedures.

Résumé

En raison de l'intérêt croissant pour l'utilisation de systèmes robotiques dans un espace de travail partagé avec des opérateurs humains, le développement de robots collaboratifs met l'interaction Homme-robot au centre des préoccupations des roboticiens. Pour cette raison, le développement de nouveaux outils de contrôle permettant la gestion des interactions est devenu un sujet de recherche important. Ainsi, la conception de solutions améliorant la dynamique d'interaction et garantissant l'intégrité de l'opérateur est d'un intérêt particulier. Dans cette thèse, plusieurs outils de contrôle pour la robotique collaborative sont proposés. Les problématiques abordées visent notamment à garantir simultanément la compliance des robots tout en gérant des contraintes, ou à modifier la dynamique d'interaction de manière sûre. L'utilisation de bio-signaux afin d'améliorer la collaboration Homme-robot est également étudiée, pour évaluer l'intention de l'utilisateur. Cet ensemble de problématiques conduit à la conception de contrôleurs dédiés. Deux preuves de concept d'applications médicales utilisant les outils proposés sont développées pour l'insertion autonome d'aiguilles en radiologie interventionnelle et pour la rééducation bimanuelle.

List of Figures

Chapter 1 Introduction

Like many revolutionary technological advances, the story of modern robotics began with the dream of a visionary. Ahead of his time, in 1920, a Czech writer Karol Čapek imagined in his science-fiction play entitled R.U.R (Rossum's Universal Robots) human like artificial workers created to take over the world's production by making it cheaper and faster [ Čapek, 1920]. In this play, the concept of robot was introduced for the first time, coming from the Czech word robota meaning forced labor. Čapek's vision originated many research efforts and has led, sixty years later, to the successful implementation of robotic devices in many areas and especially in industry. Robotic manipulators are now widespread in factories, where they replace human workers in many tasks as they are cheaper, stronger, more precise and do not get tired, making production more efficient.

Industrial robots generally evolve in close and static environments. Elementary control architectures, based on position control, without environment perception, make it dangerous for human workers to operate in a workspace shared with robots. For this reason, in most cases, a strict separation between robots and human workers is implemented in factories. The robots can execute their tasks inside a closed area, preventing any human intrusion, which would lead to an immediate shut-down of the system. Even if separating robots from humans is an efficient way to ensure security and avoid physical contacts, restricting areas in factories becomes very constraining and space-consuming. Additionally, many applications would benefit from the combination of the robot performance in terms of precision, repeatability and load capacity with the cognitive capacities of the human operator such as adaptability and fast judgment in case of unpredicted events. However, removing the separations and permitting the sharing of the same workspace by the human and the robot leads to a complete change of paradigm and requires to rethink the technology. By opening the workspace, enhanced environment awareness and dynamic control strategies are required in order to ensure safety. The recent evolution of robotics tends towards a more human oriented design of robotic mechanisms and control strategies, leading to a new class of robotic systems known as cobots [Colgate et al., 1996], for collaborative robots.

Collaborative robotics

Technological advances in mechanical system design, in the development of electronics for sensors and low level controllers, together with components cost reduction, has led to the emergence of collaborative robots in industry. These robots have the ability to combine the advantages of a robotic system with human cognitive capacities. Designed with low inertia, smooth external envelopes and power/speed limitations, as well as dedicated control architectures featuring many external information such as spacial perception or force sensing, such systems can better guarantee the operator integrity. These safety guarantees make the contact between the robot and the operator not only allowed but also desired. Driven by the concepts of Industry 4.0, often referred to as "smart factory" [START_REF] Liao | Past, present and future of industry 4.0 -a systematic literature review and research agenda proposal[END_REF], cobots become a great source of innovation and attract much attention. Some examples of commercial cobots are depicted in Figure 1.1.

Different types of application

In modern industry, as pointed out in [START_REF] Matthias | Safety of collaborative industrial robots: Certification possibilities for a collaborative assembly robot concept[END_REF], production flexibility, i.e. the ability to easily reconfigure production lines, is of particular importance, especially in areas in which products have a short lifetime and frequent production changeover is required. This flexibility is directly related to the ease of deployment and programming of robotic solutions. In fact, cobots come with simplified programming tools compared to standard solutions, which often require the intervention of specialized integrators. Cobots can be programmed using intuitive and interactive human-machine interfaces allowing the operator to teach the robot some desired task by direct interaction [START_REF] Villani | Survey on human-robot collaboration in industrial settings: Safety, intuitive interfaces and applications[END_REF]. This allows the operator to teach the robot to take over, completely or partially, the most repetitive parts of the performed task and focus on more complex ones. In this context, cobots are of particular interest, especially due to their attractive price compared to standard robotic manipulators. In [START_REF] Giraud | Robotique collaborative: Évaluation des fonctions de sécurité et retour d'expérience des travailleurs, utilisateurs et intégrateurs au Québec[END_REF], it is reported that cobotic solutions are not always deployed for the need of human-robot collaboration (HRC), but also for space saving and economical reasons. Human-robot collaborative applications in manufacturing range from object handling and transportation, welding to hand guided assembly operations [START_REF] Villani | Survey on human-robot collaboration in industrial settings: Safety, intuitive interfaces and applications[END_REF].

In the field of medical robotics, cobots attract much attention. Currently, medical cobotic solutions are focused on two main areas: robot-assisted medical interventions and functional rehabilitation [Unitec, 2018, Taylor et al., 2016, Van der Loos et al., 2016]. In the first case, cobots are used to assist the practitioner during the procedure by enhancing his capacities in terms of stability and precision [START_REF] Taylor | Medical Robotics and Computer-Integrated Surgery[END_REF]. This results in reducing the risks for the patient due to positioning errors and reduced repeatability often due to the practitioners fatigue [START_REF] Hagag | RIO: Robotic-Arm Interactive Orthopedic System MAKOplasty: User Interactive Haptic Orthopedic Robotics[END_REF]. In this context, two examples of cobot applications coming from French companies can be cited in this field. The first is the ROSA robot from Medtech, shown in Figure 1.2a, which aim is to precisely guide the practitioner movement during minimally-invasive procedures, such as biopsies and surgical needle placement. The other is the TMS-Cobot from Axilum Robotics, shown in Figure 1.2b, which is designed to precisely position a TMS-coil for transcranial magnetic stimulation (TMS) and ensure its contact with the patient's head. This solution combines the robot strength and precision with environmental awareness and interaction management.

The second application of cobots in the medical field is functional rehabilitation. This type of procedures address patients with partial limb impairment due to strokes or spinal cord injuries [START_REF] Gassert | Rehabilitation robots for the treatment of sensorimotor deficits: a neurophysiological perspective[END_REF]. In classical rehabilitation scenarios, a therapist assists the patient in his movements, which is time consuming and often results in shortening the duration of the rehabilitation [START_REF] Richards | Driving motor recovery after stroke[END_REF]. Rehabilitation tasks are well suited for automation as they generally require many repetitions of a given movement [START_REF] Basteris | Training modalities in robot-mediated upper limb rehabilitation in stroke: a framework for classification based on a systematic review[END_REF]. Robot-assisted rehabilitation is an active research topic [Oña et al., 2018, Gassert and[START_REF] Gassert | [END_REF]. More recently, some research effort has also been dedicated to investigate the usage of industrial cobots for rehabilitation therapy in [START_REF] Kyrkjebo | Feasibility of the ur5 industrial robot for robotic rehabilitation of the upper limbs after stroke[END_REF]. The role of the robotic solution is then to assist the movements of the impaired limb when the patient is not capable of executing the desired rehabilitation task on his own and limit the extend of the movement in order to prevent injuries. In contrary to classical, task-specific robotic rehabilitation devices, cobots offer an interesting alternative, thanks to their attractive price, their versatility and the capacity to work closely with humans.

Different types of collaboration

With the growing number of applications featuring the cobot and the operator, the question of characterizing human-robot collaboration arises. In [De Luca and Flacco, 2012], the authors propose a nested framework that brings together the concepts of collaboration, coexistence and safety. In this framework, ensuring safety inherently is the most important condition for a cobot to satisfy, so that it can work close to humans. It consists in avoiding contact with the operator, and if collisions occur, reduce as much as possible the damages. Building upon this concept, coexistence gives the possibility to share a common workspace between the robot and the operator, but without any interaction. Finally, collaboration builds upon the two previous concepts and is defined as the robot ability to directly interact with the user in order to reach a common goal. It is pointed out that this collaboration can be performed in two ways, separately or simultaneously. Collaboration can be physical, meaning that there is a deliberate contact between the robot and the user, or contactless when common objectives are reached by other means of information exchange, such as intention prediction or direct communication by means of gesture/voice commands [START_REF] Ajoudani | Progress and prospects of the human-robot collaboration[END_REF].

Current norms and safety requirements

Ensuring safety for both the robot and the operator is a fundamental requirement in the design process of the collaborative workspace. In this context, the introduction of specific safety standards provides unified requirements and guidelines for the design and control of collaborative robots, as well as their integration into collaborative environments. In addition to general requirements that apply to industrial robots (ISO 10218 [ISO 10218-1:2011, 2011, ISO 10218-2:2011, 2011]), cobots need to satisfy safety countermeasures that are given by the technical specification ISO/TS 15066 [ISO/TS 15066:2016, 2016], which is dedicated to collaborative operations. The key idea behind the technical specification is to give guidelines to ensure that any contact between the robot and the operator will not result in any pain or injury. The information given in ISO 10218 and ISO/TS 15066 introduce four modes of collaboration that ensure security for the operator: Safety-rated monitored stop (SMS) is the basic form for collaboration, as both the robot and the operator can share a common workspace but not in the same time. This is ensured by sensing the presence of the operator inside the workspace and pausing the current task. After the operator leaves the workspace, the robot can resume the previous operation. This mode is depicted in Figure 1.3a.

Hand-guiding (HG) allows the operator to freely move the robot around the workspace. During the movement, the robot controller compensates for the gravity and limits the speed for safety. This mode is triggered by activating a specific device placed on the end-effector. HG allows to directly teach the desired position without the need of an Speed and separation monitoring (SSM) ensures operator security by adapting the speed of the robot accordingly to the distance separating the robot and the human. As shown in Figure 1.3c, when the operator is not present or is detected to be inside the green zone, the robot can work at full speed. The speed is reduced inside the yellow zone, and if the operator enters the red zone, the robot is stopped. The separation can be dynamically computed accordingly to the relative speed between the robot and the operator [START_REF] Lasota | A survey of methods for safe human-robot interaction[END_REF].

Power and force limiting (PFL) allows a physical interaction between the robot and the operator by limiting the robot force and motor power. This mode is depicted in Figure 1.3d. It requires dedicated sensors and control architectures to deal with physical human-robot interactions (pHRI), which is the key concept of this thesis and will be further discussed in the following sections.

Even though these modes can be implemented with conventional, industrial robots by means of additional sensors and enhanced controllers, cobots implement them by default. In fact, cobots are either equipped with joint-torque sensors, or estimate torques from motor currents, which allows to compute external forces that are applied to the robot. Dedicated controllers precisely monitor robot states and interactions, as well as environmental changes in order to react accordingly with respect to safety requirements.

Physical human-robot interactions

In the case of pHRI, the robot, the human operator and the environment form a coupled dynamic system in order to accomplish the desired task. For this to work, each collaborative partner has to be able to perceive and react to the contributions of the others. The most straightforward way to coordinate all partners lies in intuitive communication. One of the current challenges in pHRI consists in building communication protocols allowing the robot to predict operator intentions in order to properly react to them. As discussed in [START_REF] Ajoudani | Progress and prospects of the human-robot collaboration[END_REF], current research topics aim to establish such communication interfaces based on different technologies. The use of body gestures and vocal commands are common ways to communicate with collaborative systems. These communication protocols, well suited for high-level commands appear however to be less suited for complex tasks. Another common approach consist in implementing on the robotic system sensors that are able to measure the force input of the user. Human activity can typically be assessed by measuring the forces applied by the operator on the robotic system using mechanical force-torque sensors. The major drawback of such a method consists in the fact that the sensor not only measures the desired human-robot interaction, but also undesired components such as gravity, friction forces and unmodeled interaction forces, resulting from some external contact with an environment. Alternatively, measurable bio-signals such as electroencephalography (EEG) or electromyography (EMG) are investigated as they give a direct insight on the operator intentions [START_REF] Li | A review of the key technologies for semg-based human-robot interaction systems[END_REF].

Another important safety challenge that arises in pHRI consists in ensuring stability regardless of the interaction between the robot and the human, while achieving the desired performance. Indeed, when the robot and the user interact, the stability of the interaction highly depends on the coupled dynamics [START_REF] Buerger | Complementary stability and loop shaping for improved human-robot interaction[END_REF]. Even if both systems are stable separately, the coupling can lead to instability and the deterioration of performance. One way to achieve stability is to regulate the robot dynamic behavior by means of impedance control (IC), as defined in [Hogan, 1985]. IC allows to make the robot behave as a virtual mass-spring-damper system under the influence of external forces, hence controlling its dynamic behavior. In the classical form of IC, the virtual system parameters are considered to be constant and cannot be chosen arbitrary in order to ensure stability [START_REF] Kronander | Stability considerations for variable impedance control[END_REF]. The concept of passivity, introduced in [Colgate and Hogan, 1988] for IC, offers a solution to this issue, even though passivity guarantees can lead to very conservative constraints and thus impact the performance.

Contributions

The goal of this thesis is to propose some general tools to address current issues of collaborative robotics and implement these tools in two scenarios of robot-assisted medical procedures.

The first challenge that is addressed in this thesis is related to the safety aspects during pHRI. Current safety norms give some limitations for cobotic systems when operating close to humans. However, in most cases, it is important to be able to set some task specific constraints on robot motion such as position, velocity or acceleration constraints. In general robotic scenarios, constraints are implemented directly in the path planning process. This however is not possible in pHRI as the operator input is generally unknown. The first contribution of this thesis consists in a design method for controllers based on model predictive control (MPC) that behave as classical IC, while ensuring such practical constraints.

The second addressed challenge consists in guaranteeing robot passivity for pHRI. In fact, while IC is one of the most commonly used state-of-the-art strategies for robot interaction control, variable impedance control (VIC) is a more recent preoccupation. While designing IC with varying parameters allows increasing the system flexibility and dexterity, it is still a challenging issue as it may result in a loss of passivity of the control system. This has an important impact on the stability and therefore on the safety of the operator and the robot. The contribution consists in proposing methods to design passivity filters that guarantee passivity of the interaction. They aim at either checking whether a desired impedance profile is passive, or modifying it if required.

The third issue addressed in this manuscript is related to the perception of operator intentions for pHRI. The proposed contribution consists in establishing a model of the human arm dynamics from EMG and force measurements that allows to assess force 1.3. OUTLINE OF THE MANUSCRIPT inputs of the user and by such better adapt the robot control for improved interactions.

Finally, this thesis investigates two potential applications making use of the proposed tools in the medical context. As described previously, cobotic applications in the medical field are split into two main fields: assistance to medical interventions and rehabilitation. This thesis proposes a first application in the field of bimanual robot-assisted rehabilitation. This proof of concept offers a technical solution to a rather novel approach of functional rehabilitation. The second application is a proof of concept for autonomous needle insertion for interventional radiology that exploits the proposed tools.

Outline of the manuscript

This manuscript first introduces some general concepts of interaction management with IC and constrained robot control in Chapter 2. Then, in Chapter 3, the design method of a MPC that behaves as IC and deals with constraints is addressed, with some practical validation experiments. Chapter 4 deals with issues of passivity guarantees when variations of the IC dynamics are allowed. In Chapter 5, the problem of assessing operator intentions based on a human EMG to force model is explored and used to enhance the interaction. The proposed tools are then used in Chapter 6, with the two previously introduced medical tasks. Chapter 7 concludes the work by summarizing the contributions and providing further perspectives.

Scientific communication 1.4.1 Related publications

Parts of the presented work were published in the following conferences :

❼ [START_REF] Pesenti | Linear parameter-varying identification of the EMG-force relationship of the human arm[END_REF]] Pesenti, M., Alkhoury, Z., [START_REF] Pesenti | Linear parameter-varying identification of the EMG-force relationship of the human arm[END_REF]. Linear parameter-varying identification of the EMG-force relationship of the human arm. Chapter 2

Robotic background

In robotics, when no interaction with the environment is needed, motion control strategies are used. At the contrary, in the presence of interactions, direct force control strategies are preferred for fine force tracking. However, they require good models of the interaction and of the environment, and are not compatible with unpredictable interactions [START_REF] Siciliano | Springer Handbook of Robotics[END_REF]. For this reason, interaction control offers a compromise in order to deal with both force and motion during interactions. Impedance control (IC) [Hogan, 1985] is a widely used interaction control method, particularly efficient for human-robot interactions. It consists in imposing an impedance model for the relationship between the robot motion and external forces. Since its introduction by Hogan, IC has attracted researchers' attention, leading to improvements in handling robots flexibility [Ott, 2008] or variable impedance parameters [START_REF] Kronander | Stability considerations for variable impedance control[END_REF], to mention just a few.

In addition to the interaction management, the robot controller should handle several constraints, such as virtual fixtures [START_REF] Bowyer | Active constraints/virtual fixtures: A survey[END_REF], speed and energy limits [START_REF] Lasota | A survey of methods for safe human-robot interaction[END_REF], limited jerk [START_REF] Macfarlane | Jerk-bounded manipulator trajectory planning: design for real-time applications[END_REF] or actuators saturation. These limits allow ensuring safety [Haddadin, 2015] as well as good ergonomics of the human-robot interaction [START_REF] Kim | Anticipatory robot assistance for the prevention of human static joint overloading in human-robot collaboration[END_REF]. However, it is not straightforward to respect these constraints using classical IC schemes. Even if some methods such as anti-windup for input saturation [ Åström and Rundqwist, 1989] can be integrated into the controller architecture, their implementation complexity in the case of multiple constraints can lead to performance issues and some difficulties to address all kinds of constraints.

In this background chapter, after introducing some general concepts of robot control, the main focus is oriented towards IC and its ability to deal with interactions in robotics. Also, the control of robotic systems under constraints is introduced.

Rigid body dynamics and linearized model

The dynamical model of the robotic system, under the form of the equation of motion, links the forces applied on the robot with its acceleration [START_REF] Siciliano | Springer Handbook of Robotics[END_REF]. The model of an n-joint serial robotic manipulator operating in an m-dimensional task space can be written as

M (q)q ¨+ c(q ˙, q) + g(q) = τ c -J(q) T f ext (2.1)
where q ¨, q ˙, q ∈ R n are joint accelerations, velocities and positions, respectively, τ c ∈ R n are the commanded joint torques, f ext ∈ R m is the end-effector wrench, J(q) ∈ R m×n is the robot Jacobian matrix and M (q) ∈ R n×n , c(q ˙, q) ∈ R n and g(q) ∈ R n are the inertia matrix and the Coriolis and gravity terms, respectively. In the following development it is considered that m ≤ n. A classical solution [START_REF] Siciliano | Springer Handbook of Robotics[END_REF] to linearize equation (2.1) is to apply the control law τ c = M (q)ν + c(q ˙, q) + g(q) + J(q) T f ext with ν the new control input, resulting in the double integrator model in joint space

q ¨= ν (2.2)
The robot dynamics can then be expressed in task-space by derivating the forward differential kinematic model, which leads to p ¨= J(q)q ¨+ J ˙(q, q ˙)q ˙(2.3) where p ¨, p ˙, p ∈ R m denote the robot end-effector Cartesian acceleration, velocity and pose, respectively, and J ˙(q, q ˙) ∈ R m×n is the time derivative of the Jacobian matrix. Using (2.2) and (2.3) one can express a new input u such that ν = J(q) + u -J(q) + J ˙(q, q ˙)q ˙(2.4) with J(q) + ∈ R n×m the pseudo-inverse of J(q), resulting in a double integrator model in task-space

p ¨= u (2.5)
The two linearization steps are depicted in Figure 2.1.

Figure 2.1: Robot dynamic model linearization.

Robot compliance by control

In traditional robotics, the environment can be considered as known, which is not the case when robots share the same space as humans and when pHRI is desired. Service robots that are designed to assist humans need to be able to perform human daily-living tasks. Inspired by human arm design, service robots often feature anthropomorphic arms. Safety requirements and practical implementation issues lead to a lightweight design of these manipulators with an almost 1:1 load to weight ratio. In fact, stiff but lightweight robot links allow to preserve the accuracy of the robot while reducing the robot inertia. Additionally, as cobots need to work alongside humans, one of the possibilities to reduce injuries in case of undesired contact consists in ensuring compliance of the robot. Structure compliance can be achieved using compliant elements that can absorb the energy generated at impact. This, however, can lead to undesired vibrations and steady state displacements during robot motion and thus the degradation of performance. For this reason, much research effort has been made to be able to ensure compliance by means of control for better interaction management. Among these methods, IC [Hogan, 1985] is the most popular control strategy to address this issue.

Impedance control

2.2.1.1 Task-space and joint-space IC

In contrast to direct force and motion control methods, the goal of IC is not to regulate forces and positions, but the overall dynamic behavior of the system given by the relationship between force and position. This virtual mechanical relationship imposes the dynamical behavior of a mass-spring-damper system to the controlled robot and can be applied either at joint-space or task-space level. The behavior of the IC controlled robot is depicted in both joint-space and task-space in Figure 2.2. In the case of joint-space IC (Figure 2.2a), each link of the robot is attracted by the virtual mass-spring-damper system to its reference, whereas in task-space IC (Figure 2.2b) only the robot end-effector is attracted. In the following development, it is considered that the robot dynamics have been linearized, as explained in the previous section, resulting in a double integrator model of the robot both in joint-space and task-space. For this reason, even though the following considerations are done in task-space, the exact same strategy can be employed to control the robot using joint-space IC.

In order to dynamically link the system positions, velocities and accelerations with the external forces, the input u of system (2.5) is computed in order to obtain

M v (p ¨r -p ¨) + D v (p ˙r -p ˙) + K v (p r -p) = f ext (2.6)
The system subject to the external forces f ext is characterized by its impedance with an apparent virtual mass M v ∈ R m×m , a desired damping D v ∈ R m×m and a desired stiffness K v ∈ R m×m , tracking a reference motion p r . Note that M v , D v and K v are symmetric positive definite matrices. The desired behavior (2.6) can be obtained by imposing

u = p ¨r + M -1 v (D v (p ˙r -p ˙) + K v (p r -p) -f ext ) (2.7)
to the linearized system (2.5).

Null-space IC

The robot is redundant when there are more joints than required to perform the task (m < n). In this case, it is possible to use the extra degrees of freedom (DoFs) to perform a secondary task [Chiaverini et al., 2016, Liégeois, 1977]. In the case of IC, it is possible to combine the task-space IC with joint-space IC as secondary task to control at the same time the interaction dynamics of the robot end-effector and of its redundancy. This can be achieved by adding an additional term into the control law (2.4), such that ν = J(q) + (u -J ˙(q, q ˙)q ˙) + (I n -J(q) + J(q))µ

where the term (I n -J(q) + J(q)) ∈ R n×n is a projector on the Jacobian kernel, I n ∈ R n×n is the identity matrix and µ ∈ R n is a joint-space acceleration, that is computed as the joint-space version of equation (2.7) and is used to obtain a compliant behavior of the robot redundancy.

IC as state-feedback control

IC can also be represented using a state-space representation. In the following, we will consider external forces with slow variations so that f ˙ext ≈ 0. This simplification corresponds to many practical cases but could seem restrictive in the case of unmodeled contact. However, if a better model of the external force is available, it can be included at this stage into the state-space model. With this assumption, the dynamic model of equation (2.5) can be written as

x ˙= A c x + B c u (2.8) with x =    p ṗ f ext    ∈ R 3m , A c =    0 0 0 I m 0 0 0 0 0    ∈ R 3m×3m , B c =    I m 0 0    ∈ R 3m×m
The external force has been integrated into the state variable in order to better deal with this disturbance, inspired by the Internal Model Principle [START_REF] Francis | The internal model principle of control theory[END_REF]. Note however that the objective is not to cancel the disturbance, but rather to have an appropriate response with respect to it.

IC can now be expressed using a state-space formalism. To do so, note that (2.7) can be written as

u = u r + K(x r -x) (2.9) with u r = p ¨r ∈ R m , x r = p ˙T r p T r 0 T ∈ R 3m , K = M -1 v D v K v I m ∈ R m×3m
(2.10)

Vector r is some reference vector that can be chosen according to the desired task and produces a reference acceleration p ¨r, velocity p ˙r and position p r . By substituting (2.9) into (2.8), the state-space model of the impedance controlled system can be written as

x ˙= (A c -B c K)x + B c (u r + Kx r ) (2.11)
IC can therefore be seen as a state-feedback tracking problem, as represented in Figure 2.3. Since the controller is implemented numerically, it is necessary to find a discrete model of the system. The discretization of equation (2.8) with zero order hold yields

x k+1 = Ax k + Bu k (2.12)
where k represents the current step, x k is the discrete state at step k and A, B and C are the discrete state-space matrices. The sampling period T s is selected short enough to accurately approximate the continuous controller (2.9), which gives

u k = u r k + K(x r k -x k ) (2.13)
with x r k and u r k the discrete forms of x r and u r .

For some applications, when a specific interaction force is required, force tracking ability can be added to the classical IC scheme. In this case, the desired dynamic model becomes

M v (p ¨r -p ¨) + D v (p ˙r -p ˙) + K v (p r -p) = f r -f ext (2.14)
that can also be expressed as a state-space feedback (2.9) with

u r = p ¨r ∈ R m , x r = p ˙T r p T r f T r T ∈ R 3m , K = M -1 v D v K v -I m ∈ R m×3m
(2.15)

For joint-space, the same strategy can be used to express joint-space IC as a state-feedback controller.

The choice of appropriate impedance parameters leads to a suitable compliant behavior of the system during interaction. The higher the stiffness parameter is chosen, the better the controlled robot rejects disturbances such as modeling errors and unmodeled friction forces, resulting in a better tracking of the desired trajectory. However, high stiffness values lead to important interaction forces during contact. This issue results in the necessity for IC to find a good compromise between tracking performance in free motion and interaction dynamics. Also, in many robotic solutions and especially those designed for industrial purposes, only motion commands can be tracked by the internal controllers, which limits the possibility of implementing classical IC that requires torque commands. Admittance control (AC) is then a classical solution to overcome this issue.

Admittance Control

As highlighted in the previous section, classical IC implementations are not performant in rejecting disturbances and require implementing joint-level torque control. AC is a solution that deals with these problems by separating motion control from IC. In this control strategy, the motion controller is designed to ensure a good performance in tracking of the desired reference, while rejecting disturbances. The motion tracked by the motion controller is given by the impedance model

M v (p ¨r -p ¨c) + D v (p ˙r -p ˙c) + K v (p r -p c ) = f ext (2.16)
where p ¨c, p ˙c, p c ∈ R m define the so-called compliant frame [Villani and De Schutter, 2016].

The compliant frame is thus computed by the model (2.16) from the values of the reference frame, given by p ¨r, p ˙r, p r and the measured external wrench applied on the end-effector f ext , such that

p ¨c = p ¨r + M -1 v (D v (p ˙r -p ˙c) + K v (p r -p c ) -f ext ) (2.17)
The principle of computing the compliant frame from the impedance model is depicted in Figure 2.4. Hence, the compliant frame, and by such the resulting robot dynamic behavior, is only modified when interaction occurs. In free motion, no contact forces are applied, f ext = 0, thus the compliant frame fully coincides with the reference frame and the trajectory is well tracked by the robot. The principle of AC is depicted in Figure 2.5.

Figure 2.5: Admittance control in task-space.

Even though the described strategy is expressed in task-space, it can also be formulated in joint-space. The reference and compliant frames are then expressed in terms of joint accelerations, velocities and positions and the external torques applied on the joints are used to compute the compliant frame.

Robot control using constrained optimization

A common way to deal with constraints consists in the resolution of constrained optimization problems. Constrained optimization is particularly advantageous to deal with multi-variable systems with multiple constraints. This strategy deals with the problem of finding an optimal control law that satisfies at the same time some objective criterion and a certain number of constraints.

Constrained QP problems

Constrained optimization problems consist in minimizing some objective function with respect to an optimization variable z ∈ R nz with n z the dimension of z, subject to a set of constraints. The optimal solution is computed as the value of z that leads to the smallest value of the objective function while satisfying the constraints. In robotics, the cost functions are often designed such that z is a control variable that allows the tracking of some desired objective given by Ez = d, with E ∈ R nz×nz and d ∈ R nz . In cost functions, the square Euclidean norm of the desired objective is often used, such that the objective function is equal to ||Ez -d|| 2 2 . Such cost functions can be written in quadratic form and solved using a quadratic programming method that seeks the optimal solution of the problem

min z 1 2 z T Hz + f T z (2.18) s.t.    Ω z z ≤ ω Π z z = π
with H ∈ R nz×nz and f ∈ R nz the Hessian matrix and the gradient vector of the objective function, respectively. The matrices Ω z ∈ R n i ×nz and Π z ∈ R ne×nz and the vectors ω ∈ R n i and π ∈ R ne are associated with desired constraints, with n i and n e the number of inequality and equality constraints, respectively. Constrained QP problems require the use of iterative methods to find the optimal solution that satisfies the given constraints.

QP solvers require high performance algorithms capable of finding optimal solutions in a short time, especially in contexts that require on-line computation. It is important to keep in mind that constrained problems can lead to situations in which there is no solution and which need to be dealt with.

Problem resolution methods

In the literature two methods are known to offer the best performance in dealing with constrained QP problems: the Interior Point method [Wright, 1997] and the Active Set method [Fletcher, 2000].

The idea behind the Active Set method consists in determining a set of the constraints that are active and use them as equality constraints to solve the resulting QP problem. In fact, when a constraint is reached, a subset of the optimization problem lies on this constraint. The challenge of Active Set algorithms lies in the ability to detect active constraints in a minimum of iterations.

In the Interior Point method, the constrained QP is transformed into an optimization problem in which the inequality constraints are included into the objective function. Problem (3.11) is rewritten as

min z 1 2 z T Hz + f T z -λ i log(ω i -ω T z,i z) (2.19) s.t. Π z z = π
where λ ≥ 0 is the weight associated with the constraints, ω T z,i is the i-th row of Ω z and ω i is the i-th element of ω. The additional term in the cost function is referred to as the barrier function and prevents the solution from going outside the feasible region. In fact, when approaching the limit, the weight associated with the barrier function increases and becomes infinite at the boundary. Convergence to the original QP problem solution is ensured by modifying the weight of the barrier function in such a way that the iterations follow a path from a strictly feasible point towards the solution.

Constraints softening

In practical implementations of constrained QP solving strategies, the imposed constraints can lead to infeasible problems. Infeasibility is a critical issue especially when dealing with safety critical online controllers. Several strategies can be applied to handle this issue, ranging from simple solutions such as re-implementing the previously computed input to more complex ones that search for the least important constraint that can be relaxed in order to regain feasibility [Rawlings andMuske, 1993, Tyler and[START_REF] Tyler | [END_REF]. A more systematic way consists in imposing soft constraints that allow some additional flexibility in critical situations rather than hard constraints that must not be crossed. This can be achieved by introducing slack variables into the QP problem [Scokaert andRawlings, 1999, Maciejowski, 2002]. These variables allow constraints violations, but are heavily penalized, forcing the optimizer to keep them as low as possible. One way to include a slack variable ǫ into the constrained QP problem (3.11) results in min z,ǫ

1 2 z T Hz + f T z + ρǫ (2.20) s.t.      Ω z z ≤ ω + ǫ ǫ ≥ 0 Π z z = π
where ǫ ∈ R n i is a vector of ǫ and ρ is a gain penalizing any increase of the slack variable. If ρ is chosen large enough, constraint violation will only occur if the problem becomes infeasible. Otherwise, the problem will result in the same feasible solutions as the original problem.

Constrained optimization in robotics

With the fast development of new and more powerful tools for solving constrained optimization problems, much research effort has been done to use this technology for robot control to handle physical and/or computational constraints. Constrained optimization can be applied in various situations such as in [START_REF] Kapoor | Constrained control for surgical assistant robots[END_REF], where the authors use constrained optimization to define virtual fixtures applied during surgical procedures, or in [Kanoun et al., 2011, Bouyarmane andKheddar, 2018] where constrained optimization is used to deal with redundant humanoid robots. Alternatively, [START_REF] Joseph | Towards X-Ray Medical Imaging with Robots in the Open: Safety Without Compromising Performances[END_REF], Joseph et al., 2018a] impose energetic constraints to ensure safety in a workspace shared by humans and robots. Also, in [START_REF] Rubrecht | Motion safety and constraints compatibility for multibody robots[END_REF], constraints are used to take into account the robot actuation limits. The combination of constrained optimization and IC was proposed in [START_REF] Zhao | Quadratic programming and impedance control for transfemoral prosthesis[END_REF] to ensure stability and robustness in the control of a transfemoral prosthesis and in [Bouyarmane et al., 2018] in a context of task-space force control. Also, in [START_REF] Hoffman | Multi-priority cartesian impedance control based on quadratic programming optimization[END_REF], the combination of IC and QP is used to deal with robot redundancy and in [START_REF] Lutscher | Hierarchical force and positioning task specification for indirect force controlled robots[END_REF] it is used to ensure stability during interaction.

In the case of constrained optimization under the form of a QP, the designed controller efficiency highly depends on the solver performance and its ability to solve the QP problem online and produce the optimal closed-loop response. Online QP resolution has proven particularly efficient when handling constraint control problems, especially for position controlled systems for which control update rates are not demanding. Modern QP solvers, such as [START_REF] Ferreau | qpOASES: a parametric active-set algorithm for quadratic programming[END_REF], are performant enough to handle online optimization at relatively high frequencies (200Hz -1kHz), as required for low level control strategies, e.g. joint-torque control. However, as the reactivity of the controlled system highly depends on the control periodicity, a system controlled with small time steps has a shorter reaction time. This can lead to important braking torques that are required to prevent the system from stepping outside the imposed limits, which can require torques that are not achievable by the actuators. One way to overcome this issue was proposed in [START_REF] Meguenani | Control of robots sharing their workspace with humans: An energetic approach to safety[END_REF], Meguenani et al., 2017] by integrating the braking capacity of the robot into the optimization problem. Even though the proposed method has proven its efficiency in limiting the system energy, it also shows a major drawback as it requires a rather high control periodicity to give the system sufficient time for braking. Thus, this issue has an important impact on the performance of the system, especially when controlling the robot at joint-torque level. Furthermore, controlling robotic systems with low periodicity using constrained optimization methods can lead to another issue that is caused by the discretization of the constraints. In fact, the computed solution satisfies constraints only for the discretized system, what can possibly lead to the violation of the constraints in between the samples. To overcome this issue, the authors of [START_REF] Lengagne | Guaranteed computation of constraints for safe path planning[END_REF], Lengagne et al., 2011] propose a method using interval analysis to guarantee that constraints are always satisfied. An alternative solution to deal with the issue of the system reactivity is to increase the ability of the controller to anticipate its future actions and the activation of constraints.

A common way to deal with constrained multi-variable control problems consists in using a predictive approach, such as Model Predictive Control (MPC). The advantage of MPC compared to the general constrained optimization methods lies in its ability to predict how and when the controlled system meets the imposed constraints. By doing so, MPC is capable of reacting accordingly, resulting in smoother transitions from free to constrained motions. MPC has been only recently used for interaction control and is further explored in the following chapter.

Chapter 3

Constrained impedance control

The use of robotic devices in environments occupied by humans has recently become a standard, in contrary to previous setups where humans and robot could only evolve in split environments. An important factor in the integration of robots in daily living activities lies in the interaction between the human and the robot, especially regarding their physical integrity. Handling even unexpected interactions can easily be achieved using classical control schemes, such as IC [Hogan, 1985]. Also, in many cases it is important to be able to constrain the robot to some task specific limits that either ensure security or allow a better execution of the desired task, which can be achieved by dedicated control strategies. The challenge, however, arises when the application requires at the same time the system to impose an impedance behavior and satisfy some given constraints. This chapter describes why MPC is particularly suited to deal with the necessity of designing controllers that are capable of ensuring robot compliance during interactions while satisfying constraints resulting from safety requirements or practical limits specific to the designed application. It is described how to design an appropriate MPC formulation to obtain a prescribed compliant behavior, while satisfying the set of imposed constraints. Finally, some simulations and experimental results are shown to validate the proposed controller.

MPC and interaction management

With the fast technological development of computer systems allowing the implementation of computationally intensive control strategies, the interest in efficient design strategies for controllers capable of dealing with complex and constrained systems has grown.

Modern MPC is the result of major advances in the field of mathematics and computer science which have greatly improved the speed and reliability of the algorithms used for optimization. First proposed for process control in the petrochemical industry [START_REF] Richalet | Algorithmic control of industrial processes[END_REF], MPC has successfully replaced classical methods in many applications [START_REF] Qin | [END_REF]Badgwell, 2003, Raković andLevine, 2018].

Recently, MPC has attracted researchers' attention for dealing with interactions in robotics. In [START_REF] Killpack | Model predictive control for fast reaching in clutter[END_REF], MPC is used as a high-level controller for path set-point generation, taking into account contact forces due to interactions with the environment. The authors of [START_REF] Matschek | Force feedback and path following using predictive control: Concept and application to a lightweight robot[END_REF] investigate a direct force control strategy based on predictive path following using MPC. In [START_REF] Wahrburg | MPC-based admittance control for robotic manipulators[END_REF], MPC is combined with admittance control for dealing with stability issues while interacting with a stiff environment. Even more recently, non-linear MPC is formulated by adding admittance dynamics into a path following problem [START_REF] Kazim | Combined predictive path following and admittance control[END_REF]. However, in all these applications, the choice of the cost functions is not directly related to the desired compliant behavior.

Model Predictive Control

MPC is a particular case of optimal control. The main difference between standard optimal control strategies like linear quadratic regulation (LQR) [ Åström and Wittenmark, 2013] and MPC is that, instead of computing the control once over some time horizon, MPC optimizes over a receding horizon, which is shifting over time and only applying the first value of the optimized control. This allows to adapt the control for constraints and makes the system closed-looped.

MPC principles

In classical MPC formalism [START_REF] Lewis | Optimal control[END_REF], the general form of an MPC problem is given by min

x,u J = min x,u φ(x H ) + H-1 k=0 l k (x k , u k ) (3.1) s.t. x k+1 = f (x k , u k ) x 0 = measured
where x k ∈ R nx and u k ∈ R nu are respectively the system state and control at instant k, H is the control horizon, x 0 ∈ R nx and x H ∈ R nx are initial and final state, f (.) represents the plant dynamics, φ and l k are respectively the terminal and running costs, and

x =       x 1 x 2 . . . x H       ∈ R nxH u =       u 0 u 1 . . . u H-1       ∈ R nuH (3.2)
The MPC problem (3.1) can also include additional equality and inequality constraints, which are not specified here for better readability.

As described in [START_REF] Lewis | Optimal control[END_REF], the optimization criterion of this type is very flexible and can be used to formulate different control problems such as ❼ Minimum time problem, for φ = 0 and l k = 1

❼ Minimum fuel problem, for φ = 0 and l k = |u k | ❼ Minimum energy problem, for φ = 1 2 x T H P x H and l k = x T k Qx k + u T k Ru k
The algorithm employed to resolve problem (3.1) using MPC consists in computing the new control input for the plant such that at each time step, the controller:

1. Takes the measurements of the system output or state 2. Computes the sequence of inputs over the finite horizon H ❼ Predicting the future states with an internal model ❼ Optimizing the cost function of future states and control inputs 3. Implements the first optimal input u 0 and discards the rest of the sequence.

MPC stabilization and QP

In order to illustrate the use of MPC, this section introduces a pedagogical example dealing with the stabilization problem of an linear time invariant (LTI) system using MPC. In the case where the system is a LTI plant model, f (x k , u k ) = Ax k + Bu k with A ∈ R nx×nx and B ∈ R nx×nu . The following stabilization problem of an LTI using MPC with a quadratic cost is considered min

x,u J = min x,u 1 2 x T H P x H + 1 2 H-1 k=0 x T k Qx k + u T k Ru k (3.3) s.t. x k+1 = Ax k + Bu k x 0 = measured where Q ∈ R nx×nx , R ∈ R nu×nu , P ∈ R nx×nx
are the weights of the particular terms such that Q = Q T ≥ 0, P = P T ≥ 0 and R = R T > 0. At this point, it is possible to express problem (3.3) in matrix form min

x,u

J = min x,u 1 2 x T Qx + 1 2 u T Ru + 1 2 x T 0 Qx 0 (3.4) s.t. x = Āx + Bu + A 0 x 0 where Q = diag(Q, ..., Q, P ) ∈ R nxH×nxH , R = diag(R, ..., R) ∈ R nuH×nuH , B = diag(B, ..., B) ∈ R nxH×nuH and Ā = 0 0 diag(A, ..., A) 0 ∈ R nxH×nxH A 0 = A T 0 . . . 0 T ∈ R nxH×nx
with 0 a null matrix of appropriate size. As 1 2 x T 0 Qx 0 is constant, it can be removed from the cost function. The previous problem can thus be simplified in a compact form

min w 1 2 w T Q 0 0 R w (3.5) s.t. ( Ā -I) B w + A 0 x 0 = 0 with w = x T u T T ∈ R (nx+nu)H . Problem (3.5
) has the form of a QP problem with constraints and can therefore be solved numerically. However, the presence of constraints does not permits to have a straightforward analytical solution. Therefore, further simplification can be performed in order to eliminate the constraint by transforming the optimization problem (3.5). In order to do so, the first step consists in removing the state vector x from the equations. This can be done as the plant model is available. In fact, from the discrete state space equations one can generalize that for k = 1, ..., H

x k =A k x 0 + A k-1 B A k-2 B . . . AB B          u 0 u 1 . . . u k-2 u k-1          leading to x = Cu + Ax 0 (3.6) with C =       B 0 ... 0 AB B ... 0 . . . . . . . . . . . . A H-1 B A H-2 B . . . B       ∈ R nxH×nuH A =       A A 2 . . . A H       x 0 ∈ R nxH×nx
J = min u 1 2 (Cu + Ax 0 ) T Q(Cu + Ax 0 ) + 1 2 u T Ru + 1 2 x T 0 Qx 0 = min u 1 2 u T (C T QC + R)u + x T 0 A T QCu + 1 2 x T 0 (A T QA + Q)x 0 (3.7)
At this point it is worth noticing that because of the simplification given by equation (3.6), the state vector x can be removed from the optimization problem defined by the cost function (3.7). For simplification, as the term 1 2 x T 0 (A T QA + Q)x 0 is constant, it can be removed from the cost function (3.7), which can be expressed as

min u 1 2 u T Hu + x T 0 F T u (3.8) with H = C T QC + R ∈ R nuH×nuH and F T = A T QC ∈ R nx×nuH .
It is worth noticing that because of the properties of the weight matrices Q, P and R, the problem is convex, as H = H T > 0. Equation (3.8) is a quadratic programming problem under the general form of problem (2.18) without constraints and therefore has an analytical solution. The optimal control vector u is the solution of equation ∇ u J (u, x 0 ) = Hu + F x 0 = 0, where ∇ u denotes the gradient with respect to vector u, leading to u = -H -1 F x 0 .

For MPC, only the first term of the optimal control u will be taken as the new control input

u 0 = -I nu 0 ... 0 H -1 F x 0
At this point it is possible to add constraints to the QP. However the resulting solution has to be computed using a numerical solver. Constraints on u k and x k (or y k = Cx k ) can, for instance, be added as follows

u k ≤ u max u k ≥ u min ⇔ I nuH -I nuH u ≤ u max -u min (3.9) x k ≤ x max x k ≥ x min ⇔ I nxH -I nxH x ≤ x max -x min (3.10)
where u min , u max ∈ R nuH and x min , x max ∈ R nxH are vectors of the upper and lower bounds of the vector of inputs and states respectively. The constraints given by equations (3.9) and (3.10) can be expressed in terms of the optimization vector u using equation (3.6). Finally, the MPC problem can be expressed as a constrained QP problem defined by

min u 1 2 u T Hu + x T 0 F T u (3.11) s.t. Ω u u ≤ ω + Ω x x 0 with Ω u ∈ R nc×nuH , ω ∈ R nc and Ω x ∈ R nc×nx formulated based on n c constraints.

Transforming a state-space controller into MPC

MPC is often applied to systems that already have a standard LTI controller, to take advantage from the ability of MPC to predict the system evolution and deal with constraints.

For this reason, it becomes very interesting to transform the existing LTI controller into MPC. The MPC is tuned in such a way, that in the case where the system operates in between the predefined boundaries, it behaves as if controlled by the existing controller. Several controller matching methods have been proposed in the literature, which are detailed in the following sections. It is assumed here that the plant and the stabilizing controller are both expressed in state space and are both linear and time invariant.

The following methods allow to tune a MPC in such a way that, in the case where no constraints are active, it will behave as the desired LTI controller.

Matching based on QP matrices

The first method is proposed in [Di Cairano andBemporad, 2009, Di Cairano andBemporad, 2010] and uses the QP version of the cost function, such as defined in problem (3.11). If no constraints are active, the solution of the optimization problem used for MPC should be equivalent to the behavior of the desired state-feedback controller, that is u k = -Kx k , with K ∈ R nu×nx . On the other hand, as shown in Section 3.1.1.2, the solution of the unconstrained stabilization problem with MPC is given by

u 0 = -[I nu 0 ... 0]H -1 F x 0
The problem of matching the input computed by the unconstrained MPC with that given by the desired state-feedback controller results in finding the weight matrices Q, R and P such that

-[I nu 0 ... 0]H -1 F x 0 = -Kx 0 (3.12)
The complexity of this problem results from the fact that the matrix [I nu 0 ... 0] is not invertible, whereas H -1 is the inverse of a function depending on Q, R and P . The method described in [Di Cairano andBemporad, 2009, Di Cairano andBemporad, 2010] consists in reformulating (3.12) as an optimization problem in which the objective is given by the equality of equation (3.12), such that min

Q,R ||HK + F || s.t. P ≥ 0, Q ≥ 0, R ≥ σI (σ ∈ R, σ > 0)
which is a convex problem subject to linear matrix inequality (LMI) constraints with

K =       K K(A -BK) . . . K(A -BK) H-1      
When the value of the cost function is 0, the MPC controller without active constraints behaves as the desired state-feedback controller. If, however, the solution is not 0 and constraints are not active, MPC behaves as the best approximation of the targeted controller.

In addition, in order to guarantee global stability and maintain the equivalence with the desired controller, an additional LMI constraint, resulting from stability considerations based on a Lyapunov function, has to be added such that

(A -BK) T P (A -BK) -P + Q + K T RK ≤ 0
Even though using this method it is possible to prove asymptotic stability, it requires solving LMI what makes the problem computationally complex and does not guarantee to have the exact matching between MPC and the desired state-feedback controller.

Zero-value cost function

As the desired behavior consists in imposing u k = -Kx k , another way to formulate the corresponding cost function for MPC consists in making the matching of the controller as the primary objective of the optimization problem. To do so, a zero-value infinite horizon cost function [START_REF] Hartley | Initial tuning of predictive controllers by reverse engineering[END_REF] can be designed such that

J = ∞ k=0 ||u k + Kx k || R
with ||y|| R = y T Ry. This cost function is optimal for the desired control law for any R = R T > 0 with R ∈ R nu×nu [Molinari, 1973]. Furthermore, in [START_REF] Kreindler | Optimality of linear control systems[END_REF] it is pointed out, that by allowing cross-terms between input and state in the quadratic infinite horizon cost function, it is possible to reproduce any multi-variable state-feedback by making this equivalence as primary control objective. Doing so results

MPC AND INTERACTION MANAGEMENT

in the following cost function

J = ∞ k=0 x k u k T K T RK K T R RK R x k u k (3.13) = ∞ k=0 x k u k T Q S S T R x k u k = ∞ k=0 x T k Qx k + u T k Ru k + x T k Su k + u T k S T x k (3.14)
where

Q = K T RK ∈ R nx×nx and S = K T R ∈ R nx×nu .
However, practical MPC implementations can only perform optimization over a finite, receding horizon. That is why the cost function has to be transformed in such a way that it uses the cost function (3.14) over the finite horizon of length H [START_REF] Hartley | Initial tuning of predictive controllers by reverse engineering[END_REF]. One classical way to do so consists in computing the cost function over the horizon H and adding a terminal cost that uses the solution P of the discrete-time algebraic Riccati equation (DARE) as the terminal cost weight, what makes the finite and infinite horizon cots functions equivalent [Maciejowski, 2007, Hartley and[START_REF] Hartley | [END_REF]. The cost function given by equation (3.14) is then equivalent to

J = x T H P x H + H-1 k=0 x T k Qx k + u T k Ru k + x T k Su k + u T k S T x k (3.15)
with the weight associated with the terminal cost being the solution of

P = A T P A -(A T P B + S)(B T P B + R) -1 (B T P A + S T ) + Q
An interesting point to notice is that in this particular case where

Q = K T RK ∈ R nx×nx and S = K T R ∈ R nx×nu ,
the terminal cost weight is the solution of

P = A T P A -(A T P B + K T R)(B T P B + R) -1 (B T P A + RK) + K T RK (3.16)
Furthermore, the optimal static state-feedback [ Åström and Wittenmark, 2013] can be defined such that

K = (B T P B + R) -1 (B T P A + RK) (3.17)
Substituting equation (3.17) into (3.16) gives

A T P A -P -(A T P B + K T R)(B T P B + R) -1 (B T P A + RK) + K T RK = 0
This shows that P = 0 is a solution and means that in the case of controller matching, no terminal cost is needed to the finite horizon cost function [START_REF] Hartley | Initial tuning of predictive controllers by reverse engineering[END_REF], i.e.

J = H-1 k=0 (u k + Kx k ) T R(u k + Kx k ) = H-1 k=0 x T k Qx k + u T k Ru k + x T k Su k + u T k S T x k
This feature is very advantageous for practical aspects, however it causes some issues when it comes to proving stability where most classical methods use x T k P x k as the Lyapunov function. In [START_REF] Hartley | Initial tuning of predictive controllers by reverse engineering[END_REF] it is pointed out that if theoretical stability guarantees are required, a terminal constraint can be used, as defined in [START_REF] Mayne | Constrained model predictive control: Stability and optimality[END_REF]. Furthermore, the equivalence does not depend on the size of the control horizon that can be chosen H = 1 resulting in the same type of cost function, ||Ez -d|| 2 2 , as given in Section 2.3. This shows that classical objective tracking can be seen as a one step ahead prediction and that methods based on such a cost function could easily be reformulated using MPC.

Model Predictive Impedance Control

In order to deal at the same time with unpredictable interactions as does IC, and different types of limits as does MPC, a novel controller, called Model Predictive Impedance Controller (MPIC) is proposed in this section.

When no constraints are active, MPIC should be equivalent to IC, i.e. to the static state-feedback controller (2.13), which, as a reminder, is given by

u k = u r k + K(x r k -x k )
As it has been shown in the previous chapter, equation (2.13) can be employed to express IC both in joint-space and task-space, which in the case of a redundant manipulator differ in dimensions. Thus, in order to keep the proposed method as general as possible, in the following development, the dimensions of the state and input vector will be generalized to n x for the state vector and n u for the input vector.

As shown in the previous section, the controller matching method based on QP matrices can make MPC behave as the desired state-feedback controller. However, this method does not guarantee the exact matching between the MPC and the desired state-feedback controller and its complexity makes it computationally expensive. In order to obtain the controller of equation (2.13), the cost function J of the MPC problem has to be designed such that its optimum is zero when the desired behavior is obtained. To do so, the method based on zero-value cost function over a control horizon H is used because of its practical simplicity and the ability to exactly reproduce the desired behavior.

MPIC formulation

The MPC cost function can be written as a zero-value cost function such that

J = H-1 k=0 ||(u k -K(x r k -x k ) -u r k )|| R (3.18)
with ||y|| R = y T Ry. The optimal solution is obtained for the desired control law (2.13) for any R = R T > 0 ∈ R nu×nu . Furthermore, J can be written as

J = H-1 k=0 u r k -u k x r k -x k T R S T S Q u r k -u k x r k -x k (3.19)
where

Q = K T RK ∈ R nx×nx and S = K T R ∈ R nx×nu
are the weights of the particular terms satisfying the conditions Q = Q T ≥ 0 and Q-SR -1 S T ≥ 0 guaranteeing convexity, and H is the control horizon. With the discrete state-space plant dynamics (2.12), a MPC problem can be written as the tracking of some references, with feedback x r k and u r k terms min

u k J = min u k H-1 k=0 u r k -u k x r k -x k T R S T S Q u r k -u k x r k -x k (3.20) s.t. x k+1 = Ax k + Bu k x 0 = measured current state
where x 0 is the initial state measured at each control cycle. In the following, the input rate of change ∆u k = u k -u k-1 is included into the MPC problem. By doing so, the variations of the control input can also be easily limited. In the case of IC, this strategy allows to bound the jerk of the system in a straightforward manner, which is particularly interesting for increasing the smoothness of the interaction. Consider the following notations

∆u =[∆u T 0 ... ∆u T H-1 ] T ∈ R nuH u =[u T -1 u T 0 ... u T H-1 ] T ∈ R nu(H+1) (3.21) x =[x T 0 ... x T H-1 ] T ∈ R nxH ρ =[u r 0 T ... u r H-1 T x r 0 T ... x r H-1 T ] T ∈ R (nu+nx)H
which define the vectors containing all the input rates of change, inputs, states and references over the control horizon H. Notation u -1 represents the last applied input. The minimization problem (3.20) can be transformed into standard QP form. To do so, it has first to be expressed in matrix form such that

min ∆u J = min ∆u       ∆u x u ρ       T Φ T R S T S Q Φ       ∆u x u ρ       (3.22) s.t. x = Ax 0 + Cu u = T ∆u + Iu -1 (3.23)
with u, ∆u, x and ρ defined by (3.21) and

Q =diag(Q, ... , Q) ∈ R nxH×nxH R =diag(R, ... , R) ∈ R nuH×nuH S =diag(S, ... , S) ∈ R nxH×nuH I =[I nu ... I nu ] T ∈ R nu(H+1)×nu
and block matrices H+1) is the H-steps state reachability matrix and A ∈ R nxH×nx the H-steps free evolution matrix [Di Cairano and Bemporad, 2010], such that

T =          0 0 . . . 0 I nu 0 . . . 0 . . . . . . . . . . . . I nu . . . . . . 0 I nu . . . . . . I nu          ∈ R nu(H+1)×nuH and Φ = -I nuH 0 -I nuH 0 I nuH 0 0 -I nxH 0 0 0 I nxH ∈ R (nu+nx)H×(3nu+2nx)H+nu Matrix C ∈ R nxH×nu(
C =          0 0 0 . . . 0 0 B 0 . . . 0 0 AB B . . . 0 . . . . . . . . . . . . . . . 0 A H-3 B A H-4 B . . . B          A =          I nx A A 2 . . . A H-2          Problem (3.22
) is a constrained QP. However, a simplification of the problem can be performed by including the constraints (3.23) into the cost function, so that

J =       ∆u x 0 u -1 ρ       T Ψ T Φ T R S T S Q ΦΨ       ∆u x 0 u -1 ρ       (3.24) with Ψ =       I nuH 0 0 0 CT A CI 0 T 0 I 0 0 0 0 I (nu+nx)H       ∈ R ((3nu+2nx)H+nu)×(nu(2H+1)+nx(H+1))
By defining matrices H+1) and

H ∈ R nuH×nuH , F ∈ R nuH×nu(2H+1)+nx(
L ∈ R nu(2H+1)+nx(H+1)×nu(2H+1)+nx(H+1) such that H F F T L = Ψ T Φ T R S T S Q ΦΨ (3.25) the optimization problem (3.22) becomes min ∆u J = min ∆u ∆u T H∆u + 2[x T 0 u T -1 ρ T ]F T ∆u (3.26)
where the term corresponding to L has been removed as it is constant. Because of the hypotheses on the weight matrices, it comes that H > 0, thus equation (3.26) is a convex QP. Note that, as the problem is convex, the solution corresponding to the zero cost function is the unique global solution. This means, that when unconstrained, the optimal input given by the MPIC controller will behave as the desired state-feedback. Without constraints, (3.26) has the following analytic solution

∆u = -H -1 F x T 0 u T -1 ρ T T (3.27)
In MPC, only the first input is applied at each iteration, i.e. u 0 = u -1 -I nu 0 ... 0 ∆u. Note that this applies only when there are no other constraints than the system model. However, in most practical cases, additional constraints have to be considered. In this case, no straightforward analytical solution exists and a numerical solver has to be used.

In order to do so, the constraints need to be written in the standard form

Ω ∆u ∆u ≤ ω + Ω x x T 0 u T -1 T (3.28)
where nu+nx) and ω ∈ R nc are formulated based on n c constraints.

Ω ∆u ∈ R nc×nuH , Ω x ∈ R nc×(

Controller matching

This section develops the optimal solution of the unconstrained MPC problem in order to prove the exact controller matching when no constraints are activated. Even though the matching of the two controllers could be directly concluded from the design of the zero-value cost function, an analytic proof can also be derived.

From equations (2.12) and (2.13) it can be shown that in the unconstrained case the input can be expressed as

u 0 = -Kx 0 + u r 0 + Kx r 0 u k = -K(A -BK) k x 0 + u r k + Kx r k + k-1 j=0 K(A -BK) k-1-j B(u r j + Kx r j ) ∀k ≥ 1 and ∆u 0 = -Kx 0 -u -1 + u r 0 + Kx r 0 ∆u 1 =(K -K(A -BK))x 0 -(I nu + KB)(u r 0 + Kx r 0 ) + u r 1 + Kx r 1 ∆u k =u k -u k-1 ∀k ≥ 2 =K((A -BK) k-1 -(A -BK) k )x 0 -(I nu + KB)(u r k-1 + Kx r k-1 ) + (u r k + Kx r k ) +K k-2 j=0 ((A -BK) k-2-j -(A -BK) k-1-j )B(u r j + Kx r j )
By defining K = diag(K, ... , K) ∈ R nuH×nxH , the input rate of change vector ∆u over the control horizon H becomes, in matrix form

∆u = V W M I nuH K x T 0 u T -1 ρ T T with V =          -K K -K(A -BK) K(A -BK) -K(A -BK) 2 . . . K((A -BK) H-3 -(A -BK) H-2 )          ∈ R nuH×nx , W =         -I nu 0 0 . . . 0         ∈ R nuH×nu and M =            I nu 0 0 . . . 0 -(I nu + KB) I nu 0 . . . 0 KB -(A -BK)B -(I nu + KB) I nu . . . 0 . . . . . . . . . . . . . . . K((A -BK) H-4 . . . . . . . . . I nu -(A -BK) H-3 )B            ∈ R nuH×nuH
In order to simplify this expression, it is possible to rewrite matrices V and W such that H+1) . The triangular matrix M is invertible and it can be shown by induction that M = (I nuH + I ˜T + KCT ) -1 which allows to define

V = -MKA W = -M(I ˜+ KC)I with I ˜= [I nuH 0] ∈ R nuH×nu(
κ = (I nuH + I ˜T + KCT ) -1 KA (I ˜+ KC)I -I nuH K ∈ R nuH×(nu+nx)(H+1) that leads to ∆u = -κ x T 0 u T -1 ρ T T
On the other side, by developing matrices H and F from (3.25) and using the definitions of S and Q from (3.19) in matrix form, i.e. S = K T R ∈ R nxH×nuH and

Q = K T RK ∈ R nxH×nxH , it comes that H =(I nuH + I ˜T + KCT ) T R(I nuH + I ˜T + KCT )
and

F =(I nH + I ˜T + KCT ) T R KA (I ˜+ KC)I -[I nH K]
From the previous statement, it can be written that Hκ = F , leading to κ = H -1 F as H is invertible. Finally, from (3.27) and the fact that only the first optimal input is used for MPC

u 0 = u -1 -I nu 0 ... 0 H -1 F x T 0 u T -1 ρ T T = u -1 -K I nu -[I nu 0 ... 0] -[K 0 ... 0] x T 0 u T -1 ρ T T = -Kx 0 + (u r 0 + Kx r 0 )
This corresponds to the desired behavior of IC in the absence of constraints, as expressed in equation (2.13).

MPIC and admittance control

In many robotic solutions and especially the ones designed for industrial applications, controlling the robot at joint-torque level is not possible and only position or velocity commands are implementable. As discussed in Chapter 2, a common way to impose compliant interaction dynamics consists in generating a compliant frame by means of admittance control and applying compliant frames to the robot position/velocity controller (Section 2.2.2). In the case of MPIC, the control input that is produced is equivalent to the desired robot acceleration. This feature can be used in order to easily adapt MPIC to robots that cannot be torque controlled efficiently. Thus, it is possible to generate position/velocity control inputs by successive integrations of the acceleration inputs computed usng MPIC. This specially applies in the case of pHRI, as robot motion is slow and it can be considered that the controller reference positions or velocities are perfectly tracked, which is a common assumption [Villani and De Schutter, 2016].

Validation

MPIC vs constrained QP optimization

In order to illustrate the prediction advantages of MPIC compared to classical constrained QP optimization, a simulation featuring both controllers is performed. The MPIC control horizon is evaluated for three values H = [3, 5, 10]. The constrained QP optimization controller is designed with a cost function given by ||u

0 -K(x r 0 -x 0 ) -u r 0 )|| R resulting in the online QP problem min u 0 u T 0 Ru 0 + 2 x T 0 u r 0 T x r 0 T    RK -R -RK    u 0 (3.29) s.t. Ω u u 0 ≤ ω + Ω x x 0
The considered system is a simple mass moving along a horizontal line, subject to a control force and an external force. The system dynamics are equivalent to a double integrator.

The desired impedance model is chosen with the following parameters:

M v = 20Ns 2 /m, K v = 1000N/m and D v = 89.
44Ns/m to have a damping ratio of 1. The system follows a smooth reference trajectory. The system is controlled using MPIC with different control horizons and the QP optimization problem (3.29) with a sampling time T s = 1ms. System velocity and position are constrained to ±0.5m/s and [0, 1]m, respectively. A constant external force, f ext = 1200N, is applied to the system between t = 1s and t = 1.18s. The simulation results are plotted in Figure 3.1. One can see that both controllers track the same way the reference, react identically to the external force and allow to satisfy the imposed constraints. The difference appears when approaching the position constraint. The capacity of MPIC to predict the activation of the constraint allows to initiate a smooth transition from free to constrained motion. In contrary, the classical QP controller only assesses the constraint one step before activation of this boundary, leading to an excessive braking torque. In real systems, such high torques could not be implemented due to hardware limitations, leading to constraint violation. Also, it is important to notice that increasing the prediction horizon of the MPIC controller allows better anticipation of the approaching constraint and thus an earlier and reduced braking torque. This clearly shows the advantage of using MPIC over constrained QP for controlling systems that require greater reactivity.

Experimental validation

In this section, IC and MPIC are compared throughout five different experiments. The first experimentation illustrates the equivalence of IC and MPIC when no constrains are applied. Then, position, velocity and input constraints are successively imposed to the MPIC controller for validation of its performances, and to illustrate the advantages it offers. Finally, a user interacts with the MPIC controlled robot.

Experimental setup

The performance of the proposed controller is evaluated with a collaborative robot following a reference trajectory with external disturbances resulting from interactions with a static environment or a human operator. At first, the experiments are performed on a simple environment in order to ensure good repeatability of the task, and so to provide comparable results. In a second part, interaction with the user occur during the execution of the task which corresponds to a scenario of collaborative manipulation or unexpected contact. The reference is a smooth polynomial trajectory on a planar surface, placed slightly under the table in order to ensure contact of the end-effector with the (x,y)-plane. In addition, a rigid obstacle is placed on the path generating unmodeled disturbances. The experimental setup is depicted in Figure 3.2.

Implementation and hardware specifications

The experiments are carried out with a KUKA iiwa 14 collaborative robot. A handle, equipped with a ATI-Mini40 force-torque sensor is mounted at the robot end-effector, allowing to interact with the environment and user and measure the end-effector wrench.

The proposed controller addresses the robot at joint-torque level, which is not a classical approach for controlling this type of manipulators. The internal structure of the low level controller of the KUKA iiwa is not well known, with some identification attempts such as in [START_REF] Chawda | IEEE International and distinct time window for recovery of arm motor control early after stroke revealed with a global measure of trajectory kinematics[END_REF]Niemeyer, 2017a, Chawda andNiemeyer, 2017b]. It is assumed that the low level control architecture behaves as described in [START_REF] Albu-Schäffer | A unified passivity-based control framework for position, torque and impedance control of flexible joint robots[END_REF]. The robot software implements joint level IC and allows the tuning of the impedance parameters. The embedded IC would however interfere with the MPIC implementation. As a result, these internal joint impedance parameters are set to zero and the torques commands generated by the MPIC are directly fed to the system using KUKA's FRI (Fast Robot Interface) protocol. This protocol is also used to measure the joint positions, whereas the joint velocities are computed using a filtered derivative with a cut-off frequency at 50Hz, which is appropriate given the slow robot motions. End-effector positions and velocities are computed using the forward kinematic and forward differential kinematic models, respectively. The control and the FRI communication rates are set at 500Hz. For solving the QP problem (3.26) with constraints, the C++ solver qpOASES [START_REF] Ferreau | qpOASES: a parametric active-set algorithm for quadratic programming[END_REF] is used. The impedance parameters are chosen such that M v = diag(5, 5, 5, 2, 2, 2) [Ns 2 /m, Nms 2 /rad for translation and rotation, respectively] , 300,300,50,500,500,500) [Ns/m, Nms/rad], and 77.5, 77.5, 31.6, 63.2, 63.2, 63.2)[N/m, Nm/rad]. K v is chosen in order to ensure: 1) a good contact in the direction normal to the (x,y)-plane, 2) a relatively low stiffness in the x and y directions, which represents a trade-off in compliance and tracking performance, and 3) a high stiffness in orientation in order to have a nearly constant wrist pose, in spite of external torques. D v is selected in order to have a critical damping ratio of 1 for the impedance model. The prediction horizon for the MPIC is set to H = 5, which corresponds to a prediction over 10ms given the sampling rate.

K v = diag(
D v = diag(

Experimental results

Comparing IC and MPIC without constraints : The IC law (2.7) is implemented on the system and tested in the same conditions as the MPIC. Figure 3.3 allows comparing the two controllers. At position level (Figure 3.3a) one can observe the almost perfect matching of both trajectories. The measured position error between IC and MPIC is respectively e ¯x = 0.24mm and e ¯y = 0.12mm in the x and y directions, with associated standard deviations σ x = 1.3mm and σ y = 0.7mm. This match can also be noticed at velocity (Figure 3.3b) and acceleration (Figure 3.3c) levels where only a small difference on the x-axis can be noted.

The difference is likely to be caused by some friction on the surface. One can also notice that both controllers do not perfectly track the reference trajectory, as the impedance model has a rather low stiffness. At some points, as it is the case in the left upper corner of the trajectory in Figure 3.3a, the performance of both controllers seems to be particularly reduced, even though they behave identically. This might be caused by the internal controller of the KUKA robot that has some issues tracking torque commands as pointed out in [START_REF] Chawda | IEEE International and distinct time window for recovery of arm motor control early after stroke revealed with a global measure of trajectory kinematics[END_REF]Niemeyer, 2017a, Chawda andNiemeyer, 2017b]. To verify the impact of these small defects due to the closed-structure low level controllers of KUKA iiwa, the real impedance parameters are estimated from experimental data, and compared to the ones chosen in the impedance model. A linear least square minimization is performed on the experimental data of IC to retrieve these parameters. The resulting estimations are given in Table 3.1. It turns out that the system impedance during the experiments fits rather well the imposed impedance model, with a 8% max error on M v , a 6.6% max error on K v , and a 6.1% max error on D v , for the x and y components. This result shows that the uncertainties on the internal controller do not heavily impact the performance and do not interfere with the robot-environment interaction. Typically, a human operator is not able to discriminate such differences that are smaller than the noticeable variations mentioned in the literature for inertia, damping and stiffness perception [Jones, 2000].

Position constrained MPIC : In this experiment, position constraints are specified to the MPIC. Figure 3.4 shows the resulting end-effector trajectories as well as the position constraint, which is set such that p x ≥ -0.006m. One can observe that the two controllers have the same behavior when working in the allowed task space. When approaching the limit, MPIC stops behaving as IC to meet the constraint. Velocity constrained MPIC : Velocity constraints can also be applied to the MPIC controller. Figure 3.5 compares the behavior of MPIC and IC at position and velocity level. In Figure 3.5b, one can see that given velocity constraints v x ≤ 0.02m/s, v y ≤ 0.02m/s, the MPIC is able to limit the system velocity below the desired threshold. As in the previous paragraph, this also results in the modification of the desired trajectory, as shown in Figure 3.5a. Input constrained MPIC : MPIC is also tested in the case of input constraints resulting in the limitation of the commanded acceleration of the system. Figure 3.6a allows comparing the accelerations obtained via IC and MPIC for a limit a y ≥ -0.4m/s 2 . Both controllers have an equivalent behavior when inside the allowed acceleration range and MPIC limits acceleration when exceeding the limit. It is worth noticing here that limiting the acceleration affects the contact force between the robot and the obstacle, as shown in Figure 3.6b. This feature is of particular interest, as it acts as a saturation that can limit the input to some maximally allowed threshold. In addition, it also acts directly on the system maximum acceleration and the resulting force applied by the robot, which is an important feature for interaction tasks. MPIC with position constraints and pHRI : Finally, the MPIC is tested in interaction with a human operator under position constraints. Figure 3.7 shows the resulting end-effector trajectories, the position constraint, which is set such that p x ≥ 0m, as well as the measured force applied on the robot end-effector and the input signal generated by the MPIC. One can see that while inside the allowed area, the system is under IC and thus reacts accordingly to the wrench applied by the operator and rejoins the reference trajectory when the operator stops interacting. It is however not possible to exit the allowed area, even when pushing against the boundary as the MPIC controlled system generates opposing forces (blue dashed areas in Figure 3.7).

The performed experiments shows that MPIC combines the features of IC with the advantages of MPC. MPIC behaves exactly as IC when no constraints are active and it satisfies position, velocity and acceleration constraints when such constraints are imposed, even in the case of pHRI. Also, the prediction capacity of MPIC is an advantage when smooth transitions between free and constrained motion are required and to increase the reaction time of the system when controlled with small sampling rates. MPIC can thus be applied to many applications, in which the desired interaction dynamics are considered constant. However, in some cases, the ability to modify the compliance of the system is required. This issue is further explored in the following chapter.

Chapter 4

Changing interaction dynamics

The development of collaborative robots designed to safely support operators during the execution of tasks has certainly played an important role in this evolution, with applications such as lifting of heavy tools and pieces, hand-guiding or collaborative assembly. In this context, the question of handling physical contact in the best possible way has become crucial. This can be obtained typically by controlling the system compliance during interaction between the robot and its environment, resulting in the classical IC scheme.

One of the main challenges in interaction management using IC lies in the proper selection of the impedance parameters, not only according to the task but also to ensure stability for all possible parameters variations. It is worth mentioning that using fixed impedance parameters makes the system passive and hence stable when interacting with a passive environment [Colgate and Hogan, 1988, Hogan, 1988], this property no longer holds for arbitrarily varying parameters [START_REF] Kronander | Stability considerations for variable impedance control[END_REF]. Not being able to modify the system impedance during the task is a notable drawback that may limit application ranges. For instance, human-robot interaction may benefit from the ability to adapt to human's force or fatigue, or to dynamically modifying the interaction [START_REF] Ikeura | Variable impedance control of a robot for cooperation with a human[END_REF].

In this chapter, after introducing the concept of variable impedance control (VIC) and reviewing some of its applications, the passivity issues that can result from the modification of the interaction dynamics are explored. Starting from state-of-the-art methods used to ensure passivity, a design method for passivity filters based on the combination of passivity conditions with an adaptation law on the impedance profile is proposed. The described method based on VIC with passivity guarantees is tested both in simulation and experimentally.

Variable impedance control

Meanwhile in classical IC the impedance profiles are kept constant, in VIC variations of these parameters are allowed resulting in a dynamical link between the system positions, velocities, accelerations and external forces applied on the system. Because of its ability to change the interaction dynamics during the task, VIC has recently attracted more and more attention. This control strategy has been explored for challenging objectives, for example to deal with explosive movements as described in [START_REF] Braun | Optimal variable stiffness control: formulation and application to explosive movement tasks[END_REF], or to optimize the performance of hammering tasks, as shown in [START_REF] Garabini | Optimality principles in variable stiffness control: The vsa hammer[END_REF]. In [START_REF] Medina | Risk-sensitive interaction control in uncertain manipulation tasks[END_REF], VIC was employed to allow risk-sensitive interactions and in [START_REF] Haddadin | Optimal control for maximizing link velocity of robotic variable stiffness joints[END_REF] it was used to maximize robot links velocities. VIC has also been implemented using reinforcement learning [START_REF] Buchli | Learning variable impedance control[END_REF], as well as adaptive approaches for human-robot collaboration based on the estimation of the human arm stiffness, from the derivatives of force and position [START_REF] Tsumugiwa | Variable impedance control based on estimation of human arm stiffness for human-robot cooperative calligraphic task[END_REF], or from the measurement of muscle activity with electromyography [START_REF] Grafakos | Variable admittance control in phri using emg-based arm muscles co-activation[END_REF]. However, in all these contributions, the resulting impedance of the robot is identified or learned during task execution. The issue of guaranteeing stability when modifying the interaction dynamics is not directly treated or worked around by means of stabilizing properties of the environment, such as the stabilizing properties of the human arm. Therefore these methods cannot be used to explicitly check before and during execution, whether a given impedance profile will lead to stable execution.

As in this thesis the passivity properties of a given system are used, the definitions of dissipativity and passivity, as defined in [Kugi, 2001, Schaft, 2017], are described here for further reference.

Consider a state-space system with inputs and outputs of the general form

x ˙=f (x, u) (4.1) y =h(x, u) (4.2)
where x ∈ R nx is the system state, u ∈ R nu the input and y ∈ R ny the output.

Definition 4.1. The system given by equations (4.1) and (4.2) is said to be dissipative with respect to the supply rate v(u, y) if there exist a non-negative function V (x), called storage function such that for all initial condition x(t 0 ) = x 0 at any time t 0 and for all t ≥ t 0 the following inequality holds

V (x(t)) -V (x 0 ) ≤ t t 0 v(u(τ ), x(τ ))dτ
Definition 4.2. The system given by equations (4.1) and (4.2) is said to be passive if it is dissipative with respect to the supply rate v(u, y) = u T y. That is if there exist a non-negative function V (x), called storage function such that for all initial condition x(t 0 ) = x 0 at any time t 0 and for all t ≥ t 0 the following inequality holds

V (x(t)) -V (x 0 ) ≤ t t 0 u(τ ) T y(τ )dτ

Passivity issues when modifying interaction dynamics

For VIC, a similar impedance behavior as for IC is imposed to the robotic system, however, the case where D v (t) and K v (t) are time-varying is considered. For notation simplicity, from now on, the impedance parameters M v , D v (t) and K v (t) are referred to as M , D(t) and K(t). The resulting system dynamics under IC are given by M e ¨p + D(t)e ˙p + K(t)e p = f ext (4.3) with e p = p r -p ∈ R m , where p is the robot end-effector pose and p r is a reference motion.

The interaction is then characterized by the impedance resulting from the apparent virtual mass M ∈ R m×m , the desired damping D(t) ∈ R m×m and the desired stiffness K(t) ∈ R m×m . In this work, M is assumed to be constant. Note that M , D and K are user defied symmetric positive definite matrices and therefor can be chosen diagonal to simplify developments. Also from now on, time dependency of the stiffness and damping will generally not be mentioned in the notations for simplicity. Note here that even though the impedance model (4.3) and further considerations are made in task-space and given in terms of the robot end-effector pose p, the same results apply to joint-space.

In order to investigate the stability of VIC given by (4.3) with varying stiffness and damping terms, the following non-negative storage function is considered

V 1 = 1 2 e ˙T p M e ˙p + 1 2 e T p Ke p (4.4)
In this case, it is worth noticing that the storage function has a physical interpretation, namely it consists of the kinetic energy of the system combined with potential energy of the spring characterized by its stiffness K. As K and M are symmetric matrices, differentiating V 1 leads to

V ˙1 = e ˙T p M
V 1 (t) -V 1 (0) ≤ t 0 e ˙T p f ext dτ (4.6)
which can be interpreted as a passivity condition, as defined in [Schaft, 2017].

However, because of stiffness time-dependency in equation (4.5), the sign of the term in brackets may change, possibly leading to a violation of the passivity condition (4.6) in case of increasing stiffness. Passivity can only be guaranteed if the stiffness is either constant or decreasing. In this case, the storage function can be used as a Lyapunov-like function to show stability. Also, it is well known that the interaction with a passive environment yields a passive, thus stable, interaction dynamics. Though in practice varying the stiffness does not necessarily lead to an unstable behavior, guaranteeing stability is a critical issue for ensuring user integrity. Thus, several contributions addressed this issue of guaranteeing a stable execution of VIC, with a given impedance profile. In [START_REF] Ficuciello | Variable impedance control of redundant manipulators for intuitive human-robot physical interaction[END_REF], the authors showed that exploiting kinematic redundancies can ensure passivity for a larger panel of impedance profiles. Gain scheduling control has also been used to address VIC by interpolating feedback gains between operating points [START_REF] Lee | Force tracking impedance control with variable target stiffness[END_REF]. Virtual energy-storing tanks [START_REF] Ferraguti | A tank-based approach to impedance control with variable stiffness[END_REF] are another solution to modify the impedance model and guarantee passivity, by storing energy dissipated by the system. This virtually stored energy can then be re-injected to simultaneously implement stiffness variations and to guarantee a stable execution. This approach has been successfully implemented in several applications, as in [START_REF] Ferraguti | An energy tank-based interactive control architecture for autonomous and teleoperated robotic surgery[END_REF] where it was used to handle control transitions in teleoperated robotic surgery. Also, in [START_REF] Schindlbeck | Unified passivity-based cartesian force/impedance control for rigid and flexible joint robots via task-energy tanks[END_REF] the authors introduce taskenergy tanks in a hybrid force/impedance control strategy for robustly handling the loss of contact which was then used for polishing tasks. In [START_REF] Ferraguti | A tank-based approach to impedance control with variable stiffness[END_REF], an energystoring tank was used to deal with the interactions resulting from tissues and membranes penetration during a robotic needle insertion. The tank-based strategy has been shown very well suited for VIC, in spite of some difficulties to tune its parameters. In fact, the tuning of the tank parameters heavily impacts the performance of the control, which can make the system stuck if not enough energy is available. Also, in [START_REF] Lutscher | Hierarchical force and positioning task specification for indirect force controlled robots[END_REF], the concept of virtual energy storage was explored in order to derive velocity constraints for passive set-point modulation. Nevertheless, these methods are dependent on the states of the system, measured during task execution and so can only be applied online.

Passivity conditions

Whereas strategies based on energy tanks are well suited to guarantee passivity during the execution of the task, for some applications the impedance profiles are well defined before the execution. Thus these applications would highly benefit from a method that could verify the passivity of the profiles before execution on the real system. Some authors treated this issue and proposed an approach that can be used to check the passivity of a given time-varying impedance profile before the execution of the task [START_REF] Kronander | Stability considerations for variable impedance control[END_REF]. In their method, the authors derive passivity conditions that allow to check whether the interaction will have passive dynamics.

In fact, as shown previously, the storage candidate function (4.4) is very conservative, as it only allows a constant or decreasing stiffness, and should not be used to conclude on system passivity for increasing stiffness. Alternative storage functions need to be explored. In order to facilitate the passivity analysis, [START_REF] Behal | Lyapunov-based control of robotic systems[END_REF] propose to use a filtered tracking error-like variable r ∈ R m such that r = e ˙p + αe p (4.7)

where α ∈ R + is a constant gain. As shown in [Dawson, 2019], if r is bounded then e p , e ˙p are bounded and, in this case, if r → 0, then e p , e ˙p → 0. With this in mind, it is possible to reformulate the control objectives using r

1 2 C ˙)e p + e T p (K + αD -α 2 M -C)e ˙p}dτ (4.9)
which can be written as

t 0 r T f ext dτ = V 2 (t) -V 2 (0) + t 0 W 2 dτ (4.10) with V 2 = 1 2 r T M r + 1 2 e T p Ce p (4.11)
and

W 2 = e ˙T p (D -αM )e ˙p + e T p (αK - 1 2 C ˙)e p + e T p (K + αD -α 2 M -C)e ˙p (4.12)
It is worth noticing that if C is defined such that C = K + αD -α 2 M , the last term of (4.12) corresponding to the cross-terms of e p and e ˙p becomes zero. In this case, the candidate storage function (4.11) is the same as proposed in the work of [START_REF] Kronander | Stability considerations for variable impedance control[END_REF]. The resulting passivity conditions are then given, such as in [START_REF] Kronander | Stability considerations for variable impedance control[END_REF], by 0 ≤D -αM (4.13)

0 ≤2αK -K ˙-αD ˙(4.14)
In fact, this particular choice leads to W 2 ≥ 0, allowing to conclude on system passivity, as in this case, equation (4.10) leads to

V 2 (t) -V 2 (0) ≤ t 0 r T f ext dτ (4.15)
Passivity conditions (4.13) and (4.14) can be used in order to verify whether passivity is guaranteed for a given impedance profile. However, they do not give any information about how to modify a profile that is not passive. In the following development, a design method capable of guaranteeing passivity properties for any given impedance profile based on conditions (4.13) and (4.14) will be proposed.

Passivity filter

Such as introduced in the previous section, being able to guarantee the passivity of the interaction when controlling a robot with VIC and variable impedance profiles is of particular importance. As it was shown, the literature proposes two main strategies to ensure passivity: online guarantees resulting form energy tanks and offline guarantees given by passivity conditions. The aim of the following development is to propose a third strategy that combines advantages of both methods. The design of passivity filters is described, which make it possible to check the passivity of a given profile or modify a non passive one. Based on the combination of passivity conditions with an adaptation law on the impedance profile, the proposed method can thus be used either offline, before the execution of the profile, or online.

Guaranteeing passivity using passivity conditions

In the further development it is considered that the change in system impedance, with the stiffness (respectively the damping) is varying between K 0 and K 1 (respectively between D 0 and D 1 ). Even though this hypothesis could appear at first conservative, it is valid for most applications. In fact, the stiffness and damping profiles are physically required to be positive and finite and thus they are always bounded. In the following development, the case K 0 < K 1 is considered without a loss of generality. The impedance profiles K and D are defined such that

K = K(t) = K 0 + Γ(t)δK (4.16) D = D(t) = D 0 + Γ(t)δD (4.17)
with δK = K 1 -K 0 , δD = D 1 -D 0 and Γ a diagonal matrix where all diagonal terms γ i are differentiable gains such that 0 ≤ γ i ≤ 1. Γ is further referred to as the switching variable.

As M , K, D were chosen to be diagonal, the impedance behavior can be decoupled. In order to ensure passivity using conditions on V 2 , equations (4.13) and (4.14) need to be satisfied for the desired profiles of K and D. Condition (4.13) allows choosing α easily. For instance, in [START_REF] Kronander | Stability considerations for variable impedance control[END_REF] it is chosen as the ratio between the smallest eigenvalue of D and the largest eigenvalue of M . This allows to chose a constant value of α that always satisfies condition (4.13). Condition (4.14) implies that for ∀i ∈ {1, 2, . . . , m}

γ i ˙(αδd i + δk i ) -2α(k 0i + γ i δk i ) ≤ 0 (4.18)
with δk i and δd i the diagonal terms of δK and δD respectively, and k 0i the diagonal terms of K 0 .

Generally, an increase in the desired stiffness does not come with a decrease in the desired damping, which would be the most constraining case in terms of stability. Then, it is assume that δd i ≥ 0. As δk i > 0 and α ≥ 0, then αδd i + δk i > 0. Equation (4.18) can then be written as

γ i ˙≤ a i γ i + b i (4.19) with a i = 2αδk i αδd i + δk i and b i = 2αk 0i αδd i + δk i
Combining equation (4.18) with a low-pass filter, in order to ensure smoothness of γ i , a passivity filter can be designed in such a way that it takes as input the desired switching profile γ ¯i and generates an output profile γ i that guarantees the system passivity. The filter can then be defined by

γ ˙i = min(β(γ ¯i -γ i ), a i γ i + b i ) (4.20)
where the integration of (4.20) permits to obtain γ i that tracks γ i ¯in such a way that passivity conditions are respected. The filter structure is shown in Figure 4.1, with β the filter gain. Note also that in the less conventional case where the damping decreases while the stiffness increases, the previous reasoning can still be applied, as long as αδd i + δk i > 0, that is δd i > -δk i /α. The case αδd i + δk i = 0 is a singularity that should be avoided.

Though the usage of the proposed filter may guarantee passivity, the tracking performance highly depends on the conservatism of the passivity conditions. The constant value α plays an important role in the conditions resulting from function V 2 . In [START_REF] Kronander | Stability considerations for variable impedance control[END_REF] it is chosen in the worst case, thus reducing the variation range of the parameters. That is why it is of particular interest to explore other candidate storage function in order to find less conservative passivity conditions resulting in better tracking performance of the switching function γ ¯.

Enhanced passivity conditions

As the previously defined storage function gives conservative passivity conditions, it is of interest to find a novel candidate storage function that may allow to conclude on system passivity in a less conservative way. To do so, a first step consists in modifying the assumptions on (4.7) to replace α by a time-dependent matrix B, leading to a novel filtered tracking error variable such that

r = e ˙p + M -1 Be p (4.21)
where B ∈ R m×m is a time-dependent diagonal matrix with bounded eigenvalues. The following lemma justifies the use of the tracking error given by (4.21) and shows that by Proof. As the matrices M and B are chosen to be diagonal, the problem can be decoupled accordingly to the diagonal terms. Hence ∀i ∈ {1, 2, . . . , m}:

⇒ if r i is bounded then e p i , e ˙pi are bounded : the solution of the differential equation

r i = e ˙pi + b i m i e p i (4.22)
is given by

e p i = e -t 0 b i (τ ) m i dτ e p i (0) + t 0 e -t τ b i (σ) m i dσ r i (τ )dτ (4.23) As b i ∈ R + is considered to be bounded such that b i ≤ b i ≤ b ¯i, (4.23
) can be upper bounded by

|e p i | ≤ |e p i (0)| + t 0 e - b i m i (t-τ ) |r i (τ )|dτ ≤ |e p i (0)| + sup t (|r i (τ )|) t 0 e - b i m i (t-τ ) dτ ≤ |e p i (0)| + sup t (|r i (τ )|) m i b i (1 -e - b i m i t ) (4.24)
what shows that if r i is bounded then e p i is bounded. Additionally, using (4.22) one get that

|e ˙pi | ≤ b ¯i m i |e p i | + |r i (t)| (4.25)
that can be written using (4.24) m i is constant, following the same steps of the proof of Lemma 1.6 in [Dawson, 2019], one can show that ||e [Dawson, 2019], as e p i , e ˙pi and ||e p i || 2 are bounded, then lim t→∞ e p i = 0. Finally, from (4.25) with lim t→∞ e p i = 0 and lim t→∞ r i = 0, one can conclude that lim t→∞ e ˙pi = 0.

|e ˙pi | ≤ b ¯i m i |e p i (0)| + sup t (|r i (τ )|) b ¯i b i (1 -e - b i m i t ) + |r i (t)| what shows that if r i is bounded then e ˙pi is bounded. ⇒ if r i → 0 then e p i ,
- b i m i t e p i (0)|| 2 is bounded and that || t 0 e - b i m i (t-τ ) r i (τ )dτ || 2 is bounded as ||r i || 2 is bounded. This shows that ||e p i || 2 is bounded. Using Corollary 1.1 from
In this case, using the same strategy as previously to obtain equation (4.9) yields

t 0 r T f ext dτ = t 0 {r T M r ˙+ e T p Ce ˙p + 1 2 e T p C ˙ep }dτ + t 0 {e ˙T p (D -B)e ˙p + e T p (BM -1 (K -B ˙) - 1 2 C ˙)e p + e T p (K + BM -1 (D -B) -B ˙-C)e ˙p}dτ (4.26) that is t 0 r T f ext dτ = V 3 (t) -V 3 (0) + t 0 W 3 dτ (4.27) with V 3 = 1 2 e ˙p e p T M B B C + BM -1 B e ˙p e p
(4.28)

and

W 3 = e ˙T p (D -B)e ˙p + e T p (BM -1 (K -B ˙) - 1 2 C ˙)e p +e T p (K + BM -1 (D -B) -B ˙-C)e ˙p (4.29)
The choice C = K + BM -1 (D -B) -B ˙allows canceling out the cross-terms. The positivity of V 3 is obtained as a consequence of Schur complement Lemma [Boyd, 1994] that states that supposing Q and R are two symmetric matrices, the linear matrix in-equality (LMI) given by

R S T S Q ≥ 0 is equivalent to R ≥ 0 and Q -SR -1 S T ≥ 0.
In the case of V 3 , this results in the condition 0 ≤ C (4.30)

According to equation (4.29) and the particular choice of C, the positivity of W 3 is obtained if the following novel passivity conditions hold

0 ≤ D -B (4.31) 0 ≤ BM -1 (K -B ˙) - 1 2 C ˙(4.32)
It is worth noticing that the storage function V 2 of equation (4.11), initially defined by [START_REF] Kronander | Stability considerations for variable impedance control[END_REF], is a special case of V 3 defined in equation (4.28). It corresponds to the case B = αM , C = K + αD -α 2 M and α constant. The analogy between the scalar α in equation (4.7) and the matrix M -1 B in equation (4.21), and the relationship between αM and D resulting from the condition (4.13) suggest that the candidate storage function could be modified to limit the conservatism of the worst-case choice for α. This leads to set B = D and C = K -D ˙. A fourth storage function is then defined, such that

V 4 = 1 2 e ˙p e p T M D D K + DM -1 D -D ˙ e ˙p e p (4.33)
and the resulting passivity conditions are deduced from equations (4.30) to (4.32)

0 ≤ K -D ˙(4.34) 0 ≤ 2DM -1 (K -D ˙) -K ˙+ D ¨(4.35)

Passivity filter using enhanced passivity conditions

Based on the new passivity conditions (4.34) and (4.35), resulting from the storage function V 4 given by equation (4.33), a new passivity filter can be designed. As previously, all matrices are chosen diagonal and the problem can be decoupled, thus (4.35) is equivalent to

k ˙i + 2 m i d i d ˙i -d ¨i - 2 m i k i d i ≤ 0 resulting in -a 1i γ ¨i + a 2i γ ˙i + a 3i γ i γ ˙i -a 4i γ i -a 5i γ 2 i -a 6i ≤ 0 (4.36) with a 1i = δd i a 2i = δk i + 2 m i d 0i δd i a 3i = 2 m i δd 2 i a 4i = 2 m i (d 0i δk i + k 0i δd i ) a 5i = 2 m i δk i δd i a 6i = 2 m i k 0i d 0i
which are all positive since both δK and δD are positive. Then, as a 2i > 0 and γ i is such that 0 ≤ γ i ≤ 1, ∀γ i , a 2i + a 3i γ i > 0. This allows to write the passivity condition

γ ˙i ≤ a 1i γ ¨i + a 4i γ i + a 5i γ 2 i + a 6i a 2i + a 3i γ i h 1i (γ ¨i, γ i ) (4.37)
Additionally, in the case where δd i = 0, equation (4.34) gives a second passivity condition for γ ˙i such that

γ ˙i ≤ δk i δd i γ i + k 0i δd i h 2i (γ i ) (4.38)
Using passivity conditions (4.37) and (4.38), a new passivity filter similar to (4.20) could be defined such that

γ ˙i = min(h 1i (γ ¨i, γ i ), h 2i (γ i ), β(γ ¯i -γ i )) (4.39)
The filter structure is shown in Figure 4.2.

Condition (4.37) is however not well suited for practical implementation as it requires derivating γ ˙i, which is the output of a min-switch and therefor can be discontinuous. One could think that one way to overcome this issue would be to introduce filtered derivatives that would smooth the output. These processing steps, however, introduce additional dynamics resulting in delays, that could result in breaking the passivity conditions. For this reason, a more conservative but practically implementable constraint is preferred, such that

γ ˙i ≤ a 4i γ i + a 5i γ 2 i + a 6i a 2i + a 3i γ i h 3i (γ i ) (4.40)
In fact, when the constraint is active, it results that equation (4.40) becomes an equality and necessarily γ ˙i ≥ 0. In this case it can be shown that

γ ¨i = 2a 5 a 2 γ i + a 5 a 3 γ 2 i + (a 2 a 4 -a 3 a 6 ) (a 2 + a 3 γ i ) 2 γ ˙i
and

a 2 a 4 -a 3 a 6 = 4 m 2 i δk i ( m i 2 d 0i + δd i k 0i + d 0 2 i δd i ) ≥ 0
resulting in γ ¨i ≥ 0. This leads to the fact that if γ ˙i respects the constraint given by equation (4.40), then it also respects the passivity condition (4.37). The proposed passivity filter resulting from conditions (4.38) and (4.40) is defined by

γ ˙i = min(h 2i (γ i ), h 3i (γ i ), β(γ ¯i -γ i )) (4.41)
Finally, the filter structure is shown in Figure 4.3.

Experimental validation

For validating the proposed approach, the passivity filter is tested in simulation and on an experimental setup featuring a KUKA iiwa robot. In both cases, the controlled system is considered to be operating in task-space under VIC.

Simulations

In order to illustrate the advantages of the proposed approach, simulations analyzing passivity properties on a one-dimensional mass-spring-damper system were performed. One can see that even if the outputs of both filters first evolve similarly, the output of filter (4.41) converges much more rapidly to the reference switching variable. This illustrates the fact that the conditions given by V 4 are less conservative than these given by V 2 . Note that no external force is applied to the system, and in this case passivity is guaranteed if the derivative of the storage (Lyapunov-like) function is negative. Comparing the evolution of V ˙1, V ˙2 and V ˙4 in the case where γ is calculated using passivity filter (4.20) (based on V 2 ), as shown in Figure 4.5, one can see, that in contrary to function V 1 that is conservative, V 2 and V 4 guarantee passivity. On the other hand, using passivity filter (4.41) (based on V 4 ), one can see in Figure 4.6 that functions V 1 and V 2 cannot guarantee passivity, whereas V 4 is capable of doing so.

Experiment

The experimental validation of the passivity filter is performed on a KUKA iiwa 14 collaborative robot. In the experimental setup, an ATI-Mini 40 force-torque sensor is mounted on the robot end-effector to measure the external wrench. The external wrench is applied by means of an elastic stretch band with a stiffness of around 180N/m and fixed along the Cartesian y-axis. The system is controlled with a Cartesian VIC in the direction of the external force and all other DoFs are controlled with a classical IC with constant parameters and a high stiffness of 1000N/m and a damping ratio of 1. Redundancy is controlled by applying a null-space stiffness of 300Nm/rad and the damping ratio of 1 with the initial joint position as reference. All inertial terms are set to 5Ns 2 /m and 5Nms 2 /rad for Cartesian and joint level controllers, respectively. The generated compliant acceleration commands are transformed into position commands by double integration In order to show the advantages of the proposed strategy, the experimentation is conducted in two phases. In the first part, the VIC is directly given the impedance switching profile, as in the second part, the switching is checked and modified by the passivity filter.

Varying impedance profiles without the passivity filter

In the first part of the experiment, the system is tested without the use of the passivity filter. Experimental results are shown in rapidly looses passivity after the increase in the impedance profile around t = 5s. In fact, passivity is only guaranteed as long as the relationship V ˙4 ≤ r T f ext holds, which is not the case here. Indeed, after the injection of energy due to the increase in stiffness the variation of the stored energy does not only increase in time but also breaks the passivity of the system. The end-effector trajectory is depicted in Figure 4.8b. On can easily note that this loss of passivity results in a violent motion of the robot end-effector that is only stopped when the robot safety stop is activated. These results clearly show that in the case where the impedance parameters are chosen without the preoccupation about the system passivity, it can create serious safety issues for the robot and its environment including a potential operator.

Varying impedance profiles using the passivity filter

In the second part of the experiment, the system is tested using the passivity filter. Experimental results are shown in Figure 4.10. In the same way as in the previous experiment, the variable impedance parameters are switched according to It is worth noticing how the variation of the impedance parameters K and D is modified by the passivity filter resulting in a slower increase compared to the reference profiles K ¯and D ¯. The influence of the passivity filter on the evolution of the switching variable can be observed in Figure 4.9b. One can see that the derivative of the switching variable γ ˙is following the passivity condition h 3 when increasing the impedance profile and evolves freely as long as it is under the passivity conditions. One can also notice that, in this case, the condition h 2 is much higher than h 3 and thus is never reached. The ability of modifying the interaction dynamics in a safe way opens many possibilities for enhanced pHRI, which is further explored in the following chapter.

Chapter 5

Enhanced human-robot interactions

The evolution of robotics and recent advances in the field of pHRI tend toward human oriented design of robotic mechanisms and control strategies. In this perspective, characterizing the behavior of the human operator comes to the fore. Even if a fully automated robotic system has generally a greater performance in precision, repeatability and load capacity than a human operator when operating in a well defined environment, many tasks require the human ability of fast judgment and adaptation in case of unpredicted events. Collaborative systems are designed to combine both the robot force and precision, and the operator's judgment and flexibility. One of the main challenges for effective collaboration between the human and the robot remains in the communication between both sides and especially in the prediction of the operator intentions [START_REF] Ajoudani | Progress and prospects of the human-robot collaboration[END_REF]. It is of particular interest in robot-assisted tasks, in which the robot is able to detect the extend of the operator activity during task execution and adapt the level of assistance when needed. This, however, requires to design and establish communication standards that are simultaneously intuitive for the user and give the robotic system sufficient information on the human intentions in various phases of the collaboration [START_REF] Sebanz | Joint action: bodies and minds moving together[END_REF], Bauer et al., 2008].

In this chapter, after reviewing some state-of-the-art strategies for human intention perception for human-robot collaboration, the focus is shifted towards the use of measurable bio-signals and in particular towards electromyography (EMG). It is explained how these signals can be acquired and used for generating models of the human arm and hence estimating operator intentions. Based on this information, an EMG driven VIC strategy is proposed and tested experimentally. In the following development it is considered that the user interacts with the robot with his upper limbs.

Human intention perception

In a context where enhanced human-robot interactions are required, the understanding and transposition of communication channels from human-human interactions to robotic systems comes to the fore [START_REF] Reed | Physical collaboration of human-human and human-robot teams[END_REF]]. Among the well known communication standards allowing to interpret human intentions, visual and language communication methods appear to be most user-friendly as they are most intuitive. The use of body gesture is a way to interpret human intentions in a vision-based framework and trigger robot assistance when needed [START_REF] Hawkins | Probabilistic human action prediction and wait-sensitive planning for responsive humanrobot collaboration[END_REF]. However, even though these communication methods appear intuitive to the human user, the information they provide is not rich enough to fully assess the human intention and thus can only be used to trigger some high level action on the robot side [START_REF] Ajoudani | Progress and prospects of the human-robot collaboration[END_REF].

The observations made in [START_REF] Reed | Physical collaboration of human-human and human-robot teams[END_REF] show that during physical collaboration between two human partners, an additional haptic communication channel emerges. Based on the assessment of the force inputs of both partners, human activity can be evaluated by measuring the extend of force applied by the operator to the robotic system. This can typically be done using mechanical force-torque sensors. Due to its simplicity, this strategy is commonly used in many collaborative applications such as object manipulation [START_REF] Ikeura | Variable impedance control of a robot for cooperation with a human[END_REF], Tsumugiwa et al., 2002, Duchaine and Gosselin, 2008]. The major drawback of such a method consists in the fact that the sensor not only measures the desired human-robot interaction. Undesired components such as gravity and friction forces are also measured, as well as unmodeled interaction forces should some external contact with an uncertain environment occur [START_REF] Ajoudani | Progress and prospects of the human-robot collaboration[END_REF]. This makes the use of such sensors unsuited for complex tasks with the simultaneous need of dissociating the human intention and the environmental forces, as it is impossible to make any distinction between them.

Alternatively, measurable bio-signals such as skin surface electromyography (EMG) can also be used for human activity assessment, particularly in medical applications for which they could provide useful information. Indeed, this type of signals directly reflects the human intentions [START_REF] Bi | A review on emg-based motor intention prediction of continuous human upper limb motion for human-robot collaboration[END_REF].

Electromyography

Electromyography (EMG) is a measurement technique, usually used for diagnostic, which consists in studying the muscle function by acquiring and analyzing the electrical signals the muscle emits during contraction. Since its discovery in 1912 [Piper, 1912], EMG has found many fields of application, especially in biomechanics, bioengineering and medicine. It is used in neurophysiology [START_REF] Merletti | Surface Electromyography: Physiology, Engineering, and Applications[END_REF], ergonomics [START_REF] Shafti | Analysis of comfort and ergonomics for clinical work environments[END_REF], movement [START_REF] Türker | Surface electromyography in sports and exercise[END_REF]] and gait analysis [START_REF] Bing | Electromyography (emg) and human locomotion[END_REF] or rehabilitation [START_REF] Maciejasz | A survey on robotic devices for upper limb rehabilitation[END_REF], to mention just a few.

The principle behind EMG consists in the recording of signals, which are the result of physiological processes involved in human muscle activity [START_REF] Farina | Biophysics of the Generation of EMG Signals[END_REF]. Muscle contraction involves the shortening of the muscle fibers, which is achieved by a change in the fiber membrane potential. This variation of potential can be recorded by a proper placement of bipolar electrodes on the skin, just above the targeted muscle. Because a muscle consists of many fibers generating potential gradients, the signal acquired by the electrodes is a superposition of spatial and temporal distributions of all involved fibers. This leads to a noisy signal that requires specific post-processing in order to assess the muscle activity. The processing of EMG signals consists of three stages [START_REF] Merletti | Surface Electromyography: Physiology, Engineering, and Applications[END_REF]. Band-pass filtering is used to reject acquisition noise and signal drifts, followed by full-wave rectification, as only the signal amplitude has a physical meaning. Finally, envelope extraction is performed using a low-pass filter. The resulting signal directly reflects the muscle activation. An example of EMG data acquired on the biceps brachii, as well as the post-processed muscle activation signal are plotted in Figure 5.1. A detailed description on the principle of EMG and signal processing steps can be found in Appendix A.

EMG for human-robot interactions

As EMG signals are the result of muscle cell bio-electrical activity, they are often used in control strategies of assistive robotic systems since they directly reflect user's muscle voluntary contractions and give an insight on his intentions.

In many applications involving EMG-controlled robots, EMG signal magnitude is used as a trigger to enable some level of assistance of the robotic system [START_REF] Krebs | Rehabilitation robotics: Performance-based progressive robot-assisted therapy[END_REF]]. An alternative to this strategy is to use bio-signals directly for robot control. In recent years, much research work has been done to use bio-signals for robot control [START_REF] Bi | A review on emg-based motor intention prediction of continuous human upper limb motion for human-robot collaboration[END_REF], Li et al., 2020]. Such a robot control strategy requires however the design of appropriate controllers to achieve good performance while guaranteeing stability and security of the human-robot collaboration. Thus, controllers based on impedance control are very well suited for these types of applications.

Recent implementations of collaborative control strategies based on VIC include information about muscle activity using EMG. In [START_REF] Grafakos | Variable admittance control in phri using emg-based arm muscles co-activation[END_REF] the operator muscle activity is used to switch the robot damping between predefined values. The authors of [START_REF] Kiguchi | An emg-based control for an upper-limb power-assist exoskeleton robot[END_REF]] also use a VIC strategy to control a power-assisted exoskeleton robot by adapting the impedance parameters using the simultaneous contraction of antagonist muscles. In [START_REF] Peternel | Teaching robots to cooperate with humans in dynamic manipulation tasks based on multi-modal human-in-the-loop approach[END_REF], the concept of tele-impedance is used in a robot teaching approach. This approach is expanded in the work of [START_REF] Peternel | A humanrobot co-manipulation approach based on human sensorimotor information[END_REF], which introduces a framework based on EMG and manipulability assessment of the human arm to provide information about the human behavior and task requirements.

In these applications, the muscle activation information is directly used in the proposed control strategy. Alternatively, muscle activation signals can also be used to estimate other characteristics of the human limb, such as its motion or applied forces. This however requires to build EMG-driven models of the human limb.

EMG-based human modeling

In order to be able to properly describe the dynamics of the human limb based on the acquired EMG signals, it is important to establish a mapping relationship between the wrench generated by the arm muscles and the corresponding EMG signals. In [START_REF] Artemiadis | Emg-based control of a robot arm using low-dimensional embeddings[END_REF] it is pointed out that estimation methods play an important role in controlling the compliance of the robotic system. One of the problems of EMGdriven models consists in determining such EMG to force models. Much research effort has been made in investigating the EMG-Force relationship of the human arm and two main types of models were proposed: phenomenological models as the Hill model [Hill, 1938] and dynamic models.

Phenomenological models describe the musculoskeletal system at the muscle and tendon tissue level by the means of nonlinear springs and contractile elements [Hill, 1938].

The advantage of using phenomenological models is that they allow to better understand the underlying physiological processes. The complexity of these type of models consists in the high variability of the subject dependent parameters describing the system behavior. Moreover, in many cases, the parameters cannot be correctly measured on the studied subject and thus need to be estimated. However, despite their complexity, EMG-driven phenomenological models have successfully been used for estimating the forces produced by upper [START_REF] Holzbaur | A model of the upper extremity for simulating musculoskeletal surgery and analyzing neuromuscular control[END_REF] and lower [START_REF] Lloyd | An emg-driven musculoskeletal model to estimate muscle forces and knee joint moments in vivo[END_REF] limbs. Nevertheless, because of the high complexity of the analyzed system, using the model resulting in expensive computation. Hence, other less complex modeling methods are preferred for real-time applications.

Less complex models linking EMG with produced torques and forces were proposed. In [Clancy et al., 2012], polynomial models are used to map the EMG/torque relationship about the elbow. In their work, the extension and flexion EMG data is used in four dynamic model structures : linear, polynomial, Hammerstein and Wiener models. Hammerstein and Wiener models consist of a series of static nonlinear and dynamic linear models [Wiener, 1966]. Moreover, [START_REF] Hashemi | Emg-force modeling using parallel cascade identification[END_REF] showed that parallel-cascade Wiener models have a good performance in estimating the arm end-point force. Alternatively, artificial intelligence algorithms can also be used, such as in [START_REF] Loconsole | An emg-based approach for on-line predicted torque control in robotic-assisted rehabilitation[END_REF] where time-delayed artificial neural networks are used to compute the human joint-torque form EMG. These models, however, turn out to be highly nonlinear and thus computationally complex. In [START_REF] Pesenti | Linear parameter-varying identification of the EMG-force relationship of the human arm[END_REF] a study on modeling the EMG to force relationship is proposed using Linear Parameter Varying (LPV) and Linear Time Invariant (LTI) models. It is pointed out that these types of models are a fair compromise in terms of accuracy and complexity. As an application example, the work presented in [START_REF] Ajoudani | Teleimpedance: Teleoperation with impedance regulation using a body-machine interface[END_REF] uses LTI models to estimate human arm endpoint force and stiffness in the context of transferring human impedance characteristics to a robotic system. Also, in the field of assistive robots, [START_REF] Teramae | EMG-based model predictive control for physical human-robot interaction: Application for assistas-needed control[END_REF] proposed to use EMG signals in combination with MPC and a linear human arm model to provide robotic assistance in the case of a disabled patient unable to follow a desired trajectory.

However, in all these applications the model parameters are considered constant over time. For a more reliable prediction, the model needs to account for some variability of the human arm dynamics. One of the main reasons is because EMG signal parameters highly depend on processing methods and are often subject to changing conditions, such as electrode displacement, muscle fatigue or sweating [START_REF] Farina | The extraction of neural information from the surface emg for the control of upper-limb prostheses: Emerging avenues and challenges[END_REF]. For these reasons, EMG based linear models of human limbs are not deterministic and thus need to be recalibrated frequently. One way to deal with such variability is addressed in [START_REF] Hahne | Concurrent adaptation of human and machine improves simultaneous and proportional myoelectric control[END_REF], where an adaptive regression is proposed in order to improve the online performance of a linear EMG to joint-angle model for EMG-based prosthetic control.

EMG-based VIC

In many robotic applications it is important to be able to identify the operator participation into some human-robot shared task. In the previous sections, it has been shown that EMG-based robot control strategies can be used to address this issue. It has also been shown that performance and stability are an important requirement that can be assessed by the use of control strategies such as VIC. Even though classical VIC cannot guarantee interaction stability when varying the impedance parameters, the method proposed in Chapter 4 can be used to overcome this issue. The importance of human arm models has also been addressed. While there exist many types of models, the computational complexity of many of them does not make them suitable for real-time human force assessment. Therefore, LTI models seem to represent the best compromise between complexity and performance. Another issue that has been addressed is the variability of the model that needs to be accounted for, and it has been shown in the literature that adaptive regression could be a solution to deal with model variability.

For all these reasons, in this section, a robot control strategy is proposed that is based on VIC featuring a passivity filter and a recursively estimated linear human arm model.

Recursively estimated linear human arm model

In this section, the use of a linear time-varying model that describes the EMG-force relationship is explored. In order to account for the subject and time variability of the EMG-based human arm dynamics, a recursive model estimation method is described. The model is then evaluated in two validation setups, which use measurements of shoulder and elbow movements.

Recursive model estimation

The EMG-force relationship, expressed as a linear, discrete time, input/output model is given by

y ˆk = - na i=1 A i,k y k-i + n b i=1 B i,k u k-i-n k +1
(5.1) with y ˆk the predicted output, y k and u k the measured output and input at time step k, respectively and A i,k and B i,k the time varying model parameter matrices. The number of past inputs and outputs used in the model are given by n b and n a , n k is an input to output delay. This model structure implies that the current output y ˆk can be predicted using a weighted sum of past and current inputs as well as past outputs.

As the model given by equation (5.1) is linear, an identification of the model parameters can be performed using least square regression min

Θ ˆ||Φ T k Θ ˆk -y k || 2 2 (5.2) with Θ ˆk = [-A 1,k , ..., -A na,k , B 1,k , ..., B n b ,k
] the vector of the estimated parameters at time step k and Φ k the regressor.

The minimization problem (5.2) can be recursively solved using the RLS method [ Åström and Wittenmark, 2013] with a forgetting factor 0 < λ < 1, such that the estimation parameter update law can be written as

Θ ˆk = Θ ˆk-1 + L k e k (5.3)
with the estimation error

e k = y k -Φ T k Θ ˆk-1 (5.4)
and the gain matrix

L k = P k-1 Φ k λ + Φ T k P k-1 Φ k (5.5)
where P k is the covariance matrix of the parameters such that

P k = 1 λ (P k-1 -L k Φ T k P k-1 ) (5.6)
The estimated model (5.1) using update law (5.3) can be then expressed in terms of the predicted output y ˆk and the estimation error e k such that

y ˆk = - na i=1 A i,k (y ˆk-i + e k-i ) + n b i=1 B i,k µ k-i-n k +1
(5.7)

By doing so, in the case where the output can be measured, the error is computed and used to update the model. However, when the output is not measurable, the error is set to 0 and the model uses the inputs and previous estimates of the output to do the estimation.

Model validation

In order to validate the proposed linear model estimation strategy, two datasets are produced. In the first dataset, the subject is applying forces on an ATI-Mini 40 forcetorque sensor by only moving his elbow, as depicted in Figure 5.2a. In the second one, both shoulder and elbow movements are combined, as shown in Figure 5.2b. The robotic system is used to guide the movement along the desired direction. The EMG signals, which represent the input of the estimator, are acquired using the Delsys Trigno Wireless System (Delsys Inc) with a sampling rate of 1kHz. The muscular activation signal is derived from the raw EMG data accordingly to the processing steps previously described. The acquisition system applies a first band-pass filter between 20 and 450 Hz, in order to filter out offsets and baseline drifts as well as to avoid aliasing. Then, during the processing, the bandwidth of the signal is further reduced using a Butterworth filter that selects the frequencies between 20 and 350 Hz, implemented by means of the series of a fifth order high-pass and a second order low-pass. Full-Wave Rectification is applied to the resulting signal by taking its absolute value. The envelope extraction, is performed using a second order Butterworth low-pass filter with a cut-off frequency of 1.775 Hz. The parameters of the filtering stages are tuned the same way as in [START_REF] Pesenti | Linear parameter-varying identification of the EMG-force relationship of the human arm[END_REF] in order to be able to compare the results between the proposed method and estimation using all data.

The measured and the estimated force outputs are compared in terms of goodness of fit (FIT) and variance accounted for (VAF), such that

F IT =100 × (1 - ||y -y ˆ|| ||y -y ¯|| ) V AF =100 × σ 2 y -y ˆ σ 2 y
where ||.|| is the Euclidean norm, y ¯is the mean value of y and σ y its variance.

Elbow EMG to force model: In this experimental setup, the model is tested on EMG data acquired on four muscles of the arm, namely biceps brachii, triceps brachii, flexor carpi radialis and brachioradialis. The EMG sensor placement is shown in Figure 5.3a. The aim is to test the recursive modeling of the EMG-force relationship on a single DoF, which is the elbow. During the acquisition, the user was pushing and pulling against the force sensor, as depicted in Figure 5.2a. The resulting muscle activation signals and force measured by the sensor (f h ) are shown on Figures 5.3b and 5.3c, respectively. Using the acquired EMG and force data, the recursive estimation was performed and then tested. As it can be observed in Figure 5.3c, two phases can be distinguished. In the first phase (white background on the figure), the model recursively estimates the model. As the force is measured, the estimate perfectly fits the output. In the second phase (red background on the figure), the estimation is stopped and the model is used to predict the force output. These two phases were repeated in order to show that the estimation can be rerun, if output data are again available. In order to assess the performance of the estimation, the FIT and VAF values where computed on the predictions with respect to the measured force in phase two. It results in FIT = 83.54% and VAF = 97.72%. This performance of the recursive estimated model can then be compared to a linear model estimated offline with all the data as it was done in [START_REF] Pesenti | Linear parameter-varying identification of the EMG-force relationship of the human arm[END_REF], resulting in a performance of the offline model given by FIT = 83.27% and VAF = 97.41%. This shows that the recursive estimation not only gives the possibility to perform online estimation, but also gives a comparable prediction of the output. These results show that for a single DoF, a recursive linear model estimation is accurate enough to describe the EMG-force relationship and thus the feasibility of intention prediction.

Shoulder EMG to force model: In this second experimentation, the model is tested on EMG data acquired on eight muscles of the arm and shoulder: biceps brachii, triceps brachii long and lateral head, anterior and posterior deltoid, pectoralis major, infraspinatus and trapezius. The EMG sensor placement is shown in Figure 5.4a. The sensors where placed accordingly to the SENIAM recommendations [START_REF] Hermens | Development of recommendations for semg sensors and sensor placement procedures[END_REF]. Here, the aim is to test the recursive modeling of the EMG-force relationship in a more complex movement that combines both shoulder and elbow movements. During the acquisition, the user was again pushing and pulling against the force sensor, but only in the frontal plane (vertical plane dividing the body into ventral and dorsal sections), as depicted in Figure 5.2b. The resulting muscle activation signals and force measured by the sensor f h are shown on Figures 5.4b and 5.4c, respectively. As in the elbow case, the two phases can be distinguished in Figure 5.4c. At first the model is estimated and then in the second (red) phase it is used to predict the produced force. The performance was evaluated at FIT = 71.67% and VAF = 93.56%. Even though the performance is significantly lower than in the previous validation, it is still sufficient for assessing the user intention. In fact, if accounting for the much higher complexity of the movement and the greater number of acquired signals, the performance is surprisingly correct. By analyzing more in detail the estimated force, one can notice the overall shape of the force is correctly reproduced. The model, however, seems at some points to struggle with the maximum amplitude of the output. One of the possible reasons could be that during the movement the muscles in the shoulder beside exerting the force measured by the force sensor also compensate for the weight of the arm, what results in parasite muscle activation signals. In fact, one can notice that the maximum amplitude estimation error observed in Figure 5.4c at about t = 32s corresponds to the spike in the muscle activity recorded by EMG 4 in Figure 5.4b, which corresponds to the anterior deltoid, a muscle highly involved in weight compensation. In order to have a better estimation of the force amplitude, more research is required to understand the underlying physiological processes of shoulder compensation and maybe a more specific sensor placement. However, the resulting prediction can be sufficient for many applications that do not require exact values but rather tendencies to assess the user intentions.

Interaction scenario with multiple force sources

In many robotic applications with human-robot interactions, the human and the robot interact not only with each other but also with some external and unmodeled environment.

In this case, distinguishing between the forces that are applied by the operator from these applied by the environment and adapting the robot behavior accordingly becomes essential.

In this context, the following scenario is considered: when no human input is detected by the system, the robot follows some reference position p r under IC and thus can react to environmental contact wrench f e , as show in Figure 5.5a. In contrary, when the operator participates to the interaction by generating some wrench f h on the system, the impedance model used by the controller is modified as shown Figure 5.5b. The desired behavior of system in contact with the environment and without any human input, is given by the impedance model characterized by its stiffness K e and damping D e . In contrary, when only the human operator is exerting wrench on the end-effector, the impedance model is given by the stiffness K h and damping D h . When human and environment interact simultaneously with the robot, the impedance parameters change accordingly to the ratio of both wrenches. The change in the system impedance is considered to be linear and vary between the stiffness (respectively damping) values K e and K h (respectively D e and D h ) such that

K(t) = K h + Γ(t)δK (5.8) D(t) = D h + Γ(t)δD
(5.9)

with δK = K e -K h , δD = D e -D h and Γ is the diagonal switching variable matrix such that for each eigenvalue of Γ, 0 ≤ γ i ≤ 1 and

Γ = g(f h , f e ) (5.10)
where g is some function of the forces generated by the environment f e and forces generated by the user arm f h . In the case when the environment wrench f e is zero, the switching variable matrix Γ is set to be the null matrix.

In practice, the robotic system is not able to differentiate the wrench applied by the environment from the one applied by the operator. It can only measure the total external wrench applied on the end-effector being the sum of all applied wrenches such that f ext = f h + f e . For this reason, in order to assess the wrench f h , the operator is equipped with EMG sensors that can measure his muscular activity. During interaction only between the robot and the operator, the human arm dynamics produce EMG signals µ and a wrench applied on the robot end-effector f ext = f h . These signals are used by a parameter estimator block that recursively estimates the EMG to force model of the human arm. The human EMG to force relationship is given by equation (5.7) that produces estimates of the human wrench f ˆh. When contact with the environment is detected by means of a dedicated sensor, the estimation is interrupted and the previously estimated model parameters are used in the controller to compute the human wrench f ˆh from the EMG signal µ. The overall control architecture is structured as shown in Figure 5.6a. Inside the controller block, as shown in Figure 5.6b, the estimated human wrench f ˆh is used to compute the impedance model switching variable such that Γ ¯= g(f ˆh, f ext )

(5.11)

In this case, the wrench applied by the environment is assessed as f e = f ext -f ˆh. This wrench also contains modeling errors of the estimation of f h . As applying Γ ¯directly to VIC parameters K(t), D(t) could lead to the loss of passivity and by such stability issues, the switching variable is first filtered by means of the passivity filter (4.41), presented in Chapter 4, resulting in Γ that guarantees passivity. This new switching variable is then used in VIC to compute K(t) and D(t) from equations (5.8) and (5.9). Finally, the controller produces the input u that is fed to the robotic system.

Experimental validation

In order to validate the proposed multi force source interaction strategy, an experiment is designed as shown in Figure 5.7. A handle is attached by means of an ATI Mini40 force-torque sensor to the collaborative KUKA iiwa 14 robot. A stretch band of about 180N/m stiffness can be attached to the handle. It permits to emulate an unmodled external forces applied on the system. The operator is equipped with 8 EMG sensors on his shoulder. The placement is done accordingly to Figure 5.4a. The EMG processing is done the same way as described previously in Section 5.2.1.2. The robot is controlled using the KUKA FRI communication interface. For simplicity, only one Cartesian DoF is considered and chosen horizontal and normal to the human and robot, that is in the and Γ ¯= 1 if f ext = f ˆh = 0, with α = 2 a gain that increases the influence of the operator input and also reduces the sensitivity to the estimation error. In fact, as α is chosen high, the weight of the estimation error |f ext -f ˆh| becomes negligible. Equation (5.12) gives Γ ¯= 0 if only the operator force is applied, resulting in an impedance profile K = K h and D = D h from equations (5.8) and (5.9). In contrary, if only the environment exerts a force on the system, Γ ¯= 1 and K = K e and D = D e . Finally, when forces come from both sources, the resulting switching variable is a ratio of the force amplitudes that can be modulated by means of the gain α, as needed for the application.

The experimentation is conducted in four phases. The phases are illustrated in Figure 5.8 with the colored areas. First, the users applies some forces on the handle so that the model estimation can be performed (first white area on figures). In the second phase, the user continues to push-pull the handle, this time however, the force is predicted by the model (red area on figures). Meanwhile, the stretch band is not attached to the handle, so that the only forces measured by the force sensor are these applied by the operator. During these two phases the VIC is given a constant low stiffness of K h = 100N/m and a damping D h = 63.25Ns/m resulting in a damping ratio of 1. In the third phase (second white area on figures), the user force prediction is disconnected and the stiffness is set to K e = 500N/m and the damping to D e = 141.42Ns/m. In this phase the stretch band is attached to the system. Finally, in the fourth phase (green area of figures), the proposed multi force source interaction strategy is launched. The different forces applied on the system come from the stretch band, the operator right arm with EMG sensors and his left arm without sensors. The force sources and the resulting measured and estimated forces are shown in Figure 5.9.

The first remark that can be made concerns the second (red) phase where the estimated model is used to predict the force (Figure 5.8a). As shown in the previous section, the predicted output fits the force intention of the user, even though the estimated amplitude is not always exact. One can also notice that in the third phase, the attached band creates a constant force of -20N on the handle that simulates an external unmodeled interaction force with an unmodeled environment, which can be observed in Figure 5.9 in the first image. In order for the environmental force not to be only static, in phase four (green), the user also applies some additional pulls on the stretch band with his left arm that has no EMG sensors placed on it. This can be observed in Figure 5.9 in the upper right image. Finally, when the operator uses his right arm with attached EMG sensors, as shown in Figure 5.9 in the lower image, the controllers estimates the intentions of moving the handle what results in a decrease of stiffness and damping and a motion of the handle (Figure 5.8c). The decrease of the impedance profiles can be observed in Figure 5.8b. When the user applies some force with the EMG equipped arm, the impedance profiles decrease to K h and D h , whereas in the other cases the profiles remain close to K e and D e . The fact that the profiles are not strictly equal to K e and D e results from the fact that even at rest, the estimation produces some force output due to the gravity compensation activity of the human arm. This difference remains, however, rather low and can further be decreased by using the gain α if more precise tuning of the profiles is required for a certain application. In fact, if α is chosen smaller, it decreases the sensitiv- ity of the controller to small force inputs made by the user. Finally, one can notice that the effective stiffness and damping profiles that are applied by the robot, in Figure 5.8b, are modified by the passivity filter in order to guarantee a stable and thus safe interaction.

The experimentation shows that using the proposed strategy, it is possible to use information about the user intention based on estimating his force with EMG signals and to tune VIC that changes the robot impedance profiles accordingly. It also shows that a fairly accurate linear model of the force produced by the human arm based on shoulder muscle activity assessment can be obtained recursively, which is of particular interest as it is user independent and can deal with the variability of the EMG-force relationship. Also, the implemented impedance profile changes are only influenced by the forces applied by the EMG equipped arm even in case of variations of the environmental forces. In order to face the requirements of different applications the proposed strategy could be used for, the switching variable Γ ¯can freely be adapted as needed. The use of the passivity filter guarantees stable and thus safe execution of the impedance changes. Among the practical applications that could benefit of such an interaction strategy, collaborative polishing tasks could be of interest, as they require to ensure constant contact with the manipulated object and the capacity of interaction with an operator. An alternative use, that will be explored in the further development of this thesis, lies in the field of functional rehabilitation for post-stroke disabled patients.

Chapter 6

Medical applications

With the gradual shift of robots from closed industrial environments into more open and crowded human populated environments, medical procedures featuring robotic solutions attract more and more attention. Current trends in medical robotics aim at including cobots into clinical practice, especially in the fields of functional rehabilitation and robotassisted medical interventions.

In this chapter, practical implementations of the proposed control methods for collaborative robots are presented in the medical context. The first application is in the field of robot-assisted functional rehabilitation of impaired upper limbs in post-stroke patients. The proposed application is a proof of concept in the use of cobots for bimanual rehabilitation. The second application is in the context of practitioner assistance for treatments that require percutaneous needle insertion. A proof of concept of an autonomous robotic needle insertion strategy is proposed, featuring a collaborative robot.

Towards robot-assisted bimanual rehabilitation

In recent years, a growing interest in robot-assisted therapy for the rehabilitation of neuro-motor function disorders has been observed. Robotic rehabilitation devices offer a promising potential to assist the patients and the practitioners during rehabilitation training. In this context, the use of cobotic solutions has come to the fore, due to their attractive cost, their versatility and the ability to work in proximity with humans.

In this section, after giving some insight into neuro-motor function disorders and introducing the necessity of rehabilitation treatment, the use of robotic solutions for rehabilitation training is explored. The main focus is then shifted to a type of rehabilitation training of upper limbs that stimulates both healthy and impaired limbs simultaneously. Finally, a novel proof-of-concept rehabilitation scenario featuring a collaborative robot is explored and tested experimentally at the end of the section.

6.1.1 Robotic rehabilitation of neuro-motor functions disorder 6.1.1.1 Neuro-motor functions disorder Brain trauma, spinal cord injuries and strokes are typical causes of neuro-motor function disorders resulting from damages to the central nervous system (CNS). In France, 140 000 cases of strokes are reported every year, making it the first cause of impairment in adult patients [INSERM, 2019]. The resulting impairment affects a variety of sensory, cognitive and psychological processes. However, the most systematic impairment affects the motor capacity of the subject resulting in a loss of movement control of the face and limbs on one side of the body, which is referred to as hemiparesis. Common problems caused by the impairment result in a condition called spastic paresis, which is a combination of a neural disorder together with a muscle disorder. Such as described in [Gracies, 2005a, Gracies, 2005b, Keyvani, 2006, Gracies, 2015], the neural disorder leads to the muscle weakness for voluntary movements (paresis) and increased reflexes resulting in the overactivity of antagonist muscles opposing any stretch movements of the limb (spasticity). Additionally, a muscle disorder called spastic myopathy caused by the constant tissue contraction, results in structural changes in the muscle and leads to its physical shortening and loss of extensibility. If untreated, this condition is aggravated by the patient psychological reaction that discourages the use of the disabled limb and encourages the compensation using healthy ones, which leads to a non-use of the impaired limb. This disuse results in a deterioration of the CNS and its motor command abilities, thus reducing the cerebral plasticity that could be used to restore the lost functions. This leads to a vicious circle that results in a permanent loss of the limb motion ability and its immobilization in short position [Gracies, 2005a]. The principle of spastic paresis is schematized in Figure 6.1.

In order to break the circle and restore, even partially, the sensory-motor functions of the patient limb, rehabilitation treatment is performed through intensive physical therapy. The aim of rehabilitation is the attempt of a rapid restoration of deficits caused by the CNS damage and the regain of the patient ability to perform daily living activities. This is done by exploiting the ability of the CNS to reorganize itself (neuroplasticity) [START_REF] Dietz | Restoration of sensorimotor functions after spinal cord injury[END_REF]. However, neuroplasticity is limited and most patients only re-Figure 6.1: Principle of spastic paresis, redrawn from [Gracies, 2005b]. cover 70 -80% of their initial capacities [START_REF] Gassert | Rehabilitation robots for the treatment of sensorimotor deficits: a neurophysiological perspective[END_REF]. The timing in which the rehabilitation is performed is of particular importance as the greatest improvements in restoration of the limb motor control are assessed within the first several weeks after the CNS damage, due to an augmented plasticity of the CNS [Cortes et al., 2017]. Also, the stimulation of the neuroplasticity and thus an optimal recovery for the patient, can be achieved by applying treatment strategies that are based on high-intensity and repetitive task-specific practice [START_REF] Langhorne | Motor recovery after stroke: a systematic review[END_REF], with an active physical and cognitive engagement of the patient [Keyvani, 2006]. Thus, optimal rehabilitation training should involve functional and high intensity exercises that feature daily living activities with an active contribution of the patient.

Rehabilitation robotics for upper limbs

The usage of robotic systems to assist rehabilitation therapies for neuro-motor function disorders is an active research topic in medical robotics [START_REF] Maciejasz | A survey on robotic devices for upper limb rehabilitation[END_REF], Gassert and Dietz, 2018, Weber and Stein, 2018]. Robotic rehabilitation devices show great potential when treating neurological injuries such as strokes or spinal cord injuries, as the treatment involves repetitive, high-effort movements with adaptive difficulty levels to restore basic motor-functions [START_REF] Basteris | Training modalities in robot-mediated upper limb rehabilitation in stroke: a framework for classification based on a systematic review[END_REF]. Even though no study has yet proved the superiority of robot-assisted rehabilitation over classical treatment, clinical trials suggest that it provides at least a similar level of recovery as classical therapist-assisted training [START_REF] Kwakkel | Effects of robot-assisted therapy on upper limb recovery after stroke: A systematic review[END_REF]. As for economical reasons the duration of practitioner-assisted rehabilitation has become shorter and shorter, robot-assisted rehabilitation has the potential to overcome this issue by reducing the therapist supervision and increasing the patient to staff ratio [START_REF] Richards | Driving motor recovery after stroke[END_REF]. Robotic devices show other advantages, as they can easily be used as measuring device to assess both the impairment and the progress made during therapy [Gracies, 2015, Cortes et al., 2017]. They can also be integrated into more complex scenarios featuring video games and virtual reality applications, to motivate the patient even further by stimulating his cognitive engagement to perform the task [START_REF] Gassert | Rehabilitation robots for the treatment of sensorimotor deficits: a neurophysiological perspective[END_REF]].

An extensive survey listing currently available robotic solutions for upper-limb rehabilitation was performed by the authors of [START_REF] Maciejasz | A survey on robotic devices for upper limb rehabilitation[END_REF]. This survey highlights two main categories of devices, which are exoskeleton robots and end-effector systems. Exoskeleton devices are far beyond the scope of this thesis. For this reason, only endeffector systems are considered in the further development of this work. As pointed out in [START_REF] Gassert | Rehabilitation robots for the treatment of sensorimotor deficits: a neurophysiological perspective[END_REF], these systems have an important advantage over exoskeleton robots, as they allow achieving higher motion dynamics together with a wider range of interaction characteristics. In fact, as pointed out in [START_REF] Metzger | Performance comparison of interaction control strategies on a hand rehabilitation robot[END_REF], the ability to adapt the robot output impedance is very important for rehabilitation training in order to limit interference with the movements performed by the patient. The survey [START_REF] Maciejasz | A survey on robotic devices for upper limb rehabilitation[END_REF]] also shows that almost all proposed systems are task specific devices that are custom-designed. This drastically increases their cost, making them less available for everyday therapy. Among well known commercial end-effector systems, two most popular examples are the InMotionARM from Bionik, Canada (Figure 6.2a) and the ArmeoPower from Hocoma, Switzerland (Figure 6.2b). More recently, motivated by the economical limitations of robot-assisted rehabilitation, some authors explored the usage of industrial cobots for therapy, which appear as an interesting alternative to high-cost dedicated robotic solutions. In their work, the authors of [START_REF] Kyrkjebo | Feasibility of the ur5 industrial robot for robotic rehabilitation of the upper limbs after stroke[END_REF] suggest the feasibility of using a standard collaborative robot for rehabilitation, in terms of safety requirements and control. Also, even more recently, a commercial rehabilitation device based on a collaborative KUKA iiwa robot was brought to the market by Life Science Robotics, Denmark. This robot, Robert (Figure 6.3), provides support in the rehabilitation of immobilized patients and currently undergoes clinical studies to evaluate its efficiency [START_REF] Bertelsen | Use of robot technology in passive mobilization of acute hospitalized geriatric medicine patients: a pilot test and feasibility study[END_REF]. These examples clearly show the high potential and the growing interest in the use of collaborative robots in the context of robot-assisted rehabilitation.

Bimanual rehabilitation

The complexity and dexterity of the human arms makes them particularly impacted by CNS damage. Arm impairments greatly limit the patient ability to perform daily living tasks and especially unimanual and bimanual reaching and grasping, thus requiring a particular attentions in research [START_REF] Gassert | Rehabilitation robots for the treatment of sensorimotor deficits: a neurophysiological perspective[END_REF]. In contrary to classical rehabilitation strategies that focus only on the impaired limb, bimanual rehabilitation consists in using both healthy and impaired arm simultaneously, in a common bimanual motion. Such a strategy is of particular interest, as many upper limb motions performed in daily living tasks require the use of both hands simultaneously. This kind of rehabilitation relies on the proprioceptive feedback from each limb and can help to retrain the motor pathways of the impaired limb [Wolf et al., 1989, Rose andWinstein, 2004]. Even though there is no clear clinical evidence of the superiority of bimanual rehabilitation over other strategies, some clinical studies suggest its efficiency in restoring some types movements such as reaching and grasping tasks [START_REF] Mudie | Can simultaneous bilateral movement involve the undamaged hemisphere in reconstruction of neural networks damaged by stroke?[END_REF]. As pointed out in [START_REF] Gassert | [END_REF]Dietz, 2018, Mutalib et al., 2019], bimanual rehabilitation could facilitate the communication between both affected and healthy hemispheres when dealing with tasks that require the action of the two limbs performing coordinated force inputs. This could lead to an active support of the healthy hemisphere in the movement of the impaired arm [START_REF] Schrafl-Altermatt | Cooperative hand movements in post-stroke subjects: Neural reorganization[END_REF]. The authors of [START_REF] Gassert | Rehabilitation robots for the treatment of sensorimotor deficits: a neurophysiological perspective[END_REF] point out the advantages of robot assistance for bimanual rehabilitation therapy and believe that its potential is currently not sufficiently explored, including the understanding of the underlying processes.

In fact, in the field of robot-assisted bimanual rehabilitation, only a few solutions have been proposed. In their pioneer study, the authors of [START_REF] Lum | Robotic assist devices for bimanual physical therapy: preliminary experiments[END_REF]] designed a robotic system to help the patient with transportation and squeezing tasks. In [START_REF] Gassert | Rehabilitation robots for the treatment of sensorimotor deficits: a neurophysiological perspective[END_REF] it is pointed out that combining two unimanual devices in a mirrored configuration could be sufficient for bimanual rehabilitation. One of the main arising challenges is the choice of the physical connection between the two hands of the patient. If the coupling stiffness is too rigid, the patient is more likely to apply minimal contribution with the impaired arm, as the healthy one dominates the motion. In contrary, if the connection is too soft, severely impaired patients would not be able to perform any bimanual motion as the strength of the weaker limb is not sufficient [START_REF] Rose | Bimanual training after stroke: are two hands better than one?[END_REF]. Some work has been done to deal with this issue, as in [START_REF] Mcamis | Design and analysis of a compliant bimanual rehabilitation device[END_REF] where the authors design a passive and compliant system for bimanual rehabilitation. This solution results in a cost effective, however task specific, mechanical system operating in a 2D plane.

Cobot-assisted bimanual rehabilitation

The previous literature study highlights several important points. First, it shows the high potential and the growing interest in the use of cobots for robot-assisted rehabilitation, due to their attractive price, their versatility and the capacity to work closely with humans. Also, it points out that even if not much explored yet, bimanual rehabilitation therapy could offer some promising advantages over classical therapy for upper limb rehabilitation, especially in terms of increased communication between the affected and healthy hemispheres. However, for such therapy to be efficient for the patient, the robotic system needs to allow smooth and natural motions. It also requires the ability to adapt its impedance. Also, the question of the proper coupling of the two arms performing the motion is to be considered.

In order to address the highlighted points, this section proposes a proof-of-concept application of cobot-assisted bimanual rehabilitation. This application features a cobot that is used to stimulate the communication between the two hemispheres of the impaired patient. This is done by enforcing the coordination between the two limbs in a daily living application. This is of particular importance in tasks that require force symmetry such as object transportation. In fact, as pointed out in [START_REF] Mutalib | Bimanual coordination during a physically coupled task in unilateral spastic cerebral palsy children[END_REF], the lack of force coordination between the healthy and impaired limb introduces a force bias towards the healthy side.

In order to overcome this issue, a system has been designed with two handles seized by the patient. The system has 1 DoF in translation to emulate an object transportation task, which is a very common daily living activity. The proposed training scheme is illustrated in Figure 6.4. During the performed task, the patient arms generate possibly Figure 6.4: Concept for bimanual rehabilitation in a 1DoF task. different forces, denoted f s for the strong and healthy arm and f w for the weaker, impaired arm. The two handles are rigidly coupled to each-other and to the robotic system. The use of the cobot has the advantage of minimally interfering with the patient motion, due to the robot good transparency. To address the challenge of the coupling between the two arms of the patient, the robotic system impedance is controlled using VIC. The main issue in the coupling between the two upper limbs is related to the possible overcompensation of the impairment by the healthy limb. For this reason, in order to increase the two limbs synergy and in the same time avoid overcompensation, the robot impedance is controlled in such a way that it can only be moved when good coordination is achieved. That is, if both arms apply equivalent forces, the apparent stiffness of the robot is low and the transportation task can easily be achieved. In contrary, when a force bias is introduced, a high stiffness makes the achievement of the task hard. Such a behavior can be obtained by defining a switching function Γ for VIC that equals 0 to encourage and 1 to penalize movement. This can be achieved by defining Γ such that enforces coordination by imposing equal force input from both arms. In fact, if only one of the arms applies a force, the switching function Γ tends towards 1. In contrary, if both arms perform the movement with the same force input, Γ tends towards 0. As it is not straightforward to apply exactly the same amplitude of force with both limbs, the switching function allows some force bias with an increasing tolerance with growing values of force. The switching function is then used to modify the system impedance parameters, such that

Γ = |f w -f s | |f w | + |f s | (6.
K(t) = K 0 + Γ(f s , f w )δK (6.2) D(t) = D 0 + Γ(f s , f w )δD (6.3) with δK = K 1 -K 0 , δD = D 1 -D 0 .
The values of K 0 and D 0 are chosen low to allow the movement of the system when both forces are equal. In contrary, the values of K 1 and D 1 are chosen high and penalize movement if the applied forces are not coordinated. Also, in this case it is considered that the forces f s and f w can be measured separately, which can either be done by using two force sensor at each handle, or by using force information from EMG, as presented in Chapter 5. As it is difficult to assess activity on the impaired arm, a more practical solution consists in using EMG on the healthy arm. This allows to estimate a correct model of the EMG to force relationship, whereas it is probably not possible to do the same because of the muscle overactivity in the impaired arm that could interference with the EMG signals.

Though the use of variable impedance is required for the application, it is not the only proposition of this thesis that can be used in this context. Indeed, when working with impaired patients, it is also important to add additional limitations to the system to prevent any harm to the already impaired limb. Form this reason, the MPIC strategy presented in Chapter 3 is of interest. It can be used to impose constraints, for instance, on position to limit the workspace and decrease the risk of damage to the patient's joints.

Experimental results

In order to test the proposed bimanual scenario for functional rehabilitation, the experimental setup shown in Figure 6.6 is proposed. In this setup, a double handle is attached Figure 6.6: Experimental setup for bimanual rehabilitation. The weak (red) and strong (green) arms produce a weak f w and a strong f s force, respectively. The EMG sensors are placed on the strong arm.

to an ATI Mini40 force-torque sensor that is fixed to the flange of the KUKA iiwa robot. This force sensor measures all the forces applied to the handle. On the patient side, EMG sensors are placed on his healthy arm. In this experiment, the EMG sensors are placed on the shoulder of the right arm. The EMG placement and processing is done following the same protocol as in Chapter 5, according to Figure 5.4a. The robot is controlled using MPIC that allows to impose compliant motion in one main direction and set workspace limits, whereas the variation of the impedance parameters is ensured to result in a passive interaction by the passivity filter. All the other DoFs are controlled using a classical IC with high stiffness. The variation of the impedance parameters is defined by the low values of stiffness K 0 = 100N/m and damping D 0 = 63.25Ns/m, enabling easy motion, and the high values of stiffness K 1 = 500N/m with damping D 1 = 141.42Ns/m, which severely penalize motion. The workspace is limited to ±0.2m with respect to the central position (note that in the obtained results this limit is never reached).

As only information about the total external force f ext and EMG-estimated strong arm force f ˆs are available, the force f w applied by the weak arm results from the difference between these two forces, i.e. f w = f ext -f ˆs. The switching variable Γ, is computed from the reference function

Γ ¯= |(f ext -f ˆs) -f ˆs| |f ext -f ˆs| + |f ˆs| (6.4)
equivalent to equation (6.1), and checked by the passivity filter (4.41), presented in Chapter 4. In the same way as presented in the experiment from Chapter 5, in the first stage of this experimentation, the model of the user EMG to force relationship is obtained recursively with a performance evaluated at FIT = 71.00% and VAF = 91.97%.

The resulting experimental data is plotted in Figure 6.7. The green areas correspond to the phases of the experiment where the user applies equivalent forces with both hands. Alternatively, the white areas correspond to the phases where the user applies a greater force with the stronger hand, imitating the situations where the patient overcompensates his disability. It is important to notice that the operator is a healthy subject, thus the emulated impairment does not reflect clinical practice. The first plot (Figure 6.7a) shows both the total force exerted on the tool and the estimated strong arm force. One can notice that in the case where both hands apply a similar force (green areas), the amplitude of the estimated strong arm force f ˆs can approximately be related to the half of the measured total force f ext . In this case, one can see in Figure 6.7b that the resulting reference switching function Γ ¯, and its passivity guaranteeing counterpart Γ, result in an important decrease of the impedance parameters K ¯, D ¯and K,D for Γ ¯and Γ, respectively. Although the decrease is significant, the switching function and thus the impedance parameters never reach their minimal value, as it is very difficult to produce exactly equal forces with both arms. This modification of the interaction dynamics allows to facilitate the motion of the handle, resulting in a greater displacement, as it can be seen in Figure 6.7c. In contrary, during the phases in which the stronger arm produces more important forces that the impaired one (white areas in Figure 6.7), the decrease of the impedance parameters is much less important, making the motion of the system more difficult. This can be seen when comparing the different amplitudes of the displacement of the handles in Figure 6.7c with the similar amplitude of applied forces in Figure 6.7a. An important point to notice, which is directly related to the choice of the switching function lies in the sensitivity of the difference between the arm inputs. In fact, one can clearly see that the smaller the forces involved in the task, the faster the switching function saturates the impedance parameters to their maximum value. In contrary, when the involved force amplitudes are higher, even unbalanced forces can modify the impedance parameters. This can lead to a limitation of the proposed switching function, as estimation errors are interpreted as a force produced by the weak arm. Such a case can be observed at time t = 28s, where the force applied on the handle only comes from the strong arm with a high amplitude of around 25N (Figure 6.7a). The rather small estimation error is interpreted as the input of the weak arm and thus allows the variation of the impedance parameters (Figure 6.7b). The proposed switching function is therefore very well suited for an application in which the forces that are to be applied by the patient are desired to have an amplitude between 10 and 20N.

The obtained results demonstrate that it is feasible to perform a bimanual rehabilitation task that stimulates the coordination between a healthy and an impaired limb using the robot control tools developed in this work. In this common daily living task, which is object transportation, the robotic system allows to modulate the coupling between both arms of the patient. It also shows that the choice of the switching function is particularly important accordingly to the desired performance and to the range of the forces that are to be applied. In order to make the task more immersive and more intellectually stimulating, which is an important factor for the efficiency of the recovery process, a scenario featuring virtual reality could easily be added to make the task more realistic.

Towards autonomous robotic needle insertion

In recent years, a generalization of minimally invasive medical procedures has occurred, because of their drastic impact on the patient well-being. In fact, this type of procedures allows to significantly reduce the hospital stay but also the discomfort for the patient by reducing pain, scars and the risk of post-procedure infections and bleeding. Minimally invasive procedures can be performed by accessing areas of interest by different means, for instance through natural orifices or directly by puncturing the skin with a needle. This last strategy called percutaneous is commonly used in interventional radiology, a medical discipline combining specialized medical know-how in radiology with latest technological innovations in the fields of medical imaging and possibly robotics [START_REF] Esr | Interventional radiology in European radiology departments: a joint survey from the European Society of Radiology (ESR) and the Cardiovascular and Interventional Radiological Society of Europe (CIRSE)[END_REF].

In this section, after a short introduction into the field of interventional radiology and its use of robotic solutions, the main focus is shifted towards autonomous robotic needle insertion. At the end of the section, an experimental proof of concept of autonomous needle insertion into an anatomical phantom is presented, featuring automatic layer detection and variable impedance control of the needle. Even though autonomous applications in robotic needle insertion tasks are currently not used in clinical practice, they represents a major milestone for medical applications in the future.

Percutaneous procedures in interventional radiology

Biopsies and interventional radiology

According to the World Health Organization (WHO), cancer is the second cause of mortality in the world, responsible for an estimated 9.6 million deaths in 2018 [World Health Organization, 2018]. Nowadays, many treatments are available to cure cancer. However, their efficiency highly depends on the early detection of the pathology. Biopsy is one of the most accurate techniques in cancer diagnosis, as it provides a precise knowledge about the tumor structure and allows the prescription of the most adapted treatment. Biopsies are not only a one-time tumor detection procedure but allow to follow up its development [START_REF] Dalag | Lung and abdominal biopsies in the age of precision medicine[END_REF] in response to standard treatment, or in the context of clinical trial. This is likely to lead to an important increase in performed biopsies in the future [Cherukuri et al., 2019]. Among the different biopsy techniques, the most interesting for the patient are these performed percutaneously in interventional radiology, as they are minimally invasive. The practitioner performs the procedure by accessing the area of interest with a needle shaped tool under the guidance of a medical imaging device. Depending on the indication, the imaging modalities used for such procedures are mainly ultrasound echography, magnetic resonance imaging (MRI) and X-ray imaging. This latter is the most popular because of its good trade-off between cost and image resolution [START_REF] Veltri | Cirse guidelines on percutaneous needle biopsy (pnb)[END_REF]. This imaging technology has however some major drawbacks such as the long term X-ray exposure of the medical staff. Beyond biopsies, interventional radiology is also very efficient in pathology treatment as it allows the percutaneous destruction of tumors by thermal-ablation or radio-frequencies.

In interventional radiology the available needle shaped tools differ not only in their diameter but also in the shape of their extremity that can be either symmetric or asymmetric, what influences the behavior inside the tissue [START_REF] Abolhassani | Needle insertion into soft tissue: A survey[END_REF]. In the following development, only symmetric needles are considered.

Robot-assisted needle insertion

Because of the precision requirements for needle insertion tasks, and even more because of the impact of radiations on the practitioner, much research effort has been done and continues to assist interventional radiology procedures. In fact, because of the harmful radiations, the practitioners are forced only to use images when strictly necessary and have to reconstruct the needle insertion trajectory only based on their experience, on the haptic feedback of the needle and on some sampled position data. The requirement for protecting the practitioner and for increasing his dexterity in the case of reduced access to the area of interest make the usage of robots particularly suitable. In recent years, many robotic systems have been designed purposely, as identified in [START_REF] Kulkarni | Review of robotic needle guide systems for percutaneous intervention[END_REF]. In particular, much research effort has been directed towards the development of robotic structures that are compatible with the medical imaging technologies. Robotassisted needle insertion also allows to perform automatic or semi-automatic insertion and path planning based on sampled information acquired by the imaging device and sensors mounted on the robot. In the future, this could allow less experienced practitioners to safely perform biopsy procedures.

Even though many robotic systems are still being developed, only a few ones are currently available commercially. This shows that the current technology is still not mature enough for everyday procedures. This is confirmed by the fact that in its 2017 report on recommendations for biopsy procedures, the Cardiovascular and Interventional Radiological Society of Europe (CIRSE) does not mention any robotic solution [START_REF] Veltri | Cirse guidelines on percutaneous needle biopsy (pnb)[END_REF]. The major reason is due to the current lack of regulatory approval for autonomous systems to enter into clinical practice, in the case where the system does not replicate the practitioner manual insertion [START_REF] Kulkarni | Review of robotic needle guide systems for percutaneous intervention[END_REF]. Therefore, a great research effort is oriented towards teleoperated robotic mater-slave systems with haptic force feedback [START_REF] Abolhassani | Needle insertion into soft tissue: A survey[END_REF]. These systems allow to protect the practitioner from radiations and provide some level of assistance while giving him a similar degree of control over the insertion as in the manual procedure. However, with the fast evolution of robotic research and the growing acceptance of robots in everyday tasks, the exploration of alternative, even more ambitious solutions, becomes important. The ultimate challenge in percutaneous procedures would be to enable robots to safely perform fully autonomous needle insertion. In this perspective it is important to understand why manual insertion is currently more reliable for these types of application. In particular, the understanding of the interactions between biological tissues and the surgical needle is of great interest. that in phase 0 the needle is in contact with the tissue surface, in phase 1 it pushes against the tissue surface and deforms it. Then, when the force applied to the tissue increases over a certain value, the needle punctures the tissue. Inside the tissue, friction forces apply on the needle shaft in addition to cutting forces. Finally, in phase 3, when the needle is extracted from the tissue, friction makes the tissue stick to the needle. During phase 1, the interaction between the needle and the tissue can be approximated by means of a viscoelastic model, as long as small motions are considered [Fung, 1993]. This allows to describe the relationship between forces and displacements of the interacting bodies. Viscoelastic models can be represented by different combinations of springs and dampers such as the Maxwell model where both elements are serial or the Kelvin-Voigt model with a parallel structure, to mention just a few [START_REF] Yang | Force modeling, identification, and feedback control of robot-assisted needle insertion: A survey of the literature[END_REF]. The Kelvin-Voigt model, shown in Figure 6.9, is the most current representation of interaction with soft tissues. The interaction between the needle and the tissue can be expressed as Figure 6.9: Kelvin-Voigt model representation to describe viscoelastic interaction between a needle and a soft tissue before puncture.

f (t) =
k(p(t) -p(0)) + dp ˙(t), if p(t) > 0 0, if p(t) ≤ 0 (6.5)

where f (t) is the force exerted by the needle on the tissue, p(t) and p ˙(t) are the position and velocity of the needle tip and k, d are the stiffness and damping coefficients of the tissue model, respectively. In such a model, and in many references in the literature, the stiffness and damping are considered constant for a given tissue, which does not allow to take into account nonlinearities in the reaction force of the tissue due to membrane punctures or to the presence of different tissue layers. To overcome this issue without complexifying the model too much, [START_REF] Barbe | In vivo model estimation and haptic characterization of needle insertions[END_REF] propose to consider time dependent stiffness and damping coefficients, resulting in a model given by

f (t) = k(t)(p(t) -p(0)) + d(t)p ˙(t), if p(t) > 0 0, if p(t) ≤ 0 (6.6)
Even if the model of equation (6.6) does not have physical accuracy during the whole needle insertion, it has proven to be efficient in estimating the force exerted by the needle. Interestingly, this allows detecting events that vary during the needle insertion such as punctures [START_REF] Barbe | In vivo model estimation and haptic characterization of needle insertions[END_REF].

Autonomous insertion by imitation

During percutaneous insertion, the needle passes through different layers of biological tissue before reaching the area of interest. For instance, in the case where the target area is placed inside the liver, the needle will have to pass through four layers: skin, fat, muscle and the liver surrounded by a fibrous capsule. These biological tissues present different characteristics [START_REF] Barbe | In vivo model estimation and haptic characterization of needle insertions[END_REF]. During a manual insertion, the practitioner has to adapt the stiffness of his arm in order to control the insertion and avoid involuntary movements that may occur for instance after a large deformation followed by a puncture.

When the tissue exhibits low resistance, resulting form low stiffness and damping parameters or low friction, the practitioner inserts the needle with a low arm stiffness. Doing this allows better controlling the insertion velocity and the exerted efforts, which improves the safety of the procedure. At the contrary, when the needle enters in contact with a more rigid structure, which can be either the next layer or some unexpected inclusion, such as a vessel, the practitioner perceives a force feedback. This information is used by the practitioner to react accordingly and stiffen his arm to continue the puncture. Such an interaction scheme allows to better deal with unpredicted contacts and react accordingly.

In order to propose a similar strategy in an autonomous needle insertion and deal with stability of the interaction a variable impedance controller can be used. A few comparable applications featuring VIC have been proposed to date, such as in [Cho et al., 2015] where a VIC is used in a needle insertion task to compensate for respiratory motion. In the considered application, the controller only varies the damping profile in order to prevent high forces exerted by the needle on the tissue in the directions normal to the insertion, in order to avoid needle bending or laceration of organs. Even though it reduces the forces applied to the tissue, this strategy is not suited when the forces need to be increased, as it is the case when a more rigid tissue needs to be punctured. In fact, in this case, a critical issue consists in guaranteeing stability of the system while modifying the impedance profiles. In [START_REF] Ferraguti | A tank-based approach to impedance control with variable stiffness[END_REF], an energy tank is used to control passivity during impedance changes occurring in needle insertion tasks.

Alternatively, as the soft tissue is a passive system, the overall passivity can be enforced by means of the passivity filter proposed in Chapter 4. The goal of the following development is therefore to propose a fully autonomous robotic needle insertion application. The insertion will be performed into a succession of biological tissues using a variable impedance controller and a passivity filter allowing to modify the desired impedance profiles. An interesting advantage of the method is that it requires no knowledge about the layers position, as no a priori information on the depth of the tissue layers is available.

Automatic layer detection

To overcome the issue of undetermined depth of the different layers, the strategy of detecting layer transitions from [START_REF] Barbe | In vivo model estimation and haptic characterization of needle insertions[END_REF]] can be employed. It uses a recursive estimation of the tissue stiffness and damping parameters [START_REF] Barbe | Online robust model estimation and haptic clues detection during in vivo needle insertions[END_REF] and the prediction of the interaction force between the needle and the tissue. In the case, when the needle punctures a new layer, the measured force does not fit the estimated model anymore, resulting in an increase of the estimation error. Its detection allows determining a change in the structure of the biological tissues. The method of [START_REF] Barbe | In vivo model estimation and haptic characterization of needle insertions[END_REF] is further described in the following paragraphs.

Recursive model parameter estimation :

The recursive parametric estimation algorithm can be formulated by expressing the model (6.6) in its discrete form

y k = Φ T k Θ k-1 + w k (6.7) where Θ k = [k k d k ]
T is the vector of model parameters at time step k and Φ T k = [p k -p 0 p ˙k] is the regression vector with p k and p ˙k the measured needle tip position and velocity respectively, and p 0 the initial position of the tissue layer surface. The measured force is given by y k = f k with w k the zero mean value measurement noise. The estimated force y ˆk = f ˆk can then be expressed after the estimation of the model parameters Θ ˆk, as

y ˆk = Φ T k Θ ˆk-1 (6.8)
In order to estimate the parameter vector Θ ˆk, a classical solution consists in using the recursive least square (RLS) algorithm to identify the model parameters resulting in the following least square regression min

Θ ˆ||y k -Φ T k Θ ˆk|| 2 2 (6.9)
The minimization problem (6.9) can be recursively solved using the RLS method [ Åström and Wittenmark, 2013] with a forgetting factor 0 < λ < 1. The estimation parameter update law can be written as Θ ˆk = Θ ˆk-1 + L k e k (6.10) with the estimation error (6.11) and the gain matrix

e k = y k -Φ T k Θ ˆk-1
L k = P k-1 Φ k λ + Φ T k P k-1 Φ k (6.12)
where P k is the covariance matrix of the parameters such that

P k = 1 λ (P k-1 -L k Φ T k P k-1 ) (6.13)
Detecting layers using the estimation error : As pointed out in [START_REF] Barbe | In vivo model estimation and haptic characterization of needle insertions[END_REF], the tissue transition corresponds to large variations of the estimation error which can be easily observed by analyzing the square of the error as detection signal s k = e 2 k . The information about the layer transition can then be assessed by means of the CuSum algorithm [START_REF] Basseville | Detection of Abrupt Changes -Theory and Application[END_REF]. The aim of the CuSum algorithm is to detect the changes in the mean value of s k . A decision function g k is then built such that

g 0 = 0 g k = max(0, g k-1 + s k -ν) (6.14)
where ν is the detection threshold that is generally chosen as half of the amplitude change to be detected [START_REF] Basseville | Detection of Abrupt Changes -Theory and Application[END_REF]. The detection of an event is decided when g k ≥ µ, with µ the decision threshold. In this case, an alarm signal is triggered. If g k continues to evolve freely, variations of the decision function may prevent it from decreasing under the threshold µ making further detection impossible. In this case, it is preferable to either reset g k after each event detection, or to saturate it at some value.

The choice depends on the desired behavior of the detector, as if g k is reset, a single event will trigger many alarms. On the contrary, if the saturation is chosen too far away from the decision threshold µ, rapidly evolving events will not be detected, as the decision function will not have enough time to decrease under the threshold.

Experimental results

A proof of concept of autonomous needle insertion into an anatomical phantom has been developed in order to illustrate the advantages of VIC with passivity guarantees. The aim is to highlight the importance of the control strategy for the task rather than to propose a solution ready for clinical application. The proposed experimental setup features a KUKA iiwa collaborative robot that could be of great benefit in a crowded medical environment, but is not suited for the general requirements of medical imaging. The setup is presented in Figure 6.10.

Figure 6.10: Experimental setup.

In the following experiment, the robot inserts the needle into a silicon phantom that imitates biological tissues. The phantom has been developed according to the procedure described in [Pfeil, 2020]. It consists of four successive layers representing respectively the skin, fat, muscle and liver. The main characteristics of the phantom are listed in The forces applied by the needle to the tissue are measured by means of an ATI-mini 40 force-torque sensor mounted on the robot end-effector. The robot is controlled with a VIC as described in Chapter 4, along the direction of insertion. All the other directions are controlled using basic IC with high stiffness. As the phantom is a succession of rigid and soft layers, the reference stiffness for the controller is selected to vary between 300N/m for the fat and liver and 500N/m for the skin and muscle. The robot damping is determined by its ratio of 1 and varies between 134.2Ns/m and 173.2Ns/m. The robot apparent inertia is set to 15kg. The passivity filter gain is set to β = 50. The reference trajectory is a smooth fifth order polynomial function in position, velocity and acceleration. The position vector is defined as p = [p x p y p z ] T . For this proof of concept the trajectory end-point is an arbitrary point inside the liver layer.

During the experiment, the controller has no input about the depth of the different layers of the phantom. The only information about the phantom that the controller is aware of is the number of layers that need to be detected. The impedance profiles are switched according to the events reported by the detection algorithm. Both detection and decision threshold are set ν = µ = 0.2. The detection is only performed during the insertion phase, as during extraction the model (6.6) is not valid anymore.

The insertion forces, as well as the output of the RLS estimation during the insertion phase are plotted in Figure 6.11a. It can be observed that the estimated forces diverges at some points from the measured one. By analyzing the square of the estimation error, visible in Figure 6.11b, one can notice that it increases at some specific points. By using the CuSum algorithm to detect the rise of the estimation error, the boundaries between the layers can clearly be observed. It is also worth noticing that at the end of the insertion the estimation error varies more often, which is caused by the model error, as the needle slows down towards the target point, what results in a modified influence of the unmodeled friction forces. That is why it is necessary to give the algorithm the number of layers that need to be detected. The layer positions are detected at p = {75, 87, 108, 122}mm for skin, fat, muscle and liver, respectively. The tissue layer thickness can therefore be estimated as 12mm for the skin layer, 21mm for the fat layer and 14mm for the muscle layer, which corresponds to the data given in Table 6.1 with a position estimation error of 2mm for the skin, 1mm for the fat and 1mm for the muscle layer. These results are specially good as the tissue is deforming during insertion.

Using the position data from the layer detector, the VIC is able to modify the robot stiffness and damping profiles accordingly to the output of the passivity filter. The evolution of the switching variable γ ˙as well as the passivity conditions h 2 and h 3 are shown in Figure 6.12. It can be seen that in the case where the impedance increases, the The resulting stiffness and damping profiles for both needle insertion and extraction are plotted in Figure 6.13. These results show that inside the skin and muscle the impedance profiles increase until reaching the desired value. It also shows that in some cases, as it is the case of inside the skin during extraction, the passivity filter does not allow the impedance profile to reach its reference value before the end of the layer. This is caused by a more important needle extraction velocity, resulting in less time for the Figure 6.12: Evolution of passivity conditions and γ ˙. For better readability, condition h 2 was scaled by a factor 10. switching variable to increase to its final value in a passive way. For most applications it is not of crucial importance to reach precisely the targeted impedance. However, if for some reason it is required, slowing down the velocity would be sufficient to overcome this issue. One can also note that outside of the phantom or inside the fat and liver, the profile decreases to its reference as a first order filter.

This experimentation shows that similarly to manual insertion, it is possible to identify the tissue layers by only using forces applied to the needle. It also shows how this information can be used to adapt the interaction dynamics between the needle and the tissue. As the proposed application is a proof of concept, only the simple case with a switch between two impedance values is proposed. However, for more realistic tasks, it is possible to use a similar strategy but to define a specific impedance profile for each layer independently. Alternatively, some more complex scenarios, exceeding the scope of autonomous insertion, could be proposed by including the practitioner into the control loop and allowing him to modify the system behavior to insert the needle. Such an application would result in a collaborative needle insertion task.

Chapter 7 Conclusion

With the growing interest in the use of robotic systems in a shared workspace with human operators, the development of collaborative robots offers new possibilities for humanrobot interactions. Therefore, developing new control tools for interaction management became an important research topic. In this context, the design of control solutions that, at the same time, allow to enhance interaction dynamics and guarantee operator integrity is of particular interest.

Contributions

In this thesis, several control tools for collaborative robotics, which ensure security and enhance the dynamics of physical human-robot interactions (pHRI), are proposed. The problem of simultaneously ensuring robot compliance and managing constraints while increasing the robot reactivity is addressed and leads to the design of model predictive impedance control (MPIC). MPIC is a combination of impedance control (IC) with model predictive control (MPC). The principle of IC is used to manage unexpected interactions between the robot, the operator and the environment. Meanwhile, MPC allows to impose constraints and in the same time uses its capacity to predict the system behavior, hence increasing its reactivity. The proposed methods consists in designing the cost function of MPC such that it enforces the desired compliance of the robot while being able to satisfy practical constraints, such as position, velocity and acceleration limits. Experimental results validate the proposed controller when interacting with an unknown static environment and in the case of human-robot collaboration and illustrate the utility of the proposed control strategy.

Whereas in many applications the desired interaction dynamics can be considered con-stant, during pHRI the ability to adapt the compliance of the system is often required. An important and dangerous issue that results from the modification of the interaction dynamics is the possible loss of passivity that can lead to instability and an unexpected behavior of the robotic system. Therefore, ensuring a passive change of interaction dynamics is an important step towards enhanced pHRI that guarantees operator integrity. In this work, a design method for passivity filters is proposed. These filters guarantee the passivity of the interaction by checking the modification of the interaction dynamics and modifying it if necessary. The first developed filter is based on state-of-the-art passivity conditions which are shown to be very conservative. For this reason, alternative passivity conditions are explored leading to the development of a second passivity filter. The presented simulations show that the second filter allows to guarantee passivity for much faster changes of the interaction dynamics. Experiential results, obtained on a setup that leads to instability without the passivity filter, show the efficiency of the proposed method.

With the possibility of changing the interaction dynamics, the understanding of the operator intentions and including this information into the robot control strategy for enhanced pHRI comes to the fore. A strategy based on electromyography (EMG) is proposed to distinguish operator force inputs from interactions with the environment and modify the robot interaction dynamics accordingly. Because of their fair trade-off between accuracy and complexity, as well as their capacity of adaptation to changing conditions, linear time-varying models that are recursively identified are explored to map the relationship between EMG and the force produced by the human arm. This information is used in the robot control strategy to perform a modification of the interaction dynamics. Experimental validation shows a good performance of the proposed EMGforce models and illustrates the advantages of such a strategy in cases where human, robot and environment interact with each other.

Finally, two proof-of-concept applications featuring the proposed control tools in the medical field are described. The first application is developed in the context of bimanual rehabilitation therapy for patients with neuro-motor function disorders. Its aim is to propose a cobot-assisted rehabilitation scenario that stimulates the coordination of the forces produced by the healthy and the impaired arm of the patient in a bimanual activity. Information about the contribution to the task of both arms, together with the estimation of the force produced by the strong arm equipped with EMG sensors, are used to modify the robot compliance. High compliance of the system resulting in the possibility to perform the transportation is only achievable if both arms contribute equivalently to the task. In other cases, the motion is penalized by a low compliance of the robot. The change in robot compliance is guaranteed passive by means of the passivity filter. In order to safeguard the patient joints form injuries, workspace constraints are ensured by means of MPIC. The proposed application was validated experimentally on a healthy subject emulating impairment.

The second proposed medical application features a cobot performing an autonomous, percutaneous needle insertion task in the context of interventional radiology. The interaction between the needle and the punctured tissues is controlled by adapting the interaction dynamics accordingly to the mechanical properties of the tissue. As the exact depths of the different tissue layers are not known in advance, a layer detection algorithm is used to determine the tissue layer that is currently being punctured. This information is used in order to modify the robotized needle interaction dynamics, which is particularly important as the tissue rigidity varies significantly between the different layers. The experimental insertion is performed into an anatomical phantom emulating human tissues and displays promising results in terms of the precision of the layer detection and the ability to adapt the behavior of the needle.

Limitations and perspectives

Even though the control methods presented in this thesis display interesting properties and have proven their efficiency in the proposed applications, some limitations can be highlighted, which leave room for future improvements.

Starting with a general remark on the hypothesis globally made in this thesis, all of the proposed methods consider the reference impedance model with a constant and imposed virtual inertia of the system. Even though this hypothesis can easily be achieved, as shown in this work, it requires the use of a force-torque sensor measuring the external wrench applied on the robot, which can be an issue for practical implementation. A possible extension of the proposed methods would consist in adapting them to be used with the real, configuration dependent inertia of the robotic system, which could allow the removal of the sensor.

Even though MPIC is in practice stable when entering and leaving the constraints, no mathematical stability analysis has been performed yet. One possible extension would consist in modifying the MPIC problem to be able to provide mathematical stability guarantees, or in adding additional constraints that ensure stability. In [START_REF] Hartley | Initial tuning of predictive controllers by reverse engineering[END_REF], the authors point out that when using zero-value cost functions, adding a terminal constraint could be sufficient to prove stability. An interesting extension of the proposed method would consist in determining how to choose such a terminal constrained without altering the match between MPIC and classical IC or to adapt other controller matching methods, different from the zero-value cost function.

Passivity filters can benefit from a more detailed exploration of storage functions to derive even less conservative passivity conditions. In addition, in their current form, the passivity filters can only be computed in scalar form with the necessity of a decoupled impedance behavior accordingly to the system DoFs. Thus, exploring the formulation of the passivity conditions in matrix form would be an interesting increment of the proposed method. In this work, a specific way of changing the impedance parameters is considered which consists in varying both stiffness and damping in a proportional way. A more generic form of parameter variation could also be considered in order to investigate the design of other types of passivity filters.

The combination of MPIC and passivity filters is performed in the bimanual rehabilitation application. The passivity filter is used to guarantee passive interaction dynamics that are then applied by MPIC, which permits to impose safety workspace constraints. During experimentation, no unexpected behavior or degradation of the desired performance was noticed. However, an interesting extension would be to combing both methods into a more complete framework for interaction management. The prediction capacity of the controller could allow to anticipate constraints but also the variations of the interaction dynamics. This could permit to overcome the issue of too conservative passivity conditions for fast changes of the interaction dynamics. Alternatively, passivity considerations could possibly be used as stability guarantees for MPIC.

In the context of EMG-force modeling, the proposed method should be trialed on a greater number of subjects and tested for the influence of fatigue on the model performance. Also, a better sensor placement could be considered for specific movements, targeting muscles with more influence on the movement. This would require a more precise study of sensor placement considering the different types of movement that are to be modeled.

As the proposed medical applications are technical proofs of concept, they would clearly benefit from a more practical point of view of a clinical practitioner. For the rehabilitation application, adding a solution to compensate for the weight of the patient's arm would be of particular interest for extending the application towards clinical trials. Also, involving impaired patients into the evaluation process could contribute to a better understanding of requirements for clinical applications, as it is very difficult for a healthy subject to emulate impairment.

Finally, in the case of needle insertion, several improvements could lead to a more realistic application. Among them, the question of dealing with moving tissues is particularly interesting, as in practice biological tissues are subject to physiological movements due, for instance, to respiration. An interesting extension of the proposed application would be to test the performance of the layer detection and the changes of interaction dynamics in the case of a moving tissue.

Appendix A Electromyography

Electromyography (EMG) is a diagnostic technique which consists in studying the muscle function by acquiring and analyzing the electrical signals the muscle emits during contraction. EMG signals can either be acquired in a intramuscular way by inserting a needle-shaped electrode inside the muscle tissue, or superficial by placing a pair of electrodes on the skin on top of the targeted muscle. Even though intramuscular EMG allows to have precise information about the studied muscle, surface EMG (sEMG) is the most commonly used acquisition strategy as it is non-invasive. In the following, only sEMG will be considered, and for simplicity it will be referred to as EMG.

A.1 Principle

The description of the neural control of muscle activation starts with the motor unit (MU). A motor unit is the basic functional unit in the neuromuscular system. It consists of a motor neuron (Alpha motoneuron) with its dendrites and axon branches and multiple muscle fibers connected to the axon by endplates [START_REF] Enoka | Physiology of Muscle Activation and Force Generation[END_REF]. A schematic representation of the MU is shown in Figure A.1. In a MU, all innervated fibers react in the same way to an action potential.

In muscle physiology, the excitation of the fibers resulting in muscle contraction is induced by the electrical properties of the muscle fiber membranes. While in steady state, an ionic equilibrium between the inside and the outside of the fiber maintained by a physiological process (ion pump) ensures a constant polarization of the membrane. At rest, the electrical potential of the membrane is maintained to be negative on the inside of the cell with respect to the outside. The resting potential of the fiber membrane is -70mV to -90mV [START_REF] Farina | Biophysics of the Generation of EMG Signals[END_REF]]. An action potential from the central nervous system or reflex induce an activation of the motoneuron inside the spinal cord. The excitation is then conducted along the motoneuron axon towards the endplates. When the action potential reaches the endplate, the release of a transmitter substance induces an endplate potential that modifies the characteristics of the muscle fiber membrane, leading to its depolarization. Locally, the equilibrium is immediately restored and the membrane is repolarized. The action potential, however, spreads along the muscle fiber towards both tendon endings resulting in a depolarization wave. This excitation leads to the shortening of the muscle fiber [Winter, 2009]. plates of each muscle fibers are places on different parts of the fiber, the signals collected in the detection area differs in time, phase and amplitude. This results, typically, in a triphasic MU action potential being the superposition of the biphasic excitation of each fiber composing the MU [Winter, 2009].

Each muscle is composed of a finite number of MU that are controlled separately. The contraction process and force modulation of each muscle is controlled by two main control strategies: recruitment of MU and firing frequency. Recruitment consists in defining the MU that are involved in a contraction process and the firing frequency is the frequency of excitation of each MU. These two strategies highly influence the magnitude and density of the observed EMG signal, as the EMG signal resulting from the muscle action potential and measured on the surface of skin is a superposition of the signals of each recruited MU. In a schematic way, the EMG signal directly reflects the firing and recruitment of each MU composing the studied muscle, which is illustrated in Figure A.4.

It is worth noticing, that the EMG signal measured on the skin surface does not present the same amplitude and frequency characteristics as the original signal inside the muscle. This is caused by the physiological properties of biological tissues that act as spatial and temporal low-pass filters [START_REF] Farina | Biophysics of the Generation of EMG Signals[END_REF]. Another interesting property of EMG signal lies in the fact that the signal can typically be recorded 20 -200ms ahead can be observed. These spikes are the result of the superposition of the activation signals of MUs that are near the electrode and thus exhibit a higher EMG signal amplitude.

As the obtained raw EMG signal is noisy, it is important to process the data in order to efficiently assess the muscle activity.

A.2 Signal processing

The muscle activity information delivered by the EMG signal are noisy, thus effective noise reduction processing is required. A typical EMG data processing strategy includes three main steps [START_REF] Merletti | Surface Electromyography: Physiology, Engineering, and Applications[END_REF], as shown in The first step consists of a band-pass filter with a cut-off frequency at 20Hz and 450Hz [START_REF] Merletti | Surface Electromyography: Physiology, Engineering, and Applications[END_REF]. As the raw EMG signal is a high-frequency signal centered on its baseline, it is important to filter out all undesired low-frequency components such as baseline shifts. These phenomena are caused by different external factors such as electrode shifts due to the relative movements of the skin with respect to the targeted muscle, temperature variations or sweat. The cut-off frequency at 450Hz, on the other hand, is used in order to act as an anti-aliasing filter, as typical EMG acquisition frequencies are chosen around 1kHz. This first processing step is responsible for the elimination of the noise generated by the signal acquisition process.

After applying the first filter, the EMG signal is fed into a full-wave rectifier. In fact, as explained in the previous section, the polarity of the signal acquired by the electrodes depends on the position of the electrode with respect to the muscle fiber endplate. For this reason, the muscle activity only depends on the amplitude of the acquired signal and thus need to be rectified [START_REF] Dakin | Rectification is required to extract oscillatory envelope modulation from surface electromyographic signals[END_REF].

Envelope extraction is the final step of the proposed processing strategy. This signal processing step directly reflects the physiological processes of the muscle activation. In fact, the muscle acts as a low-pass filter taking as input a high frequency action potential and generating low frequency muscle activation. For this reason, this filtering step is typically accomplished by the mean of a low-pass filter with a low cut-off frequency. As only the modulation of the input signal is of interest, the cut-off frequency of this filter is chosen between 1 and 2Hz, which corresponds to the physiological muscle model [Winter, 2009].

In the literature, other techniques of advanced EMG processing for noise and artifact removal are proposed. Among them, the most popular are signal whitening [START_REF] Liu | Signal whitening preprocessing for improved classification accuracies in myoelectric control[END_REF], wavelet transform [START_REF] Phinyomark | A comparative study of wavelet denoising for multifunction myoelectric control[END_REF], independent component analysis [START_REF] Tapia | Emg signal filtering based on independent component analysis and empirical mode decomposition for estimation of motor activation patterns[END_REF] and empirical mode decomposition [START_REF] Andrade | Emg signal filtering based on empirical mode decomposition[END_REF], to mention just a few. Although these methods are shown to improve the signal quality, they are computationally expensive and thus not suited for real-time applications [START_REF] Bi | A review on emg-based motor intention prediction of continuous human upper limb motion for human-robot collaboration[END_REF]. Finally, it is common practice to normalize the obtained muscle activation signal with respect to the maximum isometric voluntary contraction (MVC) of the targeted muscle [De Luca, 1997]. The resulting quantity then becomes a percentage of the MVC value, which depends on the subject physiology and sensor placement, and thus can easily be compared to other EMG recordings. In fact, even though experimental procedures that regulate the acquisition process, such as the SENIAM recommendations [START_REF] Hermens | Development of recommendations for semg sensors and sensor placement procedures[END_REF] for electrode placement, tend to ensure the highest repeatability, the acquisition conditions always differ for each procedure. However, the major drawback for MVC normalization results from the need of acquiring the MVC for each session and for each targeted muscle, which is time consuming and impractical for applications that require a fast setup.

Appendix B Résumé de la thèse B.1 Introduction

Depuis quelques années, les manipulateurs robotiques sont introduits dans les usines pour effectuer de nombreuses tâches, car ils sont plus puissants, plus précis et ne se fatiguent pas, ce qui rend la production plus efficace. Les architectures de commande de robot basées sur le contrôle de position sont peu efficaces en présence de changements de l'environnement, rendant dangereux pour des opérateurs humains de travailler dans un espace partagé avec ces robots. Pour cette raison, dans la plupart des cas, une séparation stricte entre les robots et les travailleurs a été mise en place dans les usines. Même si séparer les robots des humains est un moyen efficace d'assurer la sécurité et d'éviter les contacts physiques, restreindre ainsi les aires de travail dans les usines est très contraignant et prend beaucoup de place. De nombreuses autres applications bénéficieraient de la combinaison humain/robot, i.e. des meilleures performances du robot en termes de précision, de répétabilité et de capacité de charge avec les capacités cognitives de l'opérateur humain telles que l'adaptabilité et le jugement rapide en cas d'événements imprévus. Cependant, supprimer les séparations et laisser des espaces partagés de travail humain/robot conduit à un changement complet de paradigme et nécessite de repenser la technologie. En ouvrant l'espace de travail, des stratégies de contrôle dynamiques sont nécessaires pour assurer la sécurité. Ainsi l'évolution récente de la robotique tend vers une conception des mécanismes robotiques plus orientée vers la collaboration avec l'Homme et des stratégies de contrôle, conduisant à une nouvelle classe de systèmes robotiques appelés cobots, pour robots collaboratifs.

Dans l'industrie moderne, la flexibilité de production, c'est-à-dire la capacité de reconfigurer facilement les lignes de production, est d'une importance particulière, notam-ment dans les domaines où les produits ont une durée de vie courte et s'opèrent des changements de production fréquents. Cette flexibilité est directement liée à la facilité de déploiement et de programmation des solutions robotiques. Les solutions cobotiques sont accompagnées d'outils de programmation simplifiés par rapport aux solutions standards qui nécessitent souvent l'intervention d'intégrateurs spécialisés. Les cobots peuvent être programmés à l'aide d'interfaces Homme-machine intuitives et interactives permettant à l'opérateur d'enseigner au robot une tâche souhaitée. Cela permet à l'opérateur de facilement faire en sorte que le robot prenne en charge complètement ou partiellement les tâches les plus ennuyeuses et se concentre sur les plus complexes. Dans ce contexte, les cobots présentent un intérêt particulier, notamment en raison de leur prix attractif par rapport aux manipulateurs robotiques classiques, ainsi que de leur flexibilité. Les solutions cobotiques ne sont ainsi pas toujours déployées pour le besoin de collaboration Homme-robot, mais aussi pour des raisons spatiales et économiques.

Dans le domaine de la robotique médicale, les cobots attirent beaucoup d'attention. Actuellement, les solutions médicales cobotiques se concentrent sur deux domaines principaux: les interventions médicales assistées par robot et la rééducation fonctionnelle. Dans le premier cas, les cobots sont utilisés pour assister le praticien lors de l'intervention en améliorant ses capacités en termes de stabilité et de précision. Cela se traduit par une réduction des risques pour le patient dus aux erreurs de positionnement et une amélioration de la répétabilité souvent impactée par la fatigue des praticiens. Cela permet également de réduire dans une certaine mesure la durée des procédures, et par là même les coûts.

La deuxième application des cobots dans le domaine médical est la rééducation fonctionnelle. Ce type de procédure s'adresse en particulier aux patients présentant une atteinte partielle des membres due à un accident vasculaire cérébral ou à des lésions de la moelle épinière. Dans les scénarios de rééducation classique, un thérapeute assiste le patient dans ses mouvements, ce qui prend du temps et a souvent pour effet de raccourcir la durée de la rééducation. Les tâches de rééducation sont bien adaptées à l'automatisation, car elles nécessitent principalement de nombreuses répétitions d'un mouvement donné. C'est pourquoi plusieurs dispositifs robotiques spécifiques à des tâches ont été proposés. Plus récemment, des efforts de recherche ont été consacrés à l'étude de l'utilisation des cobots industriels pour les procédures de réhabilitation. Le rôle de la solution robotique est alors d'assister les mouvements du membre atteint lorsque le patient n'est pas capable d'exécuter lui-même la tâche de rééducation souhaitée, et de limiter l'étendue du mouvement afin d'éviter les blessures.

B.2. CONTR ÔLE EN IMP ÉDANCE SOUS CONTRAINTES

L'objectif de cette thèse est de proposer des outils généraux pour aborder les enjeux actuels de la robotique collaborative et de mettre en oeuvre ces outils sur les deux principaux scénarios cobotiques en robotique médicale.

B.2 Contrôle en impédance sous contraintes

Le premier défi abordé concerne les aspects de sécurité lors de l'interaction avec un robot. En fait, les normes de sécurité actuelles donnent des limites générales aux systèmes cobotiques lors d'une interaction. Cependant, dans la plupart des cas, il est important de pouvoir définir certaines limitations spécifiques aux tâches sur le mouvement du robot, telles que les contraintes de position, de vitesse ou d'accélération. Dans les scénarios robotiques généraux, la mise en oeuvre de ces contraintes intervient directement dans la planification de la tâche.

En robotique, lorsqu'aucune interaction avec l'environnement n'est nécessaire, des stratégies de contrôle de mouvement sont utilisées. Au contraire, en présence d'interactions, les stratégies de contrôle direct des forces sont préférées pour le suivi des forces fines. Cependant, ils nécessitent de bons modèles d'interaction et de l'environnement, et ne sont pas compatibles avec des interactions imprévisibles. Ceci n'est cependant pas possible dans le cadre de la collaboration car l'influence de l'opérateur est difficile à évaluer. Pour cette raison, le contrôle d'interaction offre un compromis afin de gérer à la fois la force et le mouvement lors des interactions. Le contrôle en impédance (IC) est une méthode de contrôle d'interaction largement utilisée, particulièrement efficace pour les interactions spontanées. Il consiste à imposer un modèle d'impédance pour la relation entre le manipulateur et l'environnement. D'autre part, avec le développement rapide de nouveaux outils plus puissants pour résoudre les problèmes d'optimisation avec gestion des contraintes, beaucoup d'efforts de recherche ont été faits pour utiliser cette technologie dans le contrôle des robots. L'optimisation des contraintes peut être appliquée dans diverses situations telles que la mise en place de barrières virtuels que le robot n'est pas autorisé à franchir ou des contraintes énergétiques. De nos jours, Les solveurs modernes sont suffisamment performants pour gérer l'optimisation en ligne, ce qui est requis pour les stratégies de contrôle bas niveau, telles que le contrôle des robots au niveau du couple articulaire.

Un moyen courant de traiter les problèmes de contrôle multi-variables contraint avec la capacité de prédiction est la commande prédictive (MPC). L'avantage du MPC par rapport aux autres méthodes d'optimisation sous contraintes réside dans sa capacité à prédire comment et quand le système contrôlé répond aux contraintes imposées. Ce faisant, MPC est capable de réagir en conséquence, ce qui entraîne des transitions plus douces entre les mouvements libres et les mouvements contraints.

La contribution proposée consiste en une méthode de conception de contrôleurs basée sur une stratégie de contrôle prédictif (MPC) qui se comporte comme un contrôleur en impédance classique. Cela permet à la fois de gérer la dynamique de l'interaction tout en garantissant les contraintes pratiques. Le contrôleur proposé est ainsi appelé commande prédictive en impédance (Model Predictive Impedance Control -MPIC). Les résultats expérimentaux obtenus montrent que la méthode proposée combine les caractéristiques du contrôle en impédance avec les avantages de la commande prédictive, permettant ainsi d'imposer des contraintes sur différents paramètres du système, tels que l'état ou l'entrée. MPIC se comporte exactement comme le contrôleur en impédance classique lors qu'aucune contrainte n'est active et il respecte les contraintes de position, de vitesse et d'accélération lorsque de telles contraintes sont imposées. notamment dans le cas des systèmes collaboratifs Homme-robot. Même si un système robotique entièrement automatisé présente des performances de précision, de répétabilité et de capacité de charge supérieures à celles d'un opérateur humain lorsqu'il fonctionne dans un environnement bien défini, de nombreuses tâches nécessitent la capacité humaine de jugement rapide et d'adaptation en cas d'événements imprévus. Dans ce contexte, les systèmes collaboratifs sont conçus pour combiner à la fois la force et la précision du robot, ainsi que le jugement et la flexibilité de l'opérateur. L'un des principaux défis d'une collaboration efficace entre l'humain et le robot reste dans la communication entre les deux côtés et surtout dans la compréhension des intentions de l'opérateur. La communication présente un intérêt particulier pour les tâches d'assistance au robot dans lesquelles le robot doit être capable de détecter l'étendue de l'activité de l'opérateur lors de l'exécution de la tâche et d'adapter le niveau d'assistance en cas de besoin. Cependant, cela nécessite de concevoir et d'établir des normes de communication qui soient simultanément intuitives pour l'utilisateur et qui donnent au système robotique des informations suffisantes sur les intentions humaines dans les différentes phases de la collaboration.

Parmi les différent canaux d'échange d'information entre le robot et l'utilisateur celle basée sur l'évaluation des apports de force des deux partenaires et particulièrement intéressante. L'activité humaine peut ainsi être évaluée en mesurant l'étendue de la force appliquée par l'opérateur sur le système robotique. Cela peut généralement être fait à l'aide de capteurs de force-couple mécaniques. En raison de sa simplicité, cette stratégie est couramment utilisée dans de nombreuses applications collaboratives telles que la manipulation d'objets. L'inconvénient majeur d'une telle méthode réside dans le fait que le capteur mesure non seulement l'interaction Homme-robot souhaitée, mais également des composants indésirables tels que les forces de gravité et de frottement, mais également des forces d'interaction inédites en cas de contact externe avec un environnement incertain. Cela rend l'utilisation de tels capteurs inadaptée à des tâches complexes avec la nécessité simultanée d'évaluer l'intention humaine et les forces environnementales, car il est impossible de faire une distinction entre elles. Alternativement, des bio-signaux mesurables tels que l'électromyographie de surface peuvent également être utilisés pour l'évaluation de l'activité humaine L'électromyographie (EMG) est une technique de diagnostic qui consiste à étudier la fonction musculaire en acquérant et en analysant les signaux électriques émis par le muscle lors de la contraction. L'EMG a trouvé de nombreux domaines d'application, notamment en biomécanique, en bio-ingénierie et en médecine. Cette technologie est utilisé en neurophysiologie, en ergonomie du mouvement, ainsi que pour l'analyse de la marche B.4. AM ÉLIORATION DES INTERACTIONS HOMME-ROBOT ou pour la rééducation. Les signaux EMG peuvent être acquis de manière superficielle en plaçant une paire d'électrodes sur la peau au-dessus du muscle ciblé. Les signaux EMG sont le résultat de l'activité bioélectrique des cellules musculaires. Pour cette raison, ils sont souvent utilisés dans les stratégies de contrôle des systèmes robotiques d'assistance car ils reflètent directement les contractions musculaires volontaires de l'utilisateur et, par conséquent, donnent un aperçu de ses intentions.

La contribution proposée consiste à établir un modèle de la dynamique du bras humain à partir de l'électromyographie (EMG) et des mesures de force, qui permet ensuite de modifier la dynamique de l'interaction. Les résultats expérimentaux obtenus (comme illustré Les résultats obtenus dans le cadre de ce travail montrent également qu'un modèle linéaire assez précis de la force produite par le bras humain basé sur l'évaluation de l'activité musculaire de l'épaule peut être obtenu de manière récursive, ce qui est particulièrement intéressant car un tel modèle est indépendant de l'utilisateur et peut gérer la variabilité de la relation EMG/force. En outre, les changements de profil d'impédance mis en oeuvre ne sont influencés que par les forces appliquées par le bras équipé d'EMG, même en cas de variations des forces environnementales. Afin de faire face aux exigences des différentes applications pour lesquelles la stratégie proposée pourrait être utilisée, la modification des paramètres d'impédance de l'interaction peut être librement adaptée selon les besoins, car l'utilisation du filtre de passivité garantit une stabilité et ainsi exécution sûre des changements d'impédance. Parmi les applications pratiques qui pourraient bénéficier d'une telle stratégie d'interaction, les tâches de polissage collaboratives pourraient être intéressantes, car elles nécessitent d'assurer un contact constant avec l'objet manipulé et la capacité d'interaction avec un opérateur. Une autre utilisation, qui est explorée dans le développement de cette thèse, réside dans le domaine de la rééducation fonctionnelle pour les patients handicapés post-AVC.

B.5 Applications médicales

Avec le passage progressif des robots des environnements industriels fermés à des environnements humains plus ouverts et encombrés, les procédures médicales comportant des solutions robotiques attirent de plus en plus d'attention. Les tendances actuelles de la robotique médicale visent à inclure les cobots dans la pratique clinique, en particulier dans les domaines de la rééducation fonctionnelle et des interventions médicales assistées par robot.

Dans cette partie, des implémentations pratiques des méthodes de contrôle proposées pour les robots collaboratifs sont présentées dans le contexte médical.

B.5.1 Insertion d'aiguille assistées par robot

Dans le travail effectué, une première application dans le domaine de l'aide au praticien est proposée avec une nouvelle méthode de contrôle pour l'insertion d'aiguille, qui utilise les outils proposés. En effet, ces dernières années, on observe une généralisation des procédures médicales minimalement invasives, en raison de leur impact drastique sur le bien-être du patient. Ce type de procédures permet de réduire significativement le temps d'hospitalisation mais aussi l'inconfort pour le patient en réduisant la douleur, les cicatrices, le risque d'infections et de saignements post-intervention. Des procédures min-B.5. APPLICATIONS M ÉDICALES imalement invasives peuvent être effectuées en accédant aux zones d'intérêt par différents moyens, par exemple à travers des orifices naturels ou directement en perforant la peau avec une aiguille. Cette dernière stratégie dite percutanée est couramment utilisée en radiologie interventionnelle, discipline médicale combinant un savoir-faire médical spécialisé en radiologie avec les dernières innovations technologiques dans les domaines de l'imagerie médicale et éventuellement de la robotique.

Dans ce travail, une application d'insertion automatique d'aiguille dans une succession de tissus mous est proposée. Dans cette application, le robot est capable de détecter les différentes transitions entre les tissus et adapter sa raideur en fonction du tissu dans lequel pénètre la pointe de l'aiguille. La passivité de la modification de la dynamique de l'interaction est assurée par le filtre de passivité.

B.5.2 Rééducation bimanuelle

La deuxième application est une preuve de concept pour la rééducation bimanuelle, qui est une approche plutôt nouvelle de la rééducation fonctionnelle, qui est étudiée de manière exploratoire, plutôt du point de vue technique que clinique. En effet, les lésions tels que les traumatismes crâniens ou les AVC sont des accidents qui souvent conduisent à une perte de mobilité des membres du patient. La rééducation robotisée est particulièrement intéressante dans ce cas de situation, car elle permet d'effectuer des exercice de rééducation de grande intensité et avec de nombreuses répétitions. Cela permet de stimuler la plasticité cérébrale et ainsi favorise la récupération des capacités motrices du patient. Contrairement à la rééducation classique qui solicite uniquement les membres diminués, la rééducation bimanuelle favorise la coordination des deux bras et ainsi stimule les connections neuronales. De plus, la rééducation bimanuelle permet une plus rapide récupération de l'autonomie des patient car une grande majorité des mouvements quotients nécessitent l'utilisation des deux membres. La méthode proposée se base sur le principe que lors d'une tâche tel que le transport d'un objet, le patient avec un membre atteint va chercher à compenser sa faiblesse avec l'autre bras. Ainsi, le système proposé raidi le robot lorsqu'un déséquilibre entre l'engagement des deux membres est détecté et empêche ainsi le mouvement. Dans le cas ou les deux bras collaborent avec la même intensité, la raideur du système diminue et favorise le mouvement. Le dispositif expérimental de l'application est illustré par la Figure B.4. Dans cette application, la passivité des modifications de la dynamique de l'interaction est garantie par le filtre de passivité et MPIC est utilisé pour imposer des contraintes sur le système pour garantir que le patient ne se blesse pas lorsqu'il effectue les mouvements.

B.6 Conclusion

Avec l'intérêt croissant pour l'utilisation de systèmes robotiques dans un espace de travail partagé avec des opérateurs humains, le développement de robots collaboratifs offre de nouvelles possibilités d'interactions Homme-robot. Par conséquent, le développement de nouveaux outils de contrôle pour la gestion des interactions est devenu un sujet de recherche important. Dans ce contexte, la conception de solutions de contrôle qui, à la fois, permettent d'améliorer la dynamique d'interaction et de garantir l'intégrité des opérateurs est particulièrement intéressante.

Dans cette thèse, plusieurs outils de contrôle pour la robotique collaborative, qui assurent la sécurité et améliorent la dynamique des interactions physiques Homme-robot, Alors que dans de nombreuses applications, la dynamique d'interaction souhaitée peut être considérée comme constante, lors de l'interaction physique Homme-robot, la capacité d'adapter la dynamique de l'interaction est souvent requise. Un problème important et dangereux qui résulte de la modification de la dynamique d'interaction est la possible perte de passivité qui peut conduire à une instabilité et à un comportement inattendu du système robotique. Par conséquent, assurer un changement passif de la dynamique d'interaction est une étape importante vers une interaction physique Hommerobot améliorée qui garantie l'intégrité de l'opérateur. Dans ce travail, une méthode de conception des filtres de passivité est proposée. Ces filtres garantissent la passivité de

1. 1

 1 Examples of commercial cobots. . . . . . . . . . . . . . . . . . . . . . . . 1.2 Examples of medical cobots. . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Main techniques for collaborative operations given by ISO 10218 and ISO/TS 15066 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.1 Robot dynamic model linearization. . . . . . . . . . . . . . . . . . . . . . 2.2 Impedance control. The reference configuration is displayed in light blue. 2.3 Robot dynamic model linearization. . . . . . . . . . . . . . . . . . . . . . 2.4 Compliant frame from impedance model. . . . . . . . . . . . . . . . . . 2.5 Admittance control in task-space. . . . . . . . . . . . . . . . . . . . . . 3.1 Comparison between MPIC and constrained QP optimization. Position and velocity constraints are at the limit of the red areas. . . . . . . . . . 3.2 Experimental setup with reference trajectory and an obstacle. . . . . . . 3.3 Experimental comparison of the positions, velocities and accelerations of the robotic manipulator end-effector, following a reference trajectory. Case 1: IC and MPIC without active constraints and with an obstacle on the path. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.4 End-effector trajectories. Case 2: IC and MPIC with position constraints (x ≥ -0.006m) and with an obstacle on the path. . . . . . . . . . . . . 3.5 Experimental comparison of the trajectories and velocities of the endeffector following a reference. Case 3: IC and MPIC with active velocity constraints (v x ≤ 0.02m/s, v y ≤ 0.02m/s) and with an obstacle on the path. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.6 Experimental comparison of acceleration and measured force applied by the robotic manipulator end-effector. Case 4: IC and MPIC control with acceleration constraints (a y ≥ -0.4m/s 2 ) and with an obstacle on the path. vii 6.6 Experimental setup for bimanual rehabilitation. The weak (red) and strong (green) arms produce a weak f w and a strong f s force, respectively. The EMG sensors are placed on the strong arm. . . . . . . . . . . . . . . . . 6.7 Experimental results for bimanual rehabilitation strategy. Green areas correspond to phases where the patient pushes equally with both hands. 6.8 Phases of needle insertion into a layer of soft tissue. . . . . . . . . . . . . 6.9 Kelvin-Voigt model representation to describe viscoelastic interaction between a needle and a soft tissue before puncture. . . . . . . . . . . . . . . 6.10 Experimental setup. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.11 Layer detection results. . . . . . . . . . . . . . . . . . . . . . . . . . . . 6.12 Evolution of passivity conditions and γ ˙. For better readability, conditionh 2 was scaled by a factor 10. . . . . . . . . . . . . . . . . . . . . . . . . . 6.13 Varying stiffness and damping modified by the passivity filter. . . . . . . A.1 Motor Unit. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.2 Depolarization wave on muscle fiber membranes and resulting signal. . . A.3 Motor unit excitation signal. . . . . . . . . . . . . . . . . . . . . . . . . . A.4 Composition of surface EMG signal. . . . . . . . . . . . . . . . . . . . . . A.5 Raw EMG signal resulting from biceps brachii contractions. . . . . . . . A.6 Signal processing steps from raw EMG data to muscle activity. . . . . . . A.7 Raw EMG data and the resulting muscle activation signal. For better readability, the muscle activation signal was scaled by a factor 15. . . . B.1 Validation expérimentale de MPIC avec contraintes de position, un obstacle et interactions physiques Homme-robot. . . . . . . . . . . . . . . . . B.2 Visualisation des résultats expérimentaux. . . . . . . . . . . . . . . . . . B.3 Dispositif expérimental. . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.4 Dispositif expérimental pour la rééducation bimanuelle. Les bras faible (rouge) et fort (vert) produisent une force faible f w et forte f s , respectivement. Les capteurs EMG sont placés sur le bras fort. . . . . . . . . . . .

Figure 1 . 1 :

 11 Figure 1.1: Examples of commercial cobots.

  (a) Rosa Brain from Medtech (medtech.fr). (b) TMS-Cobot from Axilum Robotics (axilumrobotics.com).

Figure 1 . 2 :

 12 Figure 1.2: Examples of medical cobots.

❼

  Safety-rated monitored stop ❼ Hand-guiding ❼ Speed and separation monitoring ❼ Power and force limiting These modes are summarized in Figure 1.3 and further explained in the following development.

  Figure 1.3: Main techniques for collaborative operations given by ISO 10218 and ISO/TS 15066

  (a) Impedance control in joint-space.(b) Impedance control in task-space.

Figure 2

 2 Figure 2.2: Impedance control. The reference configuration is displayed in light blue.

Figure 2 . 3 :

 23 Figure 2.3: Robot dynamic model linearization.

Figure 2 . 4 :

 24 Figure 2.4: Compliant frame from impedance model.

  (a) Control input, equivalent to acceleration. The dashed red area corresponds to the response to the external disturbance. (b) System velocity and position.

Figure 3 . 1 :

 31 Figure 3.1: Comparison between MPIC and constrained QP optimization. Position and velocity constraints are at the limit of the red areas.

Figure 3 . 2 :

 32 Figure 3.2: Experimental setup with reference trajectory and an obstacle.

  (a) End-effector trajectories under IC and MPIC. (b) End-effector velocities under IC and MPIC. (c) End-effector accelerations under IC and MPIC.

Figure 3

 3 Figure 3.3: Experimental comparison of the positions, velocities and accelerations of the robotic manipulator end-effector, following a reference trajectory. Case 1: IC and MPIC without active constraints and with an obstacle on the path.

Figure 3

 3 Figure 3.4: End-effector trajectories. Case 2: IC and MPIC with position constraints (x ≥ -0.006m) and with an obstacle on the path.

  (a) End-effector velocities under IC and MPIC. (b) End-effector trajectories under IC and MPIC.

Figure 3

 3 Figure 3.5: Experimental comparison of the trajectories and velocities of the end-effector following a reference. Case 3: IC and MPIC with active velocity constraints (v x ≤ 0.02m/s, v y ≤ 0.02m/s) and with an obstacle on the path.

  (a) End-effector accelerations under IC and MPIC. (b) Force measured at the end-effector under IC and MPIC.

Figure 3

 3 Figure 3.6: Experimental comparison of acceleration and measured force applied by the robotic manipulator end-effector. Case 4: IC and MPIC control with acceleration constraints (a y ≥ -0.4m/s 2 ) and with an obstacle on the path.

  (a) Experimental setup. (b) End-effector trajectory. (c) Force measured at the end-effector. (d) MPIC-generated acceleration input.

Figure 3

 3 Figure 3.7: Experimental validation of MPIC with position constraints (x ≥ 0m), with an obstacle on the path and pHRI.

t 0 r

 0 T f ext dτ = t 0 {r T (M e ¨p + De ˙p + Ke p )}dτ = t 0 {r T (M (e ¨p + αe ˙p) -αM e ˙p + De ˙p + Ke p )}dτ Using equation (4.7) and its derivative, it follows that t 0 r T f ext dτ = t 0 {r T M r ˙+ e ˙T p (D -αM )e ˙p + e T p (K + αD -α 2 M )e ˙p + e T p (αK)e p )}dτ (4.8) By introducing a symmetric, positive semidefinite and continuously differentiable matrix C, equation (4.8) can be written t 0 r T f ext dτ = t 0 {r T M r ˙+ e T p Ce ˙p + 1 2 e T p C ˙ep }dτ + t 0 {e ˙T p (D -αM )e ˙p + e T p (αK -

Figure 4 . 1 :

 41 Figure 4.1: Passivity filter for γ i using conditions from V 2 .

Lemma 4. 1 .

 1 Considering the filtered tracking error as defined in (4.21) with M a constant diagonal matrix and B a time-dependent diagonal matrix with bounded eigenvalues, if r is bounded then e p , e ˙p are bounded and that if r → 0 then e p , e ˙p → 0.

  e ˙pi → 0 : Considering the 2-norm for a scalar function f ∈ R defined by ||f || 2 = ∞ 0 f 2 (τ )dτ , one can show that ||e p i || 2 is bounded. This is done by bounding ||e -t 0 b i (τ ) m i dτ e p i (0)|| 2 =

Figure 4 . 2 :

 42 Figure 4.2: Passivity filter for γ i using conditions from V 4 .

Figure 4 . 3 :

 43 Figure 4.3: Practical passivity filter for γ i using conditions from V 4 .

  The simulated system has a constant mass m = 10kg, a stiffness varying between k 0 = 2N/m and k 1 = 22N/m and a damping ratio of 0.1. The considered reference trajectory was set such that p r = 10 sin 0.1t. The stiffness is expected to vary according to the switching function γ ¯illustrated in Figure 4.4 in dashed gray. The passivity filters (4.20) and (4.41), associated to storage functions V 2 and V 4 , respectively, are tuned with the parameter β = 10. The resulting switching functions γ are represented in Figure 4.4.

Figure 4 . 4 :

 44 Figure 4.4: Switching variable γ resulting from filter with conditions from V 2 and V 4 and its reference.

  (a) Storage functions. (b) Derivative of storage functions.

Figure 4 . 5 :

 45 Figure 4.5: Time evolution of storage functions for γ calculated to respect passivity constraints (4.20) (based on V 2 ).

  (a) Storage functions. (b) Derivative of storage functions.

Figure 4

 4 Figure 4.6: Time evolution of storage functions for γ calculated to respect passivity constraints (4.41) (based on V 4 ).

Figure 4 . 7 :

 47 Figure 4.7: Experimental setup.

  Figure 4.8. The variable impedance parameters are switched according to Figure 4.8a between K 0 = 100N/m and K 1 = 500N/m for the stiffness. The damping is computed with a damping ratio of 1 and varies between D 0 = 44.7Ns/m and D 1 = 100Ns/m. As it can be observed in Figure 4.8c, the system (a) Impedance stiffness and damping profiles. (b) End-effector task-space position and reference (dashed lines). (c) Passivity condition.

Figure 4 . 8 :

 48 Figure 4.8: Experimental results of varying impedance profiles without the passivity filter.

  Figure 4.9a between K 0 = 100N/m and K 1 = 500N/m for the stiffness. The damping is computed with a damping ratio of 1 and varies between D 0 = 44.7Ns/m and D 1 = 100Ns/m. The passivity filter gain is chosen such that β = 50.

  (a) Impedance stiffness and damping profiles.(b) Variations of γ ˙and passivity conditions. For better readability, the passivity condition h 2 was scaled by a factor 3.

Figure 4 . 9 :

 49 Figure 4.9: Varying impedance profiles and evolution of the switching variable using the passivity filter.

  The robot behavior, resulting from the impedance profile variations is shown in. In contrary to the previous experimentation, this time one can observe that the relationship V ˙4 ≤ r T f ext holds, as it can be seen in Figure4.10b, and thus passivity is ensured. This results in a stable motion of the robot end-effector, as shown in Figure4.10a.

  (a) End-effector task-space position and reference (dashed lines). (b) Passivity condition.

Figure 4 .

 4 Figure 4.10: Experimental results of varying impedance profiles using the passivity filter.

Figure 4 .

 4 Figure 4.11: Robot behavior in both experiments.

Figure 5

 5 Figure 5.1: Raw EMG data and the resulting muscle activation signal. For better readability, the muscle activation signal was scaled by a factor 15.

Figure 5 . 2 :

 52 Figure 5.2: Movement configuration for model validation.

  (a) EMG sensor placement on the elbow. (b) Muscle activity assessed from EMG sensors placed on the elbow. For better readability, the activation values where scaled to the maximal value for all channels. (c) Comparison between the estimated and the measured force generated by elbow movements. The model is adapting the model during white periods and using the model for estimation during red periods. During estimation (red) periods: FIT = 83.54% and VAF = 97.72%.

Figure 5 . 3 :

 53 Figure 5.3: Experimental model validation about the human elbow.

  (a) EMG sensor placement on the shoulder. (b) Muscle activity assessed from EMG sensors placed on the shoulder. For better readability, the activation values where scaled to the maximal value for all channels.(c) Comparison between the estimated and the measured force generated by shoulder movements. The model is adapting the model during white period and using the model for estimation during red period. During estimation (red) period: FIT = 71.67% and VAF = 93.56%.

Figure 5

 5 Figure 5.4: Experimental model validation about the human shoulder.

  (a) Environmental force only. (b) Human and environmental forces.

Figure 5 . 5 :

 55 Figure 5.5: System under VIC.

Figure 5

 5 Figure 5.6: Block-diagram illustrating the overall system architecture during humanrobot interaction.

Figure 5 . 7 :

 57 Figure 5.7: Experimental setup.

  (a) Comparison of the measured external force with the estimated human force. (b) Stiffness and damping profiles before and after the passivity filter. (c) Handle Cartesian position in the robot y-axis.

Figure 5

 5 Figure 5.8: Experimental data of EMG-based VIC.

Figure 5 . 9 :

 59 Figure 5.9: Visualization of experimental results.

Figure 6 . 2 :

 62 Figure 6.2: Examples of commercial end-effector rehabilitation robots.

Figure 6 . 3 :

 63 Figure 6.3: Robert from Life Science Robotics (lifescience-robotics.com).

  Figure 6.5: Representation of possible values for Γ(f s , f w ).

  (a) Measured external forces and estimated user input. (b) Stiffness and damping profiles.

  (c) Handle Cartesian position in the robot y-axis.

Figure 6 . 7 :

 67 Figure 6.7: Experimental results for bimanual rehabilitation strategy. Green areas correspond to phases where the patient pushes equally with both hands.

  Figure 6.8: Phases of needle insertion into a layer of soft tissue.

  (a) Needle force estimation and layer detection. (b) CuSum function and e 2 .

Figure 6 .

 6 Figure 6.11: Layer detection results.

Figure 6 .

 6 Figure 6.13: Varying stiffness and damping modified by the passivity filter.

Figure

  Figure A.1: Motor Unit.

  Figure A.2: Depolarization wave on muscle fiber membranes and resulting signal.

  Figure A.6.

Figure A. 6 :

 6 Figure A.6: Signal processing steps from raw EMG data to muscle activity.

  The resulting filtered signal reflects the effort produced by the targeted muscle and thus is referred to as muscle activation signal. A comparison of the activation signal with the original raw EMG signal is shown in FigureA.7.

Figure A. 7 :

 7 Figure A.7: Raw EMG data and the resulting muscle activation signal. For better readability, the muscle activation signal was scaled by a factor 15.

  En effet, lors d'interactions avec un opérateur humain sous contraintes de position, comme illustré par la Figure B.1, on peut voir qu'à l'intérieur de la zone autorisée, le système se comporte comme contrôlé en impédance et réagit ainsi en conséquence à la force appliquée par l'opérateur et rejoint la trajectoire de référence lorsque l'opérateur cesse d'interagir. Il n'est cependant pas possible de sortir de la zone autorisée, même en poussant contre la limite car le système contrôlé par MPIC génère des forces opposées.

Figure B. 1 :

 1 Figure B.1: Validation expérimentale de MPIC avec contraintes de position, un obstacle et interactions physiques Homme-robot.

  Figure B.2) montrent qu'en utilisant la stratégie proposée, il est possible d'utiliser des informations sur l'influence de l'utilisateur basées sur l'estimation de sa force avec des signaux EMG et de régler un contrôleur d'impédance variable qui modifie le comportement dynamique du robot en conséquence. Ainsi, comme illustré dans le cas de l'expérience de la Figure B.2, ce changement de comportement permet le déplacement du robot.

Figure B. 2 :

 2 Figure B.2: Visualisation des résultats expérimentaux.

Figure B. 3 :

 3 Figure B.3: Dispositif expérimental.

Figure B. 4 :

 4 Figure B.4: Dispositif expérimental pour la rééducation bimanuelle. Les bras faible (rouge) et fort (vert) produisent une force faible f w et forte f s , respectivement. Les capteurs EMG sont placés sur le bras fort.

Table 3 .

 3 1: Impedance parameters retrieved from data.

	M v	K v	D v
	x 5.40 290.04 76.72
	y 5.39 280.27 72.80

Table 6 .

 6

	1.
	Layer	Approximate stiffness [N/m] Thickness [mm]
	Skin	150	14
	Fat	50	20
	Muscle	300	15
	Liver	25	35
	Table 6.1: Anatomical phantom characteristics.

Proceedings Volumes, 44(1):6863-6871. Cited on page 50.

Acknowledgements

List of Tables

of the effective muscle contraction [Winter, 2009].

In order to understand the characteristics of the acquired EMG signal, an experimental raw EMG signal recording of several contraction bursts of the biceps brachii is shown in Figure A.5. In the acquired signal, two events can be distinguished: the rest periods and the contraction bursts. During the rest periods, when the muscle is relaxed, a rather low noise EMG baseline can be observed. Alternatively, during contraction bursts, the raw EMG signal consists of non-reproducible activation patterns recorded from the activity of the MU. This non-reproducibility of the activation shape is caused by the constantly changing recruitment of available MUs, which is also the reason why at some points spikes

B.3 Modification de la dynamique de l'interaction

Le deuxième défi abordé consiste à modifier la dynamique de l'interaction tout en garantissant la passivité du robot. En fait, alors que le contrôle en impédance est l'une des stratégies les plus couramment utilisées pour le contrôle d'interaction des robots, le contrôle d'impédance variable est une préoccupation plus récente. Si la conception d'un contrôle d'impédance avec des paramètres variables permet d'augmenter la flexibilité et la dextérité du système, cela reste un problème difficile, car il peut entraîner une perte de passivité du système contrôlé. Ceci a un impact important sur la stabilité et donc sur la sécurité de l'interaction.

L'un des principaux défis de la gestion des interactions avec IC réside dans la sélection appropriée des paramètres d'impédance, non seulement en fonction de la tâche, mais également pour assurer la stabilité de toutes les variations de paramètres possibles. En effet, si l'utilisation de paramètres d'impédance fixe rend le système passif et donc stable lorsqu'il interagit avec un environnement passif, cette propriété ne tient plus pour des paramètres variant arbitrairement. Ne pas pouvoir modifier l'impédance du système pendant la tâche est un inconvénient notable qui peut limiter les plages d'application. Par exemple, l'interaction Homme-robot peut bénéficier de la capacité d'adaptation du robot à la force appliquée par un utilisateur ou à sa fatigue et ainsi modifier dynamiquement l'interaction.

Dans cette partie, la contribution consiste à proposer des méthodes pour concevoir des filtres de passivité garantissant la passivité de l'interaction. Les filtres proposés sont basés sur des conditions de passivité issus de l'analyse de fonctions de stockage permettant de conclure sur la passivité du système après la modification de la dynamique de l'interaction. Leur rôle est de vérifier si un profil d'impédance souhaité est passif, ou de le modifier si nécessaire. Les résultats expérimentaux montrent que l'utilisation des filtres proposés permet de suivre des profiles de paramètres d'impédance tout en garantissant la passivité, ce qui n'est pas vrai dans le cas contraire.

B.4 Amélioration des interactions Homme-robot

Le troisième problème abordé est lié à l'estimation de l'intention de l'opérateur durant une tâche de collaboration. En effet, l'évolution de la robotique et les progrès récents dans le domaine de l'interaction physique Homme-robot tendent vers la conception orientée humaine de mécanismes robotiques et de stratégies de contrôle. Dans ce cadre, la caractérisation du comportement de l'opérateur humain est au centre de l'attention, l'interaction en vérifiant la modification de la dynamique d'interaction et en la modifiant si nécessaire. Les simulations présentées montrent que le second filtre permet de garantir la passivité pour des changements beaucoup plus rapides de la dynamique d'interaction. Les résultats expérimentaux, obtenus sur un montage conduisant à une instabilité sans filtre de passivité, montrent l'efficacité de la méthode proposée.

Avec la possibilité de changer la dynamique d'interaction, la compréhension des intentions de l'opérateur et l'inclusion de ces informations dans la stratégie de commande du robot pour une interaction physique Homme-robot améliorée passe au premier plan. Une stratégie basée sur l'électromyographie (EMG) est proposée pour distinguer les entrées de force de l'opérateur des interactions avec l'environnement et modifier la dynamique d'interaction du robot en conséquence. La validation expérimentale montre une bonne performance des modèles EMG-force proposés et illustre les avantages d'une telle stratégie dans les cas où l'Homme, le robot et l'environnement interagissent les uns avec les autres.

Enfin, deux preuves de concept d'applications, mettent en avant les outils de contrôle proposés dans le domaine médical. La première application médicale proposée présente un cobot effectuant une tâche d'insertion percutanée d'aiguille de menière autonome dans le cadre de la radiologie interventionnelle. L'insertion expérimentale est réalisée dans un fantôme anatomique émulant des tissus humains et affiche des résultats prometteurs en termes de précision de détection de couche et de capacité à adapter le comportement de l'aiguille. La deuxième application est développée dans le cadre d'une thérapie de rééducation bimanuelle pour les patients atteints de troubles de la fonction neuro-motrice. L'application proposée a été validée expérimentalement sur un sujet sain émulant une déficience.

Conference on Robotics and Automation (ICRA), pages 1808-1814. Cited on pages 42 and 43. [Chawda and Niemeyer, 2017b]