Model Based Automatic Tuning and Control of a Three Axis Camera Gimbal

Institutionen för systemteknik

Department of Electrical Engineering

Examensarbete

Model-based automatic tuning and control of a

three-axis camera gimbal

Examensarbete utfört i Reglerteknik vid Tekniska högskolan vid Linköpings universitet

av Henric Edlund LiTH-ISY-EX–15/4853–SE

Linköping 2015

Department of Electrical Engineering Linköpings tekniska högskola

Linköpings universitet Linköpings universitet

Model-based automatic tuning and control of a

three-axis camera gimbal

Examensarbete utfört i Reglerteknik

vid Tekniska högskolan vid Linköpings universitet

av

Henric Edlund LiTH-ISY-EX–15/4853–SE

Handledare: Daniel Simon

isy_{, Linköpings universitet} Mårten Svanfeldt

Intuitive Aerial Examinator: Martin Enqvist

isy, Linköpings universitet

Avdelning, Institution Division, Department

Automatic control

Department of Electrical Engineering SE-581 83 Linköping Datum Date 2015-06-11 Språk Language  Svenska/Swedish  Engelska/English   Rapporttyp Report category  Licentiatavhandling  Examensarbete  C-uppsats  D-uppsats  Övrig rapport  

URL för elektronisk version

http://urn.kb.se/resolve?urn=urn:nbn:se:liu:diva-XXXXX ISBN

— ISRN

LiTH-ISY-EX–15/4853–SE Serietitel och serienummer Title of series, numbering

ISSN —

Titel Title

Modellbaserad automatisk inställning och reglering av en treaxlig kameragimbal Model-based automatic tuning and control of a three-axis camera gimbal

Författare Author

Henric Edlund

Sammanfattning Abstract

En gimbal är en pivoterande anordning som frikopplar rörelser hos en plattform från dess last. Lasten är i detta fall en filmkamera som måste vara frikopplad från yttre rörelser för att kunna spela in videomaterial utan oönskade rörelser. En utmaning med en kameragim-bal är att den kan användas med många olika kombinationer av kameror och objektiv, som påverkar systemets dynamik. Detta innebär att systemet måste ställas in på nytt varje gång en ny kamerauppsättning ska användas. Då de flesta inom filmbranschen är obekanta med reglerteknik är denna justering ganska så problematisk.

Detta examensarbete undersöker vilka möjligheter som finns för att göra en automatisk inställning av ett reglersystem. Detta görs genom att excitera systemet och samla in data. Denna data används sedan för att beräkna en modell för systemet. Modellen används sedan för att reglera systemet och på så sätt behöver inte systemet justeras manuellt.

Som grund används en fysikalisk modell av systemet framtagen med hjälp av La-grangeekvationen. Denna modell innehåller okända parametrar som tröghetsmoment, tyn-gdpunktslägen och friktionskonstanter. För att beräkna dessa parametrar används metoder från systemidentifiering. Ett centralt problem vid systemidentifiering är hur systemet ska exciteras för att få fram så mycket information som möjligt om systemets dynamik. Detta problem angrips genom att formulera ett optimeringsproblem som används för att ta fram optimala trajektorior.

Den identifierade modellen används sedan för att reglera gimbalen. Olika metoder för att göra detta diskuteras. Genom att använda metoden återkopplingslinjärisering kan man kompensera för effekter som modelleras. Detta gör att systemet blir oberoende av modell-parameterarna och dessutom linjärt och frikopplat. För att styra det nya linjära systemet används en pid-regulator. Effekten av modellfel i återkopplingslinjäriseringen analyseras och för att garantera robusthet av det slutna systemet används metoden Lyapunov redesign. Slutligen presenteras experimentella resultat. Modellkvaliteten hos den skattade mod-ellen undersöks. Sedan används modmod-ellen för att reglera gimbalen och referensföljning tes-tas. Resultaten visar att det finns problem med den skattade modellens prestanda.

Nyckelord

Sammanfattning

En gimbal är en pivoterande anordning som frikopplar rörelser hos en plattform från dess last. Lasten är i detta fall en filmkamera som måste vara frikopplad från yttre rörelser för att kunna spela in videomaterial utan oönskade rörelser. En utmaning med en kameragimbal är att den kan användas med många olika kombinationer av kameror och objektiv, som påverkar systemets dynamik. Detta innebär att systemet måste ställas in på nytt varje gång en ny kamerauppsättning ska användas. Då de flesta inom filmbranschen är obekanta med reglerteknik är denna justering ganska så problematisk.

Detta examensarbete undersöker vilka möjligheter som finns för att göra en automatisk inställning av ett reglersystem. Detta görs genom att excitera syste-met och samla in data. Denna data används sedan för att beräkna en modell för systemet. Modellen används sedan för att reglera systemet och på så sätt behöver inte systemet justeras manuellt.

Som grund används en fysikalisk modell av systemet framtagen med hjälp av Lagrangeekvationen. Denna modell innehåller okända parametrar som trög-hetsmoment, tyngdpunktslägen och friktionskonstanter. För att beräkna dessa parametrar används metoder från systemidentifiering. Ett centralt problem vid systemidentifiering är hur systemet ska exciteras för att få fram så mycket infor-mation som möjligt om systemets dynamik. Detta problem angrips genom att formulera ett optimeringsproblem som används för att ta fram optimala trajekto-rior.

Den identifierade modellen används sedan för att reglera gimbalen. Olika me-toder för att göra detta diskuteras. Genom att använda metoden återkopplingslin-järisering kan man kompensera för effekter som modelleras. Detta gör att syste-met blir oberoende av modellparasyste-meterarna och dessutom linjärt och frikopplat. För att styra det nya linjära systemet används en pid-regulator. Effekten av mo-dellfel i återkopplingslinjäriseringen analyseras och för att garantera robusthet av det slutna systemet används metoden Lyapunov redesign.

Slutligen presenteras experimentella resultat. Modellkvaliteten hos den skat-tade modellen undersöks. Sedan används modellen för att reglera gimbalen och referensföljning testas. Resultaten visar att det finns problem med den skattade modellens prestanda.

Abstract

A gimbal is a pivoted device that decouples movements of a platform from its pay-load. The payload is a camera which must be stabilized to capture video without motion disturbances. A challenge with this type of gimbal is that a wide span of cameras with different sizes and weights can be used. The change of camera has significant effect on the dynamics of the gimbal and therefore the control system must be re-tuned. This tuning is inconvenient, especially for someone without knowledge of control engineering.

This thesis reviews suitable methods to perform an automatic controller tun-ing directly on the gimbal’s hardware. This tuntun-ing starts by excittun-ing the system and then using data to estimate a model. This model is then used to control the gimbal, thus removing the need for manual tuning of the system.

The foundation of this thesis is a physical model of the gimbal, derived thro-ugh the Lagrange equation. The physical model has undetermined parameters such as inertias, centre of gravity and friction constants. System identification is used to determine these parameters. A problem discussed is how the system should be excited in order to achieve data with as much information as possible about the dynamics. This problem is approached by formulating an optimization problem that can be used find suitable trajectories.

The identified model is then used to control the gimbal. Different methods for model based-control are discussed. By using a method called feedback linearisa-tion all of the parameter-dependant dynamics of the gimbal can be compensated for. Apart from being independent of model parameters the new outer system is also decoupled and linear. A pid controller is used for feedback control of the outer system. The uncertainty of the feedback linearisation is analysed to find the effects of model errors. To assure robustness of the closed loop system a Lyapunov redesign controller is used to compensate for these model errors.

Some experimental results are also presented. The quality of the estimated model is evaluated. Additionally, the reference tracking performance of the con-trol system is tested and results reveal issues with the estimated model’s perfor-mance.

Acknowledgments

I want to thank my supervisor Mårten Svanfeldt at Intuitive Aerial AB for his valuable contributions and for helping me sorting out various practical issues. I also want to thank my supervisor Daniel Simon and examiner Martin Enqvist for providing helpful feedback and support during this thesis work.

Linköping, June 2015 Henric Edlund

Contents

Notation xiii 1 Introduction 1 1.1 Background . . . 1 1.2 Purpose . . . 1 1.3 Thesis scope . . . 2 1.4 Approach . . . 2 1.5 Related work . . . 2 1.6 Thesis outline . . . 3 2 System overview 5 2.1 Introduction . . . 5 2.2 Gimbal construction . . . 5 2.3 Actuators . . . 5 2.4 Sensors . . . 5 2.5 Gimbal balancing . . . 6 3 Modelling 9 3.1 Introduction . . . 9

3.2 Mechanical multibody system theory . . . 10

3.2.1 Kinematics . . . 10

3.2.2 Differential kinematics . . . 11

3.2.3 Dynamics . . . 11

3.3 Assumptions for the gimbal . . . 13

3.4 Gimbal kinematics . . . 13

3.4.1 Coordinate representation . . . 14

3.4.2 Denavit-Hartenberg convention . . . 14

3.5 Gimbal differential kinematics . . . 16

3.5.1 Derivation of the Jacobian . . . 17

3.5.2 Jacobian for the gimbal system . . . 17

3.6 Gimbal dynamics . . . 18

3.6.1 Energy of the gimbal . . . 18

3.6.2 Lagrange equations for the gimbal . . . 20 ix

x Contents

3.6.3 Actuator modelling . . . 21

3.6.4 Friction modelling . . . 22

3.7 The resulting model for the gimbal system . . . 22

3.7.1 Complete model . . . 23

4 System identification 25 4.1 Introduction . . . 25

4.2 System identification theory . . . 26

4.2.1 Linear regression model . . . 26

4.2.2 Identification methods . . . 26

4.2.3 Identifiability and informative experiments . . . 28

4.2.4 Method for determining identifiable parameters . . . 29

4.2.5 Model structure for identification . . . 30

4.2.6 Validation . . . 31

4.3 Regression model for the gimbal . . . 31

4.3.1 Base parameters for the gimbal . . . 32

4.4 Experiment design . . . 33

4.4.1 Open loop versus closed loop identification . . . 33

4.4.2 Optimal trajectories . . . 34

4.4.3 Optimization of trajectories . . . 36

4.5 Identification for the gimbal . . . 37

4.5.1 Gimbal measurement equations . . . 37

4.5.2 Gimbal measurement filtering . . . 38

4.5.3 Gimbal predictor . . . 38

5 Control design 39 5.1 Introduction . . . 39

5.1.1 Model based control . . . 39

5.1.2 Proposed controller . . . 40

5.2 Control theory . . . 41

5.2.1 Feedback linearisation . . . 41

5.2.2 Lyapunov theory . . . 43

5.3 Feedback linearisation of the gimbal . . . 45

5.4 Outer loop control . . . 46

5.4.1 pid_{controller . . . . 46}

5.4.2 Uncertainty analysis of feedback linearisation . . . 47

5.4.3 Robust controller . . . 47 5.5 Controller tuning . . . 51 6 Experimental results 53 6.1 Introduction . . . 53 6.2 Experimental setup . . . 53 6.3 Identification results . . . 53 6.3.1 Raw measurements . . . 54 6.3.2 Data processing . . . 57 6.3.3 Estimated parameters . . . 59

Contents xi

6.3.4 Validation . . . 60

6.3.5 Conclusions . . . 65

6.4 Controller results . . . 65

6.4.1 Simultaneous axes step tracking . . . 66

6.4.2 Separate axis step tracking . . . 69

6.4.3 Conclusions . . . 71

7 Conclusions and future work 73 7.1 Results . . . 73

7.2 Discussion . . . 73

7.3 Conclusions . . . 74

7.4 Future work . . . 74

A Gimbal model 77 A.1 Components of the D matrix . . . 77

A.2 Components of the C matrix . . . 77

A.3 Gimbal regressor matrix . . . 79

Notation

Notation

Symbol Explanation

Rj_i Rotation of frame i expressed in frame j.

T_ij Homogenous transformation, rotates and translates frame i expressed in frame j.

ri/j Vector originating in point i ending in j

ω_ij Angular velocity of body i expressed in frame j Γ The model error vector of the feedback linearisation

q Gimbal axis angle vector, contains axis angles q1, q2 and q3

Abbrevations

Symbol Explanation

ia _{Intuitive Aerial}

uav _{Unmanned aerial vehicle} dh Denavit-Hartenberg convention imu Inertial measurement unit

ahrs Attitude and heading reference system. Uses the imu to compute attitude and heading measurements.

1

Introduction

1.1

Background

This thesis was performed at the company Intuitive Aerial (ia), which is special-ized at producing flying camera rigs for professional aerial photography. Their product consists of a multi-rotor platform and a gimbal that carries the camera. The gimbal is a pivoted device that allows the camera to rotate freely around its three axes. Since the gimbal is attached to a moving platform it is crucial to isolate the camera from the movements of the platform. It is also desirable to control the camera’s line of sight. This is accomplished by using electric motors and inertial sensors. The gimbal is constructed to handle cameras weighing up to 9 kilograms. Therefore, cameras with a wide span in weight and inertial prop-erties can be used with the gimbal. This poses a challenge for the control system since its performance is dependant on the dynamics of the whole system. The current control system uses pid controllers and requires cumbersome tuning for good performance. Another difficulty with the current system is the balancing of the gimbal. In order to have good controller performance the system must be carefully balanced so that the centre of mass is aligned with the rotational centre of the gimbal.

1.2

Purpose

The purpose of this thesis is to overcome the above mentioned issues by develop-ing an automatic online tundevelop-ing of the control system. This is achieved by derivdevelop-ing a mathematical model of the gimbal’s dynamics and then using methods from system identification to find the parameters of the model. Finally, this identified model can be used with a model-based control method which does not require the user to perform any tuning. Additionally, this type of controller can improve

2 1 Introduction

the overall performance, compared to a pid-based system.

1.3

Thesis scope

The system uses sensors such as gyroscopes, accelerometers, magnetometers and angular sensors at the joint axes. The processing and fusion of these signals will not be covered in this thesis. Instead, the solutions developed by ia will be used. A fundamental part of the control system is the input preprocessing which ma-nipulates the camera operator’s control inputs to smooth reference trajectories for the gimbal. This system will not be considered in this thesis.

To focus on the theoretical work some practicalities, such as the actual online tuning, are not implemented on the hardware. It is instead performed offline and the parameters are then programmed in the gimbal’s hardware. Relations for transforming sensor and reference signals into the signals required by the model have been derived but not implemented. Instead the gimbal is mounted to a rigid platform which reduces the complexity of sensor transformations for identification and control.

1.4

Approach

The base of this thesis is a physical model derived with methods from multi-body mechanics. This physical model has unknown parameters that describe the dy-namic properties of the camera and the gimbal. To find estimates of these pa-rameters, methods from system identification have been used. The performance of the controller is then evaluated on the gimbal. For ease of use the whole au-tomatic tuning procedure should be handled by the gimbal system. This means that all data collection, data preprocessing, estimation and calculations related to controller synthesis should be possible to run on the gimbal’s hardware.

1.5

Related work

Previous work on gimbals is not very common, especially system identification. Johansson [2012] is a previous master thesis on a three-axis gimbal performed at ia, where a model for the gimbal is derived and a pid controller is used. In Skoglar [2002], a model is derived and then identification and control methods are evaluated through simulations. Li et al. [2013] covers system identification of a two-axis gimbal. It results in a non-linear, demanding optimization problem which may not be suitable for this thesis. In Sangveraphunsiri and Malithong [2009] a robust feedback linearisation controller is evaluated as well as a sliding mode controller. The sliding mode controller is also evaluated in Shtessel [1994]. Kwon et al. [2007] evaluates a decentralized and linear controller for a gimbal system.

The structure of a gimbal closely resembles a robot. The field of robotics has several established approaches for modelling, identification and control.

There-1.6 Thesis outline 3

fore, many of the methods that have been used for robotics should be appropriate to apply to the gimbal.

1.6

Thesis outline

First a system overview is given in Chapter 2. This chapter is then followed by the three main parts, modeling in Chapter 3, identification in Chapter 4 and control in Chapter 5. Each main part contains a theory part going over the fundamental theory which is then applied to the gimbal. In Chapter 6, the experimental re-sults are presented and in Chapter 7 conclusions and a discussion concludes this thesis.

2

System overview

2.1

Introduction

This chapter’s purpose is to give a short overview of how the gimbal works. It briefly covers construction, sensors and actuators of the gimbal.

2.2

Gimbal construction

The mechanical structure of the gimbal consists of four parts. On the top there is the head that connects the gimbal to the multirotor platform. The connection consists of damper springs to remove vibrations. Connected to the head is the spine of the gimbal, which is also connected to the legs, see Figure 2.1. Then there is the fork and finally the cradle, where the camera is mounted. The body of the gimbal is built of carbon fibre and aluminium.

2.3

Actuators

Each axis of the gimbal is driven by an electric motor. The motors are brushless direct drive motors with high torque output. The motors are controlled by a motor control system which uses angular sensors to measure the angle of the motor. The motors are torque controlled which means that the controller takes a reference torque as input signal.

2.4

Sensors

The gimbal has two imus, one situated at the spine and one at the cradle. The imu_{s’ computer uses an attitude and heading reference system (ahrs) algorithm}

6 2 System overview

Head

Spine

Fork

Cradle _{Lower imu} Upper imu

Roll motor

Tilt motor Yaw motor

Camera mount

Figure 2.1:A schematic of the gimbal where the parts of the gimbal are pre-sented. The positions of actuators and sensors are pointed out. Illustration reproduced by courtesy of Intuitive Aerial

to obtain both attitude angle and angular velocity. The gimbal is also equipped with angular sensors at each joint. These are used in the motor controllers but are also available for use in the control unit.

2.5

Gimbal balancing

The purpose of balancing the gimbal is to adjust the positioning of the camera such that the whole system’s centre of mass is aligned with the gimbal’s rotational axes. If a perfect balancing is achieved, gravity will not interfere with the

gim-2.5 Gimbal balancing 7

bal’s dynamics thus making it easier to control. If the gimbal is not balanced the motors must compensate for static moments induced by gravity. This degrades the performance of the pid regulators. A poor balancing may result in gravity moments higher than the motors’ maximum output torque making it impossible to control the gimbal. For instance, a gimbal without a camera mounted is so badly balanced that it becomes difficult for the motors to move the gimbal.

3

Modelling

3.1

Introduction

For system identification and model-based control methods it can be advanta-geous to have a physical model of the gimbal. Since the gimbal is a mechanical system, methods used in mechanics should be a good starting point. There are many more or less systematic modelling approaches available in engineering me-chanics which can be used for derivation of a mathematical model describing the motion of the parts in the gimbal.

A camera gimbal closely resembles an industrial robot, which has been stud-ied comprehensively. The literature suggests two different modelling approaches. One is the Lagrange method and the other is the Newton-Euler method [Siciliano et al., 2009, Ch. 7]. The principal difference is that the Lagrange method is based on energy balances in the system and Newton-Euler’s method is based on balance laws of forces and torques for the system.

Before any of these two methods can be used some basic mechanical relations must be derived. These are kinematics and differential kinematics. Kinemat-ics describes the geometrical relations of the system. Differential kinematKinemat-ics is based on these relations and describes the mapping between joint velocities and link velocities. To determine the kinematical relations it is convenient to assign a coordinate system to each body of the gimbal. The coordinate assignments are made systematically using a convention called Denavit-Hartenberg [Siciliano et al., 2009, Ch. 2]. The derivation of the gimbal model is based on general n-link robot approach presented in Siciliano et al. [2009].

This approach for modelling a gimbal has been tried before in Johansson [2012]. The particular model deduced in this thesis had the problem of being difficult to work with due to its size. The main reason for this was the coordinate frame assignment which was mechanically correct but not made in the easiest

10 3 Modelling way. In this thesis, a different coordinate system assignment is used which re-sults in more simple dynamic equations.

3.2

Mechanical multibody system theory

The three fundamental parts of mechanical modelling are kinematics, differential kinematics and dynamics, as was mentioned in the introduction. We will now briefly go over the theoretical background of these.

3.2.1

Kinematics

To be able to model a mechanical system it is essential to describe the geometry of the system. Kinematics offers us a tool to describe this geometry. A comprehen-sive discussion on kinematics can be found in Siciliano et al. [2009, Ch. 2]. The geometry can be divided into rotations and translations that describe the bodies of the system. To represent the rotations, rotation matrices are useful. A rota-tion matrix relates a coordinate between two frames that have a relative rotarota-tion. The rotation in Figure 3.1 is a rotation around the z0axis by the angle θ and the resulting rotation matrix is

R0₁=         cos θ −_{sin θ 0} sin θ cos θ 0 0 0 1         (3.1)

y

0

x

0

x

1

y

1

θ

z

0

, z

1

Figure 3.1:Two coordinate frames related by a rotation matrix.

A useful property of rotation matrices is their orthogonality property. This makes it possible to invert a rotation matrix by transposing it.

To have a complete description of the relation between two coordinate sys-tems the translation must be included. This is simply done by defining the vector

3.2 Mechanical multibody system theory 11

Let us look at how rotations and translations can be used to describe points in the frames. Let p0 and p1 be the coordinate of a point p in frame 0 and 1, respectively. By using the previously established rotation and translation we can write

p0= d0/1+ R01p1 (3.2)

A more compact way to represent coordinates in different frames is to use ho-mogeneous transformations. The homegenous coordinate for an arbitrary point

p is denoted as ˜ p ="p 1 # (3.3) The complete coordinate transformation with both the rotation and translation combined can be written as

T10= " R0₁ d₁0 0 1 # (3.4) The coordinate transformation in (3.2) can now be written as

˜

p0= T₁0p˜1 (3.5)

Coordinate transformations between a sequence of n frames can be written as ˜ p0= n Y k=1 T_kk−1p˜n (3.6)

This is called a kinematic chain and is useful for describing a system’s geometry.

3.2.2

Differential kinematics

We would also like to have a description of the velocities between different frames. By using rotation matrices, velocity vectors can be related to other frames as

ω_1/00 = R0₁ω_1/01 (3.7)

Furthermore, once velocities are expressed in the same frame they can be added. These properties can be used to determine a relation between the velocities of the bodies in a system. This relation is called a Jacobian and will be determined specifically for the gimbal in Section 3.5. Differential kinematics is treated in Siciliano et al. [2009, Ch. 3].

3.2.3

Dynamics

The goal of dynamics is to find the differential equations that determine the rela-tion between inputs and system states. One method to determine the dynamics of a system is the Lagrange method. It is based on knowledge of a system’s total energy in the form of a function known as the Lagrangian which is defined as

12 3 Modelling

where T denotes the total kinetic energy and U is the total potential energy. By using the Lagrange equation [Siciliano et al., 2009]

d dt ∂L ∂ ˙q ! −∂L ∂q = Q (3.9)

where Q denotes the external force and q is the system state, the dynamics of the system is obtained.

Let us consider a simple example to illustrate this method 3.1 Example: Pendulum

θ

ℓ

p

Figure 3.2:A rigid pendulum seen from the side.

We want to model the pendulum in Figure 3.2. The pendulum has mass m and inertia I. Its centre of gravity is located at p, at the distance ` from the pivot point, which is fixated in an inertial frame. The expressions of both kinetic and potential energy must be determined first. Since the pendulum can be modelled as a rigid link, one contribution is the rotational energy, the other is the centre of mass translational energy. Hence the total kinetic energy ca be written

T ₌ 1 2I ˙θ

2₊1

2m ˙p

2 _(3.10)

The potential energy is

U _{= mg`(1 − cos(θ)) (3.11)}

To use the Lagrange method we must choose a generalized coordinate, and a simple choice is θ. Then ˙p must be expressed as a function of ˙θ. By using the

differential kinematics relation

3.3 Assumptions for the gimbal 13

we can form the Lagrangian as L₌ 1

2I ˙θ

2₊1

2m`

2_θ˙2_{−mg`(1 − cos(θ)) (3.13)} By using the Lagrange equation we have

(I + m`2) ¨θ + mgl sin(θ) = 0 (3.14)

If we would have used the Newton-Euler method instead, we would then have had to look at the forces and torques of the system to obtain the equations of motion. In this case, there is not much difference in complexity of the modelling procedures of these two methods. But as the system complexity is increased it is generally easier to compute the energy of the system rather than to calculate the torques and forces in the system. Therefore the Lagrange method has a clear advantage. A complete description of both the Lagrange method and the Newton-Euler method can be found in Siciliano et al. [2009, Ch. 7].

3.3

Assumptions for the gimbal

A fundamental dilemma when modelling a physical system is the trade-off be-tween simplicity and model accuracy. We want the model to be simple enough to be used with established methods, yet it is crucial that the model can describe the mechanical behaviour well enough. In order to achieve a balanced trade-off some assumptions are made for the gimbal.

For certain loads and velocities, the links of the system will start to deform elastically and this causes oscillations in the system. This oscillative behaviour is assumed to be insignificant compared to the overall dynamic behaviour of the gimbal. When doing the calibration it will be important not to excite these flex-ibilities for a good model estimation. Furthermore, the controller must be so robust that these flexibilities can be handled. The main reason we do not want to model these flexibilities is the increase in model complexity.

When modelling friction there are many different aspects to consider. To ob-tain a useful model many of these effects are neglected.

Since the gimbal will be attached to a multi-rotor platform some assumptions must be made about the interaction between these two. It is assumed that the movement of the gimbal does not affect the multi-rotor platform’s movement. The movement of the multi-rotor platform, however, affects the gimbal and can-not be neglected.

3.4

Gimbal kinematics

The geometry of the gimbal can be seen as the result of a sequence of rotations and (constant) translations between a fixed base and the camera mounted in the gimbal.

14 3 Modelling

3.4.1

Coordinate representation

When working with transformation matrices one has to consider how they are represented. A commonly used coordinate representation is Euler angles; it uses three parameters to describe a rotation. A drawback is the prevalence of singu-lar points. For the natural rotation order of the gimbal these singusingu-larities are located at ± 90◦

in the roll axis. These points will never be reached in any oper-ating scenario as they lie outside of the gimbal’s workspace, hence this causes no problem.

An alternative to Euler angles is quaternions. They consist of a non-minimal representation with 4 parameters and have no singularities. One the downside they lead to more complex expressions and are not as intuitive to work with. Since there is no benefit with using quaternions in this case we chose to work with Euler angles. A more comprehensive treatment of these coordinate repre-sentations can be found in Siciliano et al. [2009, Ch. 2]

3.4.2

Denavit-Hartenberg convention

The Denavit-Hartenberg (dh) convention is a procedure to assign coordinate frames to the links of a robot. The gimbal consists of 4 bodies and thus 4 co-ordinate systems are needed for a complete description of the geometry of the system. The structure can be seen in Figure 3.3. Body (0) is the base of the gim-bal which is attached to the multirotor platform. Body (1) is called the spine and connects body (0) to body (2), which is called the fork. The fork then connects to body (3) which is called the cradle, where the camera is mounted. The DH convention is described in Siciliano et al. [2009, Ch. 2] and can be used for a gen-eral link structure. The gimbal’s structure results in some simplifications when applying the DH convention, therefore a modified version is used.

Modified Denavit-Hartenberg Algoritm

1. Find the joint axes and set the corresponding z0, .., znaxes.

2. Choose frame 0 by locating the origin along the z0-axis. Axes x0and y0can be chosen arbitrarily as long as a right-handed frame is achieved.

Go through Steps 3 to 5 for i = 1, ..., n − 1

3. Locate origin 0iat the intersection of ziand zi−1

4. Assign axis xi along the common normal to axes zi and zi−1with direction from joint i to joint i + 1.

5. Choose axis yi to obtain a right-handed frame. Final steps

6. Choose the last frame n. Align znwith zn−1and choose xnaccording to step 4 and ynto obtain a right-handed system.

3.4 Gimbal kinematics 15

Denavit-Hartenberg applied to the gimbal

In Figure 3.3 a schematic view of the gimbal is presented with the coordinate systems from the dh algorithm. The rotational angles are represented by θ1, θ2 and θ3 and the axis of rotation coincides with body fixed coordinate axis z0, z1 and z2. Y0 X0 Z0 Y1 X3 Z3 Y2Z1 X2 X1 Y3 Z2

Figure 3.3:The gimbal with coordi-nate systems seen from the side.

Y0 X0 Z0 Y1 X3 Z3 Y2 Z1 X2 X1 Y3 Z2

Figure 3.4: Side view of the gimbal where coordinate systems with overlapping origins have been moved for a clearer view.

Transformation matrices for the gimbal

Now we want to derive the transformation matrices between the coordinate sys-tems defined above. We start with frame 0 and go step by step until frame 3 is reached. The rotation matrix that describes frame 1 relative to frame 0 is

R01=         −_{sin θ1 0 cos θ1} cos θ1 0 sin θ1 0 1 0         (3.15)

The translation vector is

d₁0=         0 0 h         (3.16) And finally the transformation matrix is obtained as

T₁0= " R0₁ d₁0 0 1 # (3.17)

16 3 Modelling

The rotation matrix that describes frame 2 relative to frame 1 is

R1₂=         −_{sin θ2 0 cos θ2} cos θ2 0 sin θ2 0 1 0         (3.18)

The origins of the two frames coincide, thus

d1₂=         0 0 0         (3.19)

Furthermore, the transformation matrix becomes

T₂1= " R1₂ d1₂ 0 1 # (3.20) The rotation matrix of frame 3 relative to frame 2 is

R23=         −_{sin θ3 0 cos θ3} cos θ3 0 sin θ3 0 1 0         (3.21)

Again the origin of the two frames coincides, and

d23=         0 0 0         (3.22) T₃2= " R2₃ d2₃ 0 1 # (3.23) Finally, the transformation between frame 0 and fram 3 is

T₃0=             cθ1cθ3−sθ1sθ2sθ3 −cθ2sθ1 cθ1sθ3+ cθ3sθ1sθ2 0 cθ3sθ1+ cθ1sθ2sθ3 cθ1cθ2 sθ1sθ3−cθ1cθ3sθ2 0 −_cθ 2sθ3 sθ2 cθ2cθ3 h 0 0 0 1             (3.24)

where cθi = cos(θi) and sθi = sin(θi).

3.5

Gimbal differential kinematics

To determine the equations of motion of the gimbal it is necessary to know how joint velocities map to velocities at the links of the gimbal. This mapping is called the Jacobian and it is a linear mapping of joint velocities to link velocities.

3.5 Gimbal differential kinematics 17

3.5.1

Derivation of the Jacobian

As mentioned was earlier, the gimbal motion consists of only rotations since all rotation axes intersect in a common point. Thus, when deriving the differential kinematics the parts concerning translational velocity are left out. This deriva-tion can be found in Siciliano et al. [2009, Ch. 3].

The rotation Jacobian is formulated as

ω_ij = J_ijθ˙ (3.25)

and is a mapping between the joint velocity ˙θ and the body velocity ω_ij. The vec-tor ω_ij denotes the angular velocity of body i, represented in coordinate system

j.

For a single joint it can be seen that

ωi−1_i/i−1= ˙θizi−1 (3.26) where zi−1 is the previous joint’s rotational axis. The vector zi−1 can now be represented in the base frame by using the kinematic chain (3.6) as

zi−1= i Y

k=1

Rk_k−1z0 (3.27)

Since the rotational axes are represented in the same frame, the angular velocities can be added. For each joint k, the corresponding column of the Jacobian becomes

Jk=        zk−1 k ≤ i 0 otherwise (3.28) and Jn0= h Jk Jk+1 ... i for each k = 1, ..., n (3.29)

3.5.2

Jacobian for the gimbal system

Using the briefly presented theory above, the Jacobian is

J₃0=hz2 z1 z0 i =hR0₁R1₂z0 R01z0 z0 i (3.30) where z0= h 0 0 1iT (3.31)

Note that the Jacobian is represented in the base frame. Then the Jacobian for each individual link becomes

J₁0=         0 0 0 0 0 0 1 0 0         , J₂0=         0 cos θ1 0 0 sin θ1 0 1 0 0         (3.32) and J₃0=        

0 cos θ1 −cos θ2sin θ1 0 sin θ1 cos θ1cos θ2

1 0 sin θ2         (3.33)

18 3 Modelling

3.6

Gimbal dynamics

The goal of the dynamic equations is to obtain a model that describes how the joint states depend on the actuators. As was mentioned in the introduction, there are two well-regarded methods that can be used; the Lagrange method and the Newton-Euler method. These methods are based on different principles and re-sult in different types of expressions for the dynamics. The Lagrange method results in what is known as the inverse dynamics form. This means that it gives the actuator torque as a function of the joint variables and their differentials. This form can easily be written as a state-space model.

The Newton-Euler method results in a recursive formulation of the system’s dynamics. It is possible rewrite either method’s result to the same form as the other. Since the Lagrange method requires a bit less work and gives the model on the desired form directly, it is chosen. The dynamics have been derived based on Siciliano et al. [2009, Ch. 7]. The Lagrange method is independent of the chosen coordinate frames, but it depends on joint variables qi, which must be chosen. The natural choice is the joint angles θi, thus qi = θi. To use the Lagrange method we need to find expressions for the energy of the gimbal.

3.6.1

Energy of the gimbal

Kinetic energy

The total kinetic energy of the system is T ₌ n X k=1 T_{i (3.34)} where T_{i =} 1 2 Z V_i ˙ pT_i p˙iρdV (3.35)

is the energy for link i and ˙pi is the velocity of volume element dV for link i. Let us define the vector

ri/c = pi− pc (3.36)

where pcis the rotational centre of the gimbal. The relation that describes the velocities of two points in a rigid body as it rotates is

˙

pi = ˙pc+ S(ωi)ri/c (3.37)

where S is a skew-symmetric matrix. By using this relation we have T_{i =} 1

2 Z

Vi

3.6 Gimbal dynamics 19

In Section 3.5 it was determined that ˙pc = 0 for all joints and therefore that the energy can be expressed as

T_{i =}1 2 Z Vi (S(ω_ii)ri/c)T(S(ωii)ri/c)ρdV =1 2 Z Vi (−S(ri/c)ωii)T(−S(ri/c)ωii)ρdV =1 2(ω i i)T Z Vi S(ri/c)TS(ri/c)ρdV ! ω_ii = 1 2(ω i i)TIiωii (3.39)

where we have defined the matrix

Ii =           R (r_iy2 + r_iz2)ρdV −R_rixriyρdV −R _rixrizρdV . . . R(r_ix2 + r2_iz)ρdV −R _riyrizρdV . . . . . . R(r_ix2 + r2_iy)ρdV           (3.40) =        

Iixx −Iixy −Iixz

. . . Iiyy −Iiyz . . . . . . Izz         (3.41)

Note that the matrix is symmetric and that the ". . ." elements are obtained from symmetry. Now the kinetic energy is expressed as a function of link velocity in the link’s own coordinate system. To express the kinetic energy as a function of joint variables two steps must be taken. First the velocity of each link must be expressed in the base frame. By using

ω_ii = (R0_i)Tω0_i (3.42) the expression T_{i =} 1 2(ω 0 i)TR0iIi(R0i)Tωi0 (3.43) for the kinetic energy is obtained. The next step is to express the kinetic energy as a function of the joint variables. By using the Jacobian, ω_i0can be replaced by joint variables through the following expression

ω_i0= Ji(q) ˙q (3.44)

The total kinetic energy becomes T ₌1 2 n X i=1 n X j=1 dij˙qi˙qj= 1 2q˙ T_{D(q) ˙q (3.45)} where D(q) = n X i=1 Ji(q)TR0iIi(R0i)TJi(q) (3.46) and D(q) is symmetric and positive definite [Siciliano et al., 2009, Ch. 7].

20 3 Modelling

Potential energy

The potential energy is calculated for each link as U ₌ n X i=1 U_{i (3.47)} where U_{i = −} Z Vi g₀Tpi,cρdV = −migT0rci/0 (3.48) and rci/0is the vector between the centre of gravity for link i and the base frame. Let us introduce rci/i =         `ix `iy `iz         (3.49) Then rci/0= Ti0rci/i (3.50)

where T_i0is the homogeneous transformation between frame i and 0. It is deter-mined by using the kinematic chain described by (3.6).

3.6.2

Lagrange equations for the gimbal

When the kinetic and potential energy have been computed the Lagrangian be-comes

L_{(q, ˙q) = T (q, ˙q) − U (q) (3.51)} By considering the differentiations of L according to (3.9) we get

∂L ∂ ˙qi = n X j=1 dij˙qj (3.52)

Note that U (q) does not depend on ˙q. Then we have

d dt ∂L ∂ ˙qi ! = n X j=1 dijq¨i + n X j=1 d dtdij(q) ˙qj = n X j=1 dijq¨i+ n X j=1 n X k=1 ∂dij(q) ∂qk ˙qj˙qk (3.53) and finally ∂L ∂qi = 1 2 n X j=1 n X k=1 ∂djk(q) ∂qi ˙qj˙qk− ∂U ∂qi |{z} gi (3.54)

By putting it all together, we get n X j=1 dijq¨i + n X j=1 n X k=1 hijk˙qj˙qk+ gi(q) = Qi for i = 1, ..., n (3.55)

3.6 Gimbal dynamics 21 where hijk= ∂dij ∂qk − 1 2 ∂djk ∂qi (3.56) and cij = n X k=1 hijk˙qk (3.57)

To summarize, (3.55) can be written on compact matrix form as

D(q) ¨q + C(q, ˙q) ˙q + G(q) = Q (3.58)

Physical interpretation

The Lagrange method may not be as physically intuitive as, for instance, Newton-Euler’s equations. The following description gives a physical interpretation of what the model represents.

The matrix D(q), with elements dij, describes inertial effects of the gimbal. When i = j the element describes the inertia of that particular axis. When i 6= j the element describes inertial effect on axis i by axis j.

Next we have the matrix C(q, ˙q), with elementscij. According to (3.57), these elements can be written as hijk. The element hijj describes the centrifugal ef-fects induced on joint i by joint j. The elements hijkrepresent the Coriolis effect induced on joint i by joints j and k.

Finally, the element gi in the vector G(q) describes the moment induced by gravity on link i.

3.6.3

Actuator modelling

The actuators consist of brushless motors. The torque τ is proportional to the mo-tor’s current, according to τ = Kii. The motor’s torque is not linear but saturated at some τmax. Furthermore, the current in a brushless electric motor is described by

Ld

dti + Ri + V − Vemf = 0 (3.59)

where L is the inductance and R is the resistance of the motor. Vemf = Kbω is the counter electromotive force which is proportional the velocity. The inductance of the motor is omitted by letting L = 0, thus making the motor model static. Furthermore, since the motor is directly coupled to the gimbal axis, without a gearbox, the rotational speed of the motor is going to be very low. The term Vemf is thus insignificant and removed from the equations. The torque then becomes

τ =        KiV R V < Vmax τmax V > Vmax (3.60) The actuator system is considered as a separate system and is therefore not in-cluded explicitly in the gimbal model. Therefore we consider the gimbal system input to be the actuator torque.

22 3 Modelling

3.6.4

Friction modelling

Friction is difficult to model and therefore some friction effects are omitted. By using a very simplified model, which only accounts for static and viscous friction, the complexity can be reduced. This model is expressed as

Fi = Fissign( ˙qi) + Fivq˙i (3.61) where Fisis the static friction parameter and Fivis the viscous friction parameter. This model is linear with respect to the parameters which is useful when estimat-ing friction together with the rest of the dynamics. A drawback with this model is that it is only valid for velocities over a certain threshold.

3.7

The resulting model for the gimbal system

So far, the equations have just been derived from theory with only a few assump-tions made about the gimbal system. The results from kinematics, differential kinematics and dynamics will now be put together to get the model for the gim-bal.

To obtain a more usable model, the inertia matrices are assumed to be diago-nal. This simplification may not be needed for identification and control, as the number of terms in the model should not be a problem for the hardware of the gimbal. This will be tested and discussed in detail later. Note that some elements of the inertia matrix are in fact zero due to geometric symmetries. Along with the model equations the parameters of the model will be presented as well.

Model parameters

In Table 3.1 the parameters of the model are summarized. Note that some pre-viously defined parameters have no effect due to the restriction of motion. For instance body 1 only has one degree of freedom and therefore only one inertial parameter is needed. Furthermore, the motion of body 1 is not affected by gravity and therefore no gravity parameters are needed to describe it’s motion. Body 2 has two degrees of freedom are therefore two inertial parameters are needed as well as two gravity parameters.

Table 3.1:Table with the model’s parameters.

D I1yy I2xxI2yy I2zz I3xx I3yy I3zz

C I1yy I2xxI2yy I2zz I3xx I3yy I3zz

G `2mx`2mz`3mx`3my `3mz

F F1sF1vF2sF2v F3sF3v

Calculating D

The matrix D is obtained from (3.46) as

3.7 The resulting model for the gimbal system 23

All elements of the matrix are presented in Appendix A.

Calculating C

The matrix C(q, ˙q) is obtained from D(q), by using (3.56) and then (3.57). All elements of the matrix are presented in Appendix A.

Calculating G

By using the relations from Section 3.6.1, G (q) is obtained. The notation is simpli-fied by letting one variable, l∗_m∗, describe the product of length, mass and gravity. Each element describes the gravity effect och the corresponding gimbal axis.

The first element of G is

g1= 0 (3.63)

Furthermore, for the second axis we have the following expression

g2= l3my cos(q2) + l2mz cos(q2) − l2mxsin(q2) −

l3mz cos(q3) sin(q2) + l3mx sin(q2) sin(q3)

(3.64) Finally, for the third axis

g3= − cos(q2) (l3mxcos(q3) + l3mz sin(q3)) (3.65)

Calculating Q

The modelled non-conservative torques are actuator torque and friction torque, respectively

Q_{= τ − F ( ˙q) (3.66)}

with the first element

Q_{1= τ1}−_{F1ssign( ˙q1) − F1v˙q1 (3.67)} Furthermore, the second element is

Q_{2= τ2}−_{F2ssign( ˙q2) − F2v˙q2 (3.68)} and the third element is

Q_{3= τ3}−_{F3ssign( ˙q3) − F3v˙q3 (3.69)}

3.7.1

Complete model

Finally, the complete model has the form

D(q) ¨q + C(q, ˙q) ˙q + G(q) + F ( ˙q) = τ (3.70) This structure is called inverse dynamic form since it gives a static relation be-tween the dynamic variables q, ˙q and ¨q and the torque τ. This structure will be

24 3 Modelling

It is interesting to be able to simulate the gimbal model. For simulation we want to know how a torque input affects the motion of the gimbal. This can be done by writing the model on state-space form. By solving the equation for ¨q we

get

¨

q = D(q)−1(τ − C(q, ˙q) ˙q − G(q) − F ( ˙q)) (3.71) and by introducing the state vectors

xθ = q (3.72)

xω = ˙q (3.73)

the dynamic equations can be written on state-space form as " ˙xθ ˙ xω # = " xω D(xθ)−1τ − C(xθ, xω)xω − G(xθ) − F(xω) # (3.74)

4

System identification

4.1

Introduction

The physical model derived in Chapter 3 has parameters which depend on the gimbal and the mounted camera. We want our model to describe the actual sys-tem and therefore we must find these parameter values. The accuracy of the model will largely depend on how accurately we can determine these parameters. The performance and robustness of model based controllers will also depend on the quality of the parameters.

One way of determining some of the required parameters is to use a cad-model of the gimbal to find the centre of gravity and inertia for the gimbal’s parts. Parameters as friction cannot be found by this method and the inertia of the camera varies depending on what camera is used. Therefore, we cannot rely on this method, but it can still be used to find reference values for the gimbal-related parameters.

Another method is to use system identification to find these parameters. The basis of system identification is to observe the input/output behaviour and use this data to estimate models. Estimates of the parameters can be found by min-imizing the difference between the model output and the measured output, and this is called the prediction-error method. System identification has three funda-mental aspects that must be considered.

The first aspect is to consider is the structure of the model. This concerns both the input/output dynamics and also the modelling of disturbances in the system. When a model structure has been chosen we must consider how to compute es-timates for this particular structure. Depending on properties such as linearity and disturbance modelling different estimation methods must be used.

Finally, we must consider how we acquire and use the data for estimation. For example, what sensors should be used, and how should the data be preprocessed.

26 4 System identification

A major question is how the system should be excited to acquire data with as much information as possible of the system’s dynamics. This area has been stud-ied extensively for robots and there exist a few books. One such book is Khalil and Dombre [2002] where an approach is presented for identification of robots. A good property with the model structure called inverse dynamics form, obtained in Chapter 3, is the possibility to rewrite the model as a linear regression. Before the regression model can be used to compute estimated parameters we must find out which parameters can be identified.

4.2

System identification theory

Here, the theory for general regression models as well as a method for determin-ing the identifiability of regression models will be presented. For a complete review of regression theory, see Ljung [1987].

4.2.1

Linear regression model

A linear regression is a model structure which is linear with respect to its param-eters and can be written

y(t) = ϕT(t)θ + e(t) (4.1) where, in the particular case of an arx model,

ϕ(t) =h−_{y(t − 1) . . . y(t − na) u(t) . . . u(t − nb)}iT _(4.2) is the regression vector, y is the output measurement, u is the input and θ is the parameter. The parameter na is the number of poles and nb is the number of zeros.

4.2.2

Identification methods

The prediction

ˆ

y (t|θ) = ϕT(t)θ (4.3) of y(t) makes it possible to define the prediction error

 (t, θ) = y(t) − ˆy (t|θ) (4.4) Finding the θ that results in the smallest possible prediction error is a common estimation objective. One approach is the least squares method.

Least squares

The prediction error can be used to define an optimization criterion for estimat-ing θ. One such criterion is the sum of squared prediction errors

VN(θ) = 1 N N X t=1 1 2 2_{(t, θ) (4.5)}

4.2 System identification theory 27

and by minimizing this sum with respect to θ, we get the parameter estimate ˆ

θN = arg min

θ VN (4.6)

We thus have the least squares method. The criterion in (4.6) can be minimized analytically by substituting (4.3) and (4.4) into (4.5), which gives

VN(θ) = 1 N N X t=1 1 2(y(t) − ˆy (t|θ)) 2 = 1 N N X t=1 1 2  y (t)2−_2θT_{ϕ(t)y(t) + θ}T_{ϕ ϕ}T_θ _(4.7) By identifying fN = 1 N N X t=1 ϕ(t)y(t) (4.8) RN = 1 N N X t=1 ϕ(t)ϕT(t) (4.9) (4.7) can now be written as

VN(θ) = 1 N N X t=1 1 2  y(t)2− fNTR −₁ NfN + (θ − R −₁ NfN)TRN(θ − R −₁ NfN)  (4.10) The first two terms are independent of θ and will not affect the minimum. Since

RN is assumed positive semi-definite, the last term will always be positive and the minimum is found at

ˆ

θN = R

−₁

NfN (4.11)

The covariance of the estimate is calculated as

PN = λ0        N X t=1 ϕ(t)ϕT(t)        −₁ (4.12) where λ0 is the variance of the noise term e(t). This variance can be estimated with λ0= 1 N − nθ N X t=1 (t, ˆθN)2 (4.13) where N is the number of samples and nθis the length of θ.

28 4 System identification

Matrix formulation

The summations can be substituted with matrices for a more compact represen-tation. The measurements can be written as

YN =           y(1) .. . y(N )           (4.14)

and the regression matrix as

ΦN =            ϕT(1) .. . ϕT(N )            (4.15)

The estimate becomes

ˆ

θN = (ΦNTΦN)

−₁

ΦNTYN (4.16)

Multivariable case

The multivariable case is obtained by a slight modification of the one-dimensional case. The measurement is

y(t) =           y1(t) .. . yn(t)           (4.17)

and the regression matrix is

ϕT(t) =            ϕ₁T(t) .. . ϕ_nT(t)            (4.18)

for both the matrix and summation forms. The target function becomes

VN(θ) = 1 N N X t=1 1 2k(t, θ)k 2 _(4.19)

4.2.3

Identifiability and informative experiments

Identifiability concerns which parameters of a model that can be estimated uni-quely. When doing analysis of identifiability it is important to separate between identifiability of the model structure and the excitation of the system. The input of the system may not affect some dynamics which makes the associated parame-ters practically unidentifiable.

To study identifiability, the regressor ϕ is used. For each row j of the regressor,

4.2 System identification theory 29

with corresponding regressor rows we can form Φ_NTΦN, which can be used to determine identifiability by calculating its rank.

By removing all unidentifiable parameters and replacing all linearly depen-dent parameters with new parameters we can obtain a fully idepen-dentifiable set, called base parameters. Note that the parameter introduced to replace a linear combi-nation has no physical meaning. Since the goal often is to achieve a predictive model the loss of physical interpretation may not be restrictive.

4.2.4

Method for determining identifiable parameters

Finding the identifiable parameters can be done both analytically and numeri-cally. Finding unidentifiable parameters is a simple procedure. By going through each regressor row ϕj, with respective parameter θj, if all columns of row ϕj is zero the corresponding parameter θj will not affect the model output and is therefore not identifiable.

The other parameter types are not as easy to find directly in the regressor. Finding linear combinations of identifiable parameter and independently identi-fiable parameters is cumbersome analytically, especially for a large models. There-fore a numerical method is preferred.

Numerical method for determining base parameters

A numerical method is presented in Khalil and Dombre [2002]. By using N ran-dom values for the regressor inputs and then placing them in the regressor matrix (4.15), Φ with size (N, nθ), θ has size (nθ, 1). The rank, denoted b, ofΦ is the number of identifiable parameters.Φ can be split into two parts, one which con-tains the independent columns and one with dependent columns.Φp= (Φi,Φd) whereΦi contains the linear independent columns andΦd contains the linear dependent columns. The corresponding is done for splitting θ into θi and θd.

SinceΦispansΦd we can write

Φd =Φiβ (4.20)

where β has dimension (b, nθ−b) and determines the linear combinations. We can now express the regressor model as

YN =Φiθbase+ e (4.21)

where

θbase = θi + βθd (4.22)

is the new parameter vector that will be used. Note that θbaseis ambiguous, since it depends on the order of the columns inΦdandΦi. We now have to determine Φi,Φdand β. By using a qr factorization ofΦ we get

Φ = Q(n,n) " R(c,c)

0_(n−c,c)

#

30 4 System identification

where Q is an orthogonal matrix and R is an upper triangular matrix. By study-ing the diagonal of R, the elements equal to zero indicates linearly dependant columns ofΦ, which are put intp Φd. The remaining independent columns are put intoΦi. What is left to determine is β. By defining

Φp =hΦi Φd i

(4.24) and using a qr factorization of

Φp= h Qi Qd i"R_i R_d 0 0 # (4.25) we have Φd =ΦiR −₁ i Rd (4.26)

and from (4.20) we have

β = R−_i1Rd (4.27) For reference, an alternative numerical method for finding base parameters is found in Gautier [1991]. It uses a svd factorisation instead of a qr factoriza-tion. Additionally, a method for doing this analysis analytically can be found in Gautier and Khalil [1990].

4.2.5

Model structure for identification

In the modelling chapter, two model structures were presented, the inverse dy-namics form and state space form. The inverse dydy-namics form can be written as a regression model and the state-space model results in a greybox model.

Regression model

Since the physical model (3.70) is linear with respect to its parameters, which were presented in Table 3.1, it can be formulated as a regression. The regressions has dynamic variables q, ˙q and ¨q as inputs and torque τ as output. The linearity

property reduces the complexity of the estimation. The regression model has a simple noise and it is assumed that the noise is white. If this requirement is not fulfilled there might be a bias in the estimate.

Grey-box model

The grey-box model is based on the state-space model (3.74) which is nonlinear with respect to its parameters. It has torque τ as input and dynamic variables q,

˙

q and ¨q as outputs. The nonlinearity property makes estimation more difficult

as a more complex and computationally demanding method must be used. One advantage with grey-box modelling is the capability of handling complex noise models. The dynamics of the noise model, which may be non-linear, can be added directly to the state-space model. The resulting model will then have white noise as input.

4.3 Regression model for the gimbal 31

4.2.6

Validation

One commonly used validation method is to use the model for simulation of the system response for a particular input signal. If a model can describe the behaviour of the system well it can be considered being a good model. This is done by using cross validation. It means that the model is estimated on one data sequence and then validated against another sequence. Preferably this data sequence is collected using a different trajectory than the one that is used to esti-mate the parameters, since this reduces the risk for overmodeling.

It is desired to have a measure of how good a model is for simulation, and one such measure is called model fit and is defined as

M_{f it= 100 1 −}ky − ykˆ ky − ¯yk

!

(4.28) where ¯y is the mean of y. Fit is the percentage of the output that is described by

the model.

4.3

Regression model for the gimbal

The regression approach is well proven in the mentioned literature, it is also less complicated than the grey-box model approach and its computational require-ments seems appropriate for the hardware platform. The approach with regres-sion model structure is therefore used.

The physical model is linear with respect to its parameters so it can be written on a regression form

τ = ϕT(q, ˙q, ¨q) θ (4.29) where θ contains all the parameters of the gimbal. τ is the torque from the motors and ϕ is the regression vector. The parameter vector for body i is

θi = h

Iixx Iiyy Iizz `imx `imy `imz Fis Fiv iT

(4.30) Then θifor each body can be combined to form

θ =         θ1 θ2 θ3         (4.31)

The regression matrix has the following form

ϕT(q, ˙q, ¨q) =          ϕ1,1 ϕ1,2 ϕ1,3 · · · ϕ1,3p 0 · · · _{0 ϕ2,p+1} · · · _ϕ2,3p 0 0 · · · _{0 ϕ3,2p+1} · · · _ϕ3,3p          (4.32)

It has three rows, one for each axis of the gimbal.

Since the model of the gimbal is on an inverse form the measurements of the gimbal becomes inputs to the model and the input of the gimbal becomes

32 4 System identification

the output. The regressor model does not handle input noise and therefore the measurements of the gimbal must have low noise levels in order to avoid bias in the estimates. It is therefore important to use accurate sensors and appropriate preprocessing of the measurements.

4.3.1

Base parameters for the gimbal

By using the theory of identifiability we can find the base parameters of the gim-bal. By forming the regressor matrix and doing a quick inspection the unidenti-fiable parameters can be found by looking for zero columns. This yields that the parameters that are unidentifiable are

Su = {I1xx, I1zz, `1mx, `1my, `1mz, `2my} (4.33) These parameters could easily be found by looking at how the motion of the gim-bal is constrained. For instance since the spine of the gimgim-bal has just one axis of rotation, the inertial parameters for the other two axes will not have any effect on the gimbal. Let us use the numerical method to find one set of base parameters. One vector of identifiable parameters is

θi =                                                              I1yy F1s F1v I2xx I2yy `2mx `2mz F2s F2v I3xx I3yy `3mx `3mz F3s F3v                                                             

and the matrix of linear combinations is β =

                                                            1 0 0 0 0 0 0 0 0 −_{1 1 0} 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 −_{1 0} 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0                                                            

The vector of unidentifiable parameters is

θd =         I2zz I3zz `3my        