ESO-Based Robust and High-Precision Tracking Control for Aerial Manipulation

This paper studies the tracking control problem of an aerial manipulator that consists of a quadcopter flying base and a Delta robotic arm. We propose a novel control approach that consists of extended state observers (ESOs) for dynamic coupling estimation, ESO-based flight controllers, and a cooperative trajectory planner. Compared to the state-of-the-art approaches, the proposed one has some attractive features. First, it requires much less measurement information as opposed to the full-body control approaches and hence can be implemented conveniently and efficiently in practice. Second, while the existing approaches estimate the coupling effect based on precise models, the proposed ESOs can do that based on much less information about the system model. The proposed approach is verified by four experiments on a real aerial manipulation platform. The experimental results show that the average tracking error can reach 1 cm by the proposed approach as opposed to 10 cm by the PX4 baseline controller. Although force control is not considered specifically in the approach, the system can complete aerial weaving tasks thanks to the ESOs in the presence of drag forces applied to the end-effector during manipulation.Note to Practitioners—Aerial manipulators have received increasing research attention in recent years due to their wide range of applications. In this paper, we particularly focus on the high-precision and robust control of aerial manipulators. We propose a novel control approach that consists of extended state observers (ESOs) for dynamic coupling estimation, ESO-based flight controllers, and a cooperative trajectory planner. Four experiments on a real aerial manipulation platform demonstrate the effectiveness of the approach. In future research, we will address the control problem when the aerial manipulator contacts the environment.


I. INTRODUCTION
Aerial manipulators have received increasing research attention in recent years due to their wide range of applications.An aerial manipulator combines a flying base, which is usually a multirotor aerial vehicle, with one or multiple robotic arms.Compared to conventional flying robots such as multirotors, an aerial manipulator can interact with the H. Cao, Y. Li, and S. Zhao are With the Research Center for Industries of the Future (RCIF) at Westlake University, School of Engineering at Westlake University, and Westlake Institute for Advanced Study, Hangzhou, China.
environment to complete various tasks such as pick-andplace [1], tree cavity inspection [2], contact-based inspection [3], opening/closing of a valve [4], assembly of structures [5].Compared to ground mobile manipulators, an aerial manipulator can fly to areas that are difficult to reach by humans or ground robots, hence greatly enlarging the scope of manipulation tasks.Aerial manipulation has been studied from various aspects such as platform, control, perception, teleoperation, and cooperation up to now.See [6] and [7] for recent surveys.In this paper, we particularly focus on the high-precision and robust control of aerial manipulators.Different from the control of a multirotor, control of an aerial manipulator is much more challenging due to the multi-DoF robotic arm attached to the multirotor base.The dynamics of the multirotor and the robotic arm are coupled: the movement of the robotic arm and the multirotor base mutually affect each other.
The existing control approaches of aerial manipulators can be divided into two categories: full-body control and separate control.Full-body control requires establishing the nonlinear model of the entire aerial manipulation system and then applies appropriate controller design methods [8].This kind of approach is promising in terms of control accuracy and robustness if the nonlinear dynamics can be precisely modelled.However, they are challenging to implement in practice because on the one hand the precise nonlinear model is difficult to obtain and on the other hand various measurements required by this approach are difficult to obtain.For instance, the full-body control approach proposed in [9] requires velocities and accelerations of the joints of the robotic arm, which is usually difficult to measure for normal actuators.As a comparison, separate control is easier to implement in practice because it does not require knowing the nonlinear model of the entire system and it merely has a minimal requirement on the measurements.The basic idea of separate control is to control the multirotor base and the robotic arm separately.A literature review on related works is given below.
Up to now, researchers have proposed a variety of separate controllers for aerial manipulators.In the early works, the dynamic coupling between the quadcopter and the robotic arm is ignored to simplify controller design [10].As a consequence, the control accuracy is relatively low [11].The work in [12] addressed the dynamic coupling from a mechanical point of view: a sliding battery box was used to compensate for the motion of the robotic arm.However, the performance of this method is limited by the sliding velocity of the battery box.In order to well resolve the problem of dynamic coupling, we need to consider three aspects.
The first aspect is the estimation of the dynamic coupling effects.The work in [13] introduced a variable parameter integral backstepping method to design the quadcopter controller.It can only estimate the coupling torques generated by the displacement of the center of mass from the quadrotor geometry center, which may gain more errors when the velocity of the manipulator increases.The works in [14] and [15] adopted model-based disturbance estimators to estimate the disturbance caused by the manipulators.However, the model-based disturbance estimators rely on accurate models of aerial manipulators that are usually hard to obtain.
The second aspect is the compensation of the dynamic coupling term in the control.To take the external disturbance into account, nonlinear control methods such as sliding mode control [16] and backstepping control [17] are introduced into the motion control of the quadcopter.These approaches treat the dynamic coupling term as unknown external disturbances without estimating it.To suppress the influence of the dynamic coupling term, the works in [13]- [15] integrated mode-based disturbance estimators into controllers to compensate for the dynamics coupling term, requiring the accurate models of the aerial manipulators.
The third aspect is cooperative planning of the trajectories for the multirotor base and the robotic arm.Common approaches for cooperative planning include the closedloop inverse kinematics (CLIK) method [18], [19], the null space-based method [20], and multiple task-priority inverse kinematics method [21].However, these methods did not consider the physical constraints and, therefore, the solutions of these methods may be infeasible.Recently, the work in [22] proposed a nonlinear model predictive control (NMPC) method and formulated the cooperative planning problem as a constrained optimization problem.The NMPC method is, however, computationally expensive, which may request high computational power or a low update rate.
The above analysis reveals the limitations in the three aspects of the existing control approaches for aerial manipulation.This paper aims to overcome these limitations by proposing a new framework that incorporates a coupling term estimator, a motion controller, and a cooperative trajectory planner.The proposed algorithms are verified by carefully designed experiments on an aerial manipulation platform.The novelty of the proposed algorithms is summarized below.1) We propose a partially coupled motion controller.The control of the quadcopter and the delta arm are separate, but based on extended state observers (ESOs) to estimate and compensate for the dynamic coupling between the quadcopter and the manipulator.The controller can be divided into two parts: the flight controller and the manipulator controller.The flight controller is designed by an ESO-based nonlinear control method while the manipulator is controlled by the PID method.The dynamic coupling effect in the mathematical model is divided into two known terms and unknown terms.The known terms are directly used in the flight controller design.The unknown terms are estimated by the ESOs in real time.Compared to [14], [15], the advantage of the proposed ESOs is that they rely on less measurement information for estimating the dynamic coupling.Besides, the closedloop system with the proposed flight controller converges to predesigned dynamics, which are constituted by four desired subloops that can be continently used to tune control gains.All parameters of the controller can be easily determined with the predesigned dynamics according to the performance requirement.
2) A novel cooperative planner with two modes is proposed to coordinate the motions of the quadcopter and the manipulator.The physical constraints of the aerial manipulator are transformed to reduce the computation cost of the planner.To adapt to different manipulation tasks, two control modes, P-P and E-P, are proposed.The P-P mode, designed based on the CLIK method, can be used when the end-effector and quadcopter track separate trajectories.The E-P mode can be used when the end-effector is required to track a given trajectory while the quadcopter base is not.The E-P mode formulates the planning problem as a quadratic programming (QP) problem.The two modes are applicable to different scenarios.The P-P mode is applicable to pick and place task [23], peg-in-hole task, aerial repair task [24], etc.The E-P mode is applicable to the trajectory tracking task, the pulling/pushing task, etc.
The proposed algorithms are verified by four experiments on an aerial manipulator platform.Unlike the traditional Delta arm, the Delta arm used in this paper drives the joint angles by three four-bar linkages to magnify the control forces [25].Experiments (including disturbance rejection, end-effector stabilization, and end-effector trajectory tracking) are conducted to validate the novelties.At last, we apply the proposed method and aerial manipulator to the aerial weaving (see Fig. 1), which shows the ability to complete tasks.
The remainder of this paper is structured as follows.Problem setup and preliminaries used in this paper are given in Section II.The control system overview is presented in

Σ I
inertial frame such that the z-axis is in the direction of the gravity vector Σ B body-fixed frame that is rigidly attached to the quadcopter base at its center of gravity.
Σ D Delta arm frame that is rigidly attached to the arm base at its center p F .
p, p d position and desired position of the quadcopter (p ∈ R 3 and p d ∈ R 3 ) in Σ I q, q d joint vector and desired joint vector of the Delta arm (q ∈ R 3 and q d ∈ R Section III.Section IV proposes the position control system of the quadcopter base.The attitude control of the quadcopter base is proposed in Section V. Section VI gives the proposed cooperative planner.Then, the experimental verification for the proposed methods is given in Section VII. Conclusions are drawn in Section VIII.

II. PROBLEM SETUP AND PRELIMINARIES
This section presents the problem setup and some necessary preliminary results.

A. Problem setup
The aerial manipulator considered in this paper consists of a quadcopter base and a Delta arm (see Fig. 2).The base of the Delta arm is attached underneath the quadcopter.The position of the end-effector of the Delta arm can be controlled by changing the torques applied to the three actuators fixed on the base.However, the orientation of the end-effector remains the same as the base of the Delta arm [26].Hence, the Delta arm has three translational DoFs and three control inputs (i.e., the torques applied to the three actuators).On the other hand, the quadcopter base has six DoFs and four inputs (i.e., the thrusts of the four blades).Therefore, the entire system has nine DoFs and seven control inputs.
The aerial manipulator has three reference frames: the inertial frame Σ I , the quadcopter body-fixed frame Σ B , and the Delta arm frame Σ D (see Fig. 2).Σ I is an inertial frame where the z-axis is in the direction of the gravity vector.Σ B is rigidly attached to the quadcopter base.Its origin coincides where R B D ∈ SO(3) is the rotation matrix from Σ D to Σ B .The lengths for the upper and lower arms are represented by l U ∈ R and l L ∈ R as illustrated in Fig. 2. Circumradii of the top base and the bottom end-effector base are, respectively, defined as r F ∈ R and r M ∈ R. The relationship between the end-effector position p D E ∈ R 3 and the joint vector q = [q 1 , q 2 , q 3 ] T ∈ R 3 is where On the one hand, given a joint vector q, the position p D E can be solved from (2) based on the forward kinematics.On the other hand, given a position p D E , the joint vector q can be solved from (2) by the inverse kinematics.Details can be found in [27], [28].
As can be seen from Fig. 2, the joint angles of the Delta arm are driven by planar four-bar linkages.The relationship between the joint angles and the crank position angles can be calculated by the kinematics of the planar four-bar linkage [29,Section 3.6].

C. Kinematics of the aerial manipulator
The position of the end-effector in Σ I can be calculated by p E = p + Rp B E , where R ∈ SO(3) is the rotation matrix from Σ B to Σ I .The time derivative of p E is where ω ∈ R 3 is the angular velocity vector of the quadcopter expressed in Σ B , [•] × denotes the skew-symmetric matrix.Let Then, (4) can be rewritten as where J = [I 3 , RR B D ] is the Jacobian matrix, I 3 is the 3 × 3 identity matrix.

D. Dynamics of the aerial manipulator
The dynamics of the aerial manipulator can be divided into the dynamics of the quadcopter base and the Delta arm.Inspired by [14], the dynamics of the quadcopter base are surveyed to design the flight controller of the aerial manipulator.The dynamics of the quadcopter base in the aerial manipulator are where v ∈ R 3 represents the velocity vector of the quadcopter base, g ∈ R 3 is the gravity vector, m B represents the mass of the quadcopter base, m M is the mass of the Delta arm, M ∈ R 3×3 denotes the inertia matrix of the quadcopter base, τ ∈ R 3 represents torque vector of the rotors.In addition, f = R[0, 0, f ] T ∈ R 3 is control force vector, where f ∈ R represents total force of the rotors.
Here, f c and τ c are dynamic coupling force and torque caused by the Delta arm.They are given as [14] where p B C is the center of mass of the aerial manipulator in Σ B , M B M is the inertia matrix of the Delta arm in Σ B .
Although the expressions of f c and τ c are given in ( 7) and ( 8), they cannot be measured in practice since Ṁ B M , ṗB C , and pB C are difficult to measure.The mass m M and the inertia matrix M B M of the Delta arm are contained in model ( 6).In the absence of the Delta arm, we have m M ≡ 0, M B M ≡ 0, p B C ≡ 0 and Ṁ B M ≡ 0. In this case, model ( 6) degenerates to a standard quadcopter model.

E. Extended state observer
Preliminaries of the ESO are given below.The stability analysis and other details about the ESO can be found in [30].
Consider a system with input u ∈ R and output y ∈ R in the form y (n) = ∆ + u, where ∆ ∈ R is the unknown term of the system.The output y can be measured by sensors.Denote n) .We define y n+1 as an extended state variable and y n+1 = ∆.The original plant is now described as ẏ1 = y 2 , ẏ2 = y 3 , • • • , ẏn−1 = y n , ẏn = y n+1 + u, ẏn+1 = ∆.The purpose of the ESO is to estimate ∆ by utilizing the measurable output y.
Let ŷ1 , • • • , ŷn+1 denote the estimated values of y 1 , • • • , y n+1 .Then, the ESO can be designed as follows [30]: where ) denote the estimation errors.Combining the original plant and (9), the error dynamics are given as According to [30], if ∆ is bounded, i.e., | ∆| ≤ δ, then there exists a constant c i > 0 and a finite time for some positive integer k, where ỹsum (0 11), it can be concluded that the estimation error bound of the ESO is determined by the initial estimation error, δ, and observer bandwidth w o .Since the unknown term ∆ cannot be measured, the initial estimation error and δ cannot be utilized to reduce c i .Gain w o is inversely proportional to 1/c k i .Then, a smaller value of c i can be achieved by an ESO with a larger w o .
Quadcopter Delta arm PID controller Fig. 3: The proposed control structure of the aerial manipulator system.easily determined with the predesigned dynamics according to the performance requirement.
The third component the Delta arm controller.Its input is the desired joint angles of the Delta arm.The output is the torques that each actuator should generates.The Delta arm controller is constituted by three actuator controllers.We use Dynamixel AX28 servomotors as the actuators of the Delta arm.The controllers of the servomotors are designed by the traditional PID method [28].

IV. POSITION CONTROL OF THE QUADCOPTER BASE
In this section, the position controller of the quadcopter base is designed by an ESO-based nonlinear control method.First, the coupling between the quadcopter base and the Delta arm is estimated by a position ESO.Then, the feedback linearization control method is adopted to design the position controller.
From model ( 6) and ( 7), the position dynamics can be rewritten in a compact form as where the coupling force is rewritten as ) represents the first order item of f c and can be used in controller design directly, represents the second order item of f c .∆ v is treated as a unknown item, because ω, ṗB C and pB C can not be sensed by the aerial manipulator.

A. Position ESO for coupling estimation
A third-order position ESO is designed to estimate unknown item ∆ v .Let p, v, ∆v denote the estimates of p, v, ∆, respectively.The position ESO is designed as where pi , vi , ∆v,i denote the i-th elements of p, v, ∆, respectively.In addition, w p = [w p,1 , w p,2 , w p,3 ] T ∈ R 3 is a positive constant vector and u v,i is the i-th element of u v .
According to Section II-E and (10), the estimate error is bound, i.e., lim t→∞ |∆ v,i − ∆v,i | ≤ c p,i , i = 1, 2, 3.According to (10), c p,1 , c p,2 , c p,3 are inversely proportional to w p,1 , w p,2 , w p,3 to the k-th power.Then, the position ESO with the larger gain vector w p achieves the smaller estimate error bound.Let c p denote the maximum of bounds c p,1 , c p,2 , c p,3 , i,e., c p = max{c p,1 , c p,2 , c p,3 }.Then, we have

B. Nonlinear position controller
The dynamic behavior is significant for the aerial manipulator system.To adapt more flight mode of the quadcopter base (e.g.position mode, altitude model), the closed-loop dynamics of the position loop are divided into two subloops: position and velocity subloops.Let p d ∈ R 3 and p = p−p d ∈ R 3 denote the desired position and the position error of the quadcopter base, respectively.The desired dynamics of the two subloops are given as Fig. 3: The proposed control structure of the aerial manipulator system.

III. CONTROL SYSTEM OVERVIEW
The overall control system of the aerial manipulator is decomposed into three components as illustrated in Fig. 3.
The first component is the cooperative planner.Its input is the desired trajectory of the end-effector p E,d .Its outputs are the planned trajectories p d and q d for the quadcopter base and the Delta arm.The proposed planner has two modes: P-P and E-P modes, which makes the planner suitable for different scenarios.Besides, the physical constraints are considered in the planner to ensure the solution is feasible.Compared with the existing methods, the proposed planner adapts to more scenarios since it has two modes and considers the physical constraints.
The second component is the flight controller for the quadcopter base.Its input is the desired trajectory planned by the cooperative planner.Its outputs are τ and f , which are the torque and force commands of the quadcopter.This flight controller can be further decomposed into a few subcomponents.The first is an ESO-based position controller.It aims to generate the force vector f for the quadcopter so that the desired position p d can be tracked.The closed-loop position system converges to predesigned dynamics that are constituted by position and velocity subloops.The second is an ESO-based attitude controller.It aims to generate the torque vector command for the quadcopter so that the desired attitude R d can be tracked.The closed-loop attitude system also converges to predesigned dynamics that are constituted by attitude and angular velocity subloops.All the corresponding sections introducing these subcomponents are listed in Fig. 3.All parameters of the flight controller are included in the predesigned dynamics and can be determined according to the performance requirement.
The third component is the Delta arm controller.Its input is the desired joint angles of the Delta arm.The output is the torques that each actuator should generate.The Delta arm controller is constituted of three actuator controllers.We use Dynamixel AX28 servomotors as the actuators of the Delta arm.The controllers of the servomotors are designed by the traditional PID method [31].

IV. POSITION CONTROL OF THE QUADCOPTER BASE
In this section, the position controller of the quadcopter base is designed by an ESO-based nonlinear control method.First, the coupling between the quadcopter base and the Delta arm is estimated by a position ESO.Then, the feedback linearization control method is adopted to design the position controller.
From model ( 6) and ( 7), the position dynamics can be rewritten in a compact form as where the coupling force is rewritten as ) represents the first order term of f c and can be used in controller design directly, represents the second order term of f c .∆ v is treated as an unknown term, because ω, ṗB C and pB C can not be sensed by the aerial manipulator.

A. Position ESO for coupling estimation
A third-order position ESO is designed to estimate unknown term ∆ v .Let p, v, ∆v denote the estimates of p, v, ∆, respectively.According to Section II-E, the position ESO is designed as where pi , vi , ∆v,i denote the i-th elements of p, v, ∆, respectively, u v,i is the i-th element of u v .In addition, w p = [w p,1 , w p,2 , w p,3 ] T ∈ R 3 is the bandwidth vector of the position ESO and it can be adjusted.
According to Section II-E and (11), the estimation error is bounded, i.e., lim According to (11), c p,1 , c p,2 , c p,3 are inversely proportional to w p,1 , w p,2 , w p,3 to the k-th power.Then, the position ESO with the larger gain vector w p achieves the smaller estimation error bound.Let c p denote three times the maximum of the bounds c p,1 , c p,2 , c p,3 , i,e., c p = 3 max{c p,1 , c p,2 , c p,3 }.Then, we have lim

B. Nonlinear position controller
The dynamic behavior is significant for the aerial manipulator system.To control more flight modes of the quadcopter (e.g., position mode, altitude mode), the closed-loop dynamics of the position loop are divided into two subloops: position and velocity subloops.The two-subloop structure matches the cascade structure of the autopilot used in the quadcopter.Let p d ∈ R 3 and p = p − p d ∈ R 3 denote the desired position and the position error of the quadcopter base, respectively.The desired dynamics of the two subloops are given as where Λ p ∈ R 3×3 and K v ∈ R 3×3 are positive diagonal matrices.In addition, r v = v − v r , where v r is the desired velocity reference trajectory.To achieve the desired position subloop, we design v r as The relationship between v r and ṗd can be shown in (18).Then, we have From (19), one can conclude that the desired position subloop can be achieved as r v → 0. The relationship between r v and p is described in Lemma 1.
, where tr(Λ −1 p ) denotes the trace of Λ −1 p .Proof: According to the definition of r v , we have Let pi denote the i-th element of p.Let [Λ p ] i,j denote the element in the i-th row and j-th column of Λ p .Then, we have where , where A i is the Hurwitz matrix and its eigenvalue is [Λ p ] i,i with multiplicity 2. According to bounded input bounded output (BIBO) stability [32, Section 5.2], we have To achieve the desired velocity subloop, the position control law is designed as Substituting control law ( 22) into (12) gives According to Section IV-A, ∆ v − ∆v is bounded.It means that the desired velocity subloop can be achieved boundedly.
Then, the control force vector f can be inferred from ( 13): We now prove the stability of the position dynamics with control law (24).Theorem 1: Assuming that ∆v is bounded, i.e., ∆v ≤ δ v , if the position controller is designed as (24) with position ESO (15), then the position tracking error p is bounded, i.e., p ≤ c p tr(Λ −1 p )/λ min (K v ) as t → ∞, where c p is the upper bound of the estimation error from the position ESO.
Proof: Combining the position control law (24) and position model (12), the error dynamics can be written as Define a Lyapunov function as With error dynamics (25), the time derivative of V is According to Section IV-A, the estimation error of the position ESO (15) After finite time T p , we can obtain where λ min (K v ) is the minimum eigenvalue of K v .Since K v is a positive diagonal matrix, we have λ min (K v ) > 0.
According to the BIBO stability [32, Section 5.2], we have With the definition of W , we have r v ≤ c p /λ min (K v ) as t → ∞.The relationship between p and r v are given in Lemma 1.
From the relationship, we have p ≤ c p tr(Λ −1 p )/λ min (K v ) as t → ∞.
Theorem 1 assumes that ∆v is bounded.The unknown term ∆ v is formulated as (14) and, therefore, ∆v is determined by motions of the quadcopter base and the Delta arm.This assumption is reasonable since the motions are physically limited.The control parameters in the position controller are Λ p and K v .The two parameters are contained in the desired dynamics of the two subloops.The desired dynamics ( 16) is expressed as a second-order system while the desired dynamics ( 17) is expressed as a first-order system.Therefore, the control parameters can be determined by utilizing the features of these systems and the performance requirements.
An analysis of the tracking error is given below.The upper bound of the position tracking error is proportional to c p , tr(Λ −1 p ) and inversely proportional to λ min (K v ).As a result, the controller with the larger control gains Λ p and K v achieves the lower position tracking error bound.Another way to reduce the position tracking error bound is by reducing the estimate bound c p .According to Section IV-A, the position ESO with the larger bandwidth vector w p results in a lower value of c p .Therefore, the controller with the larger gain vector w p can achieve the lower position tracking error bound.

V. ATTITUDE CONTROLLER OF THE QUADCOPTER BASE
In this section, an attitude controller of the quadcopter base is proposed.In particular, the coupling between the quadcopter base and the Delta arm is estimated by an attitude ESO.Then, the feedback linearization control method is adopted to design the attitude controller of the quadcopter base.
From model ( 6) and ( 8), the dynamics of the attitude loop can be rewritten as The first term in the second equation of ( 29) is a composite vector expressed as the coupling torque can be rewritten as τ c = τ s + M ∆ ω , where represents the first order term of τ c and can be used in controller design directly.The second term in the second equation of (29) represents the second order term of τ c and can be expressed as By definition, the dynamic influence ∆ ω is a function of ω, v, Ṁ B M , ṗB C , pB C .We treat ∆ ω as an unknown term to be estimated, since ω, v, Ṁ B M , ṗB C , pB C are not measurable.

A. Attitude ESO for coupling estimation
An attitude ESO is designed to estimate unknown term ∆ ω according to the method in Section II-E.Let ω, ∆ω denote the estimates of ω, ∆ ω , respectively.To estimate ∆ ω , the attitude ESO is designed as where u ω,i is the i-th element of u ω , w ω = [w ω,1 , w ω,2 , w ω,3 ] T ∈ R 3 is the bandwidth vector of the attitude ESO and it can be adjusted, ωi , ∆ω,i denote the i-th elements of ω and ∆ω , respectively.

B. Nonlinear attitude controller
Let 3) denote the desired rotation matrix of the quadcopter.We define a = [cos ψ, sin ψ, 0] T .Then, we have Let ω d denote the desired angular velocity of the quadcopter and it can be calculated by where [•] ∨ denotes the vee map which is the inverse operation of [•] × .The error rotation matrix is defined as According to [33], the error function used in SO( 3) is equivalent to the vector part of error quaternion Let Ṙ be the derivative of the error rotation matrix.It can be obtained by To control the flight modes of the quadcopter (e.g., stabilized mode, acro mode), the closed-loop dynamics of the attitude loop are divided into two subloops: attitude and angular velocity subloops.According to [34], the designed dynamics of the two subloops are given as where K ω ∈ R 3×3 and Λ ∈ R 3×3 in (42) are positive diagonal matrices, k β > 0 is a constant gain.Moreover, r ω = ω − ω r , where ω r is the desired angular velocity reference trajectory.
To achieve the desired attitude subloop, we design ω r as Then, we have According to [34], the desired attitude subloop can be achieved as r ω → 0.
To achieve the desired angular velocity subloop, the attitude controller of the quadcopter base is designed as follows: where where I 3 ∈ R 3×3 is the identity matrix of order 3 × 3.
Substituting the control law (42) into (29) yields According to Section V-A, ∆ ω − ∆ω is bounded.It means that the desired angular velocity subloop can be achieved boundedly under the control law (42).Then, the control torque vector τ is given as We now prove the stability of the attitude loop with control law (45).Theorem 2: Assuming that ∆ω is bounded, i.e., ∆ω ≤ δ ω , if the attitude controller is designed as (45) with the attitude ESO (33), then the attitude tracking error is bounded, i.e., βv (t) ≤ c q as t → ∞, where c q is the upper bound of the attitude tracking error.
From Theorem 2, one can conclude that the attitude tracking error is bounded, i.e., βv (t) ≤ c q as t → ∞.An analysis of the tracking error is given below.By definitions of k 1 , k 2 , k 3 , and k 4 , (54) can be rewritten as where Then, we can conclude that the upper bound of the attitude tracking error is determined by K ω , k β , λ min (Λ q ), and c ω .The control gain K ω , k β , and λ min (Λ q ) are inversely proportional to upper bound c q .As a result, the controller with the larger control gains achieves the lower attitude tracking error bound.Another way to reduce the estimation error of the attitude ESO is by reducing the estimate bound c q , since the estimation error bound is proportional to c q .According to Section V-A, the attitude ESO with the larger bandwidth vector w ω obtains a lower value of c ω .Therefore, the controller with the larger gain vector w ω can achieve the lower attitude tracking error bound.

VI. COOPERATIVE PLANNING
In this section, we design a cooperative planner for the aerial manipulator.The relationship between the motion controller and the cooperative planner is described in Section III.The references of the motion controller are calculated by the cooperative planner with the given end-effector trajectory.The motion controller tracks the result references.
The physical constraints affect the feasibility of the planner.The obtained solution may violate such constraints, resulting infeasible results for the robot without considering the physical constraints.Therefore, the physical constraints of the aerial manipulator are considered in the cooperative planning.With the constraints, a planner with two modes is designed to coordinate the motions of the quadcopter and the multirotor.The two modes adapt to different scenarios.

A. Physical constraints
The physical constraints of the aerial manipulator consist of position, velocity, and acceleration constraints and can be written as These constraints are expressed in terms of s, ṡ, s, respectively.Using these expressions will introduce 18 variables in solving the cooperative problem.To reduce the number of variables, we convert the constraints of s, s to those of ṡ.Then, the number of the variables in solving the cooperative problem is reduced to 6.
Lemma 2: The physical constraints (57)-( 59) can be approximated as where ṡ = max{s p , ṡmin , ṡv } and ṡ = min{s p , ṡmax , ṡv }.In addition, s p and s p are lower and upper approximate bounds of the actuated state position and can be calculated by (62).ṡmin and ṡmax can be determined by (58).ṡv and ṡv are lower and upper approximate bounds of the actuated state acceleration and can be calculated by (66).Proof: First, the position constraint is now converted to the constraint of ṡ.Let ∆t denote the sampling time.The actuated state in the next time step is defined as s next ∈ R 6 .According to the position constraint, s next should satisfy the position constraint, i.e., Then, the position constraint (61) can be rewritten as (62) (64) (64) (58) (58) si ṡi Fig. 3: The approximate physical condition.Blue line: approximate position constraint (41).Purple line: velocity constraint (37).Green line: approximate condition (43).Gray area: feasible area (39).line passing through the two points to linearise the viability condition.Then the viability condition (42) is linearised as From (43), we have where ṡv ∈ R 6 , ṡv ∈ R 6 are the lower and upper bound of the viability condition, ṡv,i , ṡv,i respectively represent the i-th element of ṡv and ṡv , The relationships among the approximate position constraint, velocity constraint and approximation condition are illustrated in Figure 3.

B. Cooperative Method
The proposed planner provides two modes: End-effector position-Quadcopter position (P-P) and End-effector position (E-P) modes.The P-P mode steers the manipulator and the quadcopter base with the references of the end-effector and the quadcopter, respectively.The E-P mode steers the manipulator and the quadcopter base only with the reference of the end-effector.The P-P mode uses quadcopter base for large scale movements and manipulator for compensating the tracking error of the end-effector.It performs better in tracking errors due to the less relative motion.However, the manipulation area of the P-P mode is small since the size of the manipulator's workspace is limited.Compared with the P-P mode, t area while s The two mod mode is app task, aerial r to trajectory 1) P-P m positions of (4), we obta Inspired by t end-effector designed as where p E,d R 3 is the de a positive di We then manipulator vectors that Then, we ha Thus, the des by integratin can be obta manipulator 2) E-P m copter base effector posi formulated a whre F ( ṡ) i positive diag s c is a feedb  Second, we introduce the viability condition to handle the acceleration constraint [35].Let s i denote the i-th element of s.Then, the viability condition is expressed as where smin,i , smax,i , s min,i , s max,i denote the i-th elements of smin , smax , s min , s max , respectively.However, the viability condition (63) is nonlinear.We linearise the viability condition by lines passing through the two terminal points of the curves in (63).For instance, the curve of the first equation in (63) passes through (s i = s i,max , ṡi = 0) and (s i = s max,i + ṡ2 max,i /(2s m ), ṡi = ṡmax,i ).We can use a line passing through the two points to linearise the viability condition.Then the viability condition (63) is linearised as From (64), we have where ṡv ∈ R 6 , ṡv ∈ R 6 are the lower and upper bound of the viability condition, ṡv,i , ṡv,i respectively represent the i-th element of ṡv and ṡv , The relationships among the approximate position constraint, velocity constraint, and approximation condition are illustrated in Fig. 4.

B. Cooperative method
The proposed planner provides two modes: End-effector position-Quadcopter position (P-P) and End-effector position (E-P) modes.The P-P mode steers the Delta arm and the quadcopter base with the references of the end-effector and the quadcopter, respectively.The E-P mode steers the Delta arm and the quadcopter base only with the reference of the end-effector.The P-P mode uses the quadcopter base for large-scale movements and the Delta arm for compensating for the tracking error of the end-effector.It performs better in tracking errors due to the less relative motion.However, the manipulation area of the P-P mode is small since the size of the Delta arm's workspace is limited.Compared with the P-P mode, the E-P mode has a large size of manipulation area while suffers more influence from dynamic couplings.The two modes are applicable to different scenarios.The P-P mode is applicable to pick-and-place task [23], peg-in-hole task, aerial repair task [24], etc.The E-P mode is applicable to trajectory tracking task, pulling/pushing task, etc.
1) P-P mode: The inputs of the P-P mode are desired positions of the end-effector and the quadcopter base.From (4), we obtain Inspired by the CLIK method [20, Section VA], the desired end-effector velocity in the Delta arm frame Σ D is designed as where p E,d ∈ R 3 is the desired end-effector position, ṗE,d ∈ R 3 is the desired velocity of the end-effector, K c ∈ R 3×3 is a positive diagonal matrix.
We then introduce the physical constraints of the aerial manipulator into the P-P mode.Let ṗD E and ṗD E denote the vectors that consist of the last three elements of ṡ and ṡ.Then, we have Thus, the desired end-effector position in Σ D can be obtained by integrating ṗD E,d .The desired actuated joint vector q d can be obtained by the inverse kinematics of the Delta arm with p D E,d from (69).
2) E-P mode: To calculate the references of the quadcopter base and the Delta arm with the given desired endeffector position, the cooperative problem is mathematically formulated as a QP problem: whre F ( ṡ) is the cost function of the QP problem, W is a positive diagonal matrix, J is the actuated Jacobian matrix, s c is a feedback inspired by the CLIK algorithm in [20, Section VA] and designed as The cost function shows the constraints on the energy of the cooperative planning.The weight matrix W also reflects the priority of the actuated DOFs.A smaller diagonal element means a higher priority.Therefore, we define W as follows: where w u > 0 is the weighting coefficient to determine the priorities of the quadcopter base and the Delta arm.w u with a small value means the Delta arm has a higher priority.w u with a large value means the Delta arm has a lower priority.

VII. EXPERIMENTAL VERIFICATION
This section presents experimental results to verify the effectiveness of the proposed algorithms.The experimental video is available at https://youtu.be/QMjGtBCUi-E.
First of all, we describe the experimental setup (see Fig.   II.

A. Example 1: Disturbance rejection
In the first experiment, a payload of 1.2 kg is mounted at the end of the Delta arm (Fig. 6(b)).The control objective is to maintain the position of the quadcopter base unchanged while the payload swings fast (Fig. 6(a) and (c)).Fig. 6(c) shows the tracking performance of the Delta arm.In particular, the reference of x D E is a sinusoidal signal with the period as 0.5π s.Since it is mechanical control, the tracking error is maintained on the level of 10 −3 m.
The purpose of this experiment is to verify the effectiveness of the ESO and the flight controller.A PX4 1.9.0 controller is selected as a baseline for comparison.The PX4 controller has been widely used in the motion control of the aerial manipulator [36]- [38].
Fig. 6(d)-(e) shows the position tracking performances of the proposed controller and the PX4 controller.As can be seen in Fig. 6(e), with the proposed controller, the mean position error of the quadcopter base is 0.02 m and the maximum position error reaches 0.04 m.As a comparison, with the PX4 controller, the mean position error of the quadcopter base is 0.09 m and the maximum position error reaches 0.20 m.Compared to the PX4 controller, our proposed controller can reduce the mean and maximum position errors by 78% and 80%, respectively.
In addition to the experiment on the real platform, we also conducted a simulation to verify the effectiveness of the proposed ESOs.The reason that we conduct simulation experiments is that the ground truth of the disturbances generated by the swing payload is unknown in real experiments but can be known in simulation.The simulation is conducted in Matlab Simscape without considering the noise of the aerial manipulator, a mechanical system simulation environment.The reference of the end-effector in the simulation is the same as that in the real experiment (Fig. 7(b)).
As shown in Fig. 7(a), the mean position error of the quadcopter base achieved by the proposed controller is 0.002 m, which is 1/10 of the mean error in the real experiment.It is, however, not surprising because the real experiment is more challenging than the simulation due to various uncertainties.More importantly, Fig. 7(c) shows the estimation results of the position and the attitude ESOs.The true values of the dynamic coupling terms are calculated from ( 12) and ( 29).As can be seen, the proposed ESOs can well track the dynamic coupling terms though there is a phase delay of about 0.18 s.The good performance of the ESOs underpins the overall performance of the ESO-based controller.
To further verify the robustness of the proposed method, we conducted Monte-Carlo simulation.The simulation considers measurement noises and inertia uncertainties.The details of the noises and the uncertainties are given in Table III.In particular, we consider three cases for the ranges of the uncertainties: ±5%, ±10%, and ±20% of the corresponding nominal values.The three cases are referred to as low, middle, and high uncertainties, respectively.One thousand times of stochastic simulation were conducted for each level  The mean values and the standard derivations of the position error of the quadcopter base in the Monte Carlo simulation results are shown in Table IV.With the proposed method, the mean values of the results in low, middle, and high uncertainties are 0.009 m, 0.016 m, and 0.024 m, respectively.The standard derivations of the results in the low, middle, and high uncertainties are 0.004 m, 0.008 m, and 0.018 m, respectively.The results demonstrate the robustness of the proposed method given different levels of uncertainties.

B. Example 2: End-effector stabilization
The control objective of this experiment is to stabilize the end-effector to a fixed position while the quadcopter is flying along a circle with a diameter of 0.12 m (see Fig. 8(a)).The size of the circle is limited by the size of the Delta arm's workspace to ensure the end-effector can reach the fixed position.The purpose of this experiment is to validate the P-P mode of the proposed cooperative planner together with the flight controller.This scenario considered in this experiment is important because practical tasks may require assigning separate trajectories for the quadcopter base and the end-effector.For example, the end-effector must track a desired trajectory to complete a manipulating task whereas the quadcopter base must move to avoid obstacles or observe the manipulating target from different angles using onboard cameras.
The control results are shown in

C. Example 3: End-effector trajectory tracking
The control objective in the third experiment is that the end-effector can accurately track a desired complex trajectory.The purpose of this experiment is to validate the E-P mode of the proposed cooperative planner together with the flight controller.Different from the second example, this example does not require specifying the desired trajectory of the quadcopter base.
The desired trajectory of the end-effector is set as a shape of a four-petal flower (Fig. 9(a)).In particular, let

D. Example 4: Aerial weaving
In the final example, a rope is tied at the end-effector.The control objective is that the end-effector must pass through four hangers precisely and complete an aerial weaving task (see Fig. 1).The four hangers are mounted on the top of four pillars.A passive rope dispenser is mounted on one of the pillars.The rope must be pulled out from the dispenser by the aerial manipulator during the flight.Therefore, this example is more challenging than the third one because the aerial manipulator must overcome the drag force of the rope.Hence, this example can better demonstrate the robustness of the overall system in practical manipulation tasks.
The reference of the end-effector is calculated based on the positions of the hangers.The control is archived by the flight controller together with the E-P mode of the cooperative planner.As shown in Fig. 10, the task is completed by the proposed method.In particular, the mean tracking error of the end-effector is 0.02 m, which is greater compared to the third example where the mean tracking error is 0.01 m.This is reasonable due to the drag force of the rope.Since force control is not specifically considered in our algorithms, the present performance verifies the ability of the proposed algorithm to counter external forces applied to the endeffector.

VIII. CONCLUSION
This paper proposed an ESO-based approach for endeffector tracking control of an aerial manipulator.It is verified by four experimental results.The results of the disturbance rejection experiment show the mean position error of the quadcopter base with the proposed method is 0.02 m when the Delta arm swings fast with a payload of 1.2 kg.As a comparison, the mean position error of the quadcopter base with the PX4 controller is 0.09 m.The mean tracking errors of the end-effector stabilization and the end-effector trajectory tracking experiments are 0.011 m and 0.01 m, respectively.In the aerial weaving experiment, the mean tracking error of the end-effector is 0.02 m.These results verify that the proposed approach can achieve accurate and robust control performance.The approach requires minimal measurement information and hence can be implemented conveniently in real time.Although force control is not specifically addressed in this approach, the ESOs can estimate external disturbances applied to the end-effector to a certain extent.It is verified by the aerial weaving experiment that the aerial manipulator with the proposed method can accomplish the experiment under the drag force of the rope.The entire system shows certain robustness against external forces.However, in order to handle strong force interaction between the manipulator and the environment, force sensors and force control must be introduced.This will be one important research direction for future research.
This work was supported by the Research Center for Industries of the Future (RCIF) at Westlake University (Grant WU2022C027) and the Hangzhou Key Technology Research and Development Program (Grant 20200416A16).

3 )pF
E , p E,d position and desired position of the end-effector (p E ∈ R 3 and p E,d ∈ R 3 ) in Σ I p D E , p D E,d position and desired position of the end-effector (p D E ∈ R 3 and p D E,d ∈ R 3 ) in Σ D p B E position of the end-effector (p B E ∈ R 3 ) in Σ B R rotation matrix (R ∈ SO(3)) from Σ B to Σ I R B D rotation matrix (R B D ∈ SO(3)) from Σ D to Σ B ψ, ψ d yaw angle and desired yaw angle of the quadcopter (ψ ∈ R and ψ d ∈ R) f, τ total force (f ∈ R) and torque vector (τ ∈ R 3 ) of the rotors fc, τc dynamic coupling force vector (fc ∈ R 3 ) and torque vector (τc ∈ R 3 ) caused by the Delta arm τ M torque vector of the Delta arm's actuated joints (τ M ∈ R 3 ) p B center position of the Delta arm's base (p B F ∈ R 3 ) in Σ B
where s = [p T , p D,T E ] T consists of the position of the quadcopter base in Σ I and the position of the end-effector in Σ D , [•] min denotes the lower bound of the corresponding variable, and [•] max denotes the upper bound.
5).The aerial manipulator platform used in the experiments consists of a quadcopter and a Delta arm.The wheelbase of the quadcopter is 0.65 m.The mass of the quadcopter (including a battery) is 3.60 kg.The Delta arm consists of a mounting base (0.56 kg) and a movable robotic arm (0.44 kg).The ESO-based flight controller runs on a Pixhawk 4 autopilot.The cooperative planner and the inverse kinematics of the Delta arm run on an onboard Intel NUC i7 computer with ROS (an open-source robotics middleware suite).The experiments are conducted in a Vicon system, which provides

Fig 8 .
First, it is shown in Fig 8(b) that the mean position error of the end-effector is 0.011 m.The standard deviation of the position error is 0.004 m.Second, the planned trajectories for the quadcopter base and Delta arm are shown in Fig 8(c).As can be seen, the mean tracking error of the quadcopter base can reach 0.018 m.

p
E,0 = [x E,0 , y E,0 , z E,0 ] T be the initial position.The desired trajectory of the end-effector is designed asx E,d = x E,0 + 0.3(1 − cos 0.08πt) cos 0.02πt, y E,d = y E,0 + 0.3(1 − cos 0.08πt) sin 0.02πt, z E,d = z E,0 .(73)The tracking performance of the end-effector is shown in Fig.9(b).The mean tracking error of the end-effector is as low as 0.01 m.

TABLE I :
List of important notations.
is bounded.The upper bound of the estimation error is represented as c p .Then, we have ∆ v − ∆v ≤ c p as t > T p , where T p represents a finite time.Using the comparison theorem [32, Section 9.3], we define W

TABLE II :
Physical constraints of the aerial manipulator.

TABLE III :
Measurement noises and inertia uncertainties of the Monte Carlo simulation.

TABLE IV :
Mean values and standard derivations of the Monte Carlo simulation results.