<figure><img src="image_1.png" /></figure>

HAL
open science

Controlled Shaking of Trees With an Aerial Manipulator

Antonio González-Morgado, Eugenio Cuniato, Marco Tognon, Guillermo

Heredia, Roland Siegwart, Aníbal Ollero

► To cite this version:

Antonio González-Morgado, Eugenio Cuniato, Marco Tognon, Guillermo Heredia, Roland Siegwart, et al. Controlled Shaking of Trees With an Aerial Manipulator. IEEE/ASME Transactions on Mechatronics, 2024, pp.1-12. 10.13039/501100011033. hal-04618213v2

HAL Id: hal-04618213

https://hal.science/hal-04618213v2

Submitted on 21 Jun 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.
# Controlled Shaking of Trees With an Aerial Manipulator

Antonio González-Morgado, Eugenio Cuniato, Marco Tognon, Guillermo

Heredia, Roland Siegwart, Aníbal Ollero

► To cite this version:

Antonio González-Morgado, Eugenio Cuniato, Marco Tognon, Guillermo Heredia, Roland Siegwart, et al. Controlled Shaking of Trees With an Aerial Manipulator. IEEE/ASME Transactions on Mechatronics, 2025, 30 (2), pp.1026-1037. 10.13039/501100011033. hal-04618213v2

HAL Id: hal-04618213

https://hal.science/hal-04618213v2

Submitted on 21 Jun 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>
# Controlled Shaking Of Trees With An Aerial Manipulator

A. Gonzalez-Morgado, *Graduate Student Member, IEEE*, E. Cuniato, *Graduate Student Member, IEEE*, M. Tognon, *Member, IEEE*, G. Heredia, *Member, IEEE*, R. Siegwart and A. Ollero, *Fellow, IEEE*

***Abstract***—In recent years, the fields of application of aerial manipulators have expanded, ranging from infrastructure inspection to physical interaction with flexible elements such as branches and trees. This paper presents the controlled shaking of a tree with an aerial manipulator. Our work aims at contributing to applications like the identification of tree parameters for environmental health monitoring or the collection of samples and fruits by vibration. To this end, we propose a control strategy for controlled shaking of flexible systems. We adopt a self-excited oscillation strategy that induces vibrations at the natural frequency of the system, at which the greatest amplification and therefore the greatest vibrations occur. Likewise, this work presents a simplified 1 DoF model based on the Rayleigh-Ritz method to analyze dynamic interaction between a tree and the aerial manipulator with the controlled shaking strategy. The proposed control strategy is evaluated through indoor experiments, where an aerial manipulator shakes an indoor tree made of bamboo canes. Experimental results show how the proposed model can estimate properly the amplitude of the vibration and the frequency of the vibration, depending on the grasping point and the control gain of the self-excited oscillation strategy<fnref n="1" />.

***Index Terms***—Aerial Robotics, Aerial Physical Interaction, Control for Manipulation, Control of Flexible Elements, Controlled Shaking

<figure><img src="image_2.png" /><figcaption>Fig. 1. Aerial manipulator controlled shakes a tree. Different frames during the oscillation of the system are shown.</figcaption></figure>

## I. INTRODUCTION

AERIAL manipulators bring together the flexibility and maneuverability of aerial platforms with the functionalities and capabilities of robotic manipulators [1]. They offer a chance to diminish the duration, expenses, and potential hazards associated with operations conducted at elevated, dangerous, or hard-to-reach areas. This led to the expansion of application fields for UAVs (Unmanned Aerial Vehicles) which range from support in post-disaster situations [2] to the inspection of oil and gas refineries [3], chemical plants [4], power lines [5], or civil infrastructures [6]. This shows the huge recent development of aerial manipulators.

However, some challenges need to be addressed to develop secure, robust and outdoor capable aerial platforms. Standard coplanar multirotors are underactuated systems, since they can

This work has been supported by the AEROTRAIN Marie Skłodowska-Curie project (MSCA-ITN-2020-953454), funded by the European Commission, and by the MARTIN project (PID2022-143267OB-100) funded by MCIN/AEI/ 10.13039/501100011033 and by ERDF A way of making Europe.

A. Gonzalez-Morgado (corresponding author: mantonio@us.es.), G. Heredia and A. Ollero are with the GRVC Robotics Lab, Escuela Tecnica Superior de Ingenieria, Universidad de Sevilla, 41012 Sevilla, Spain.

E. Cuniato, and R. Siegwart are with the Autonomous Systems Lab (ASL), ETH Zurich.

M. Tognon is with the Univ Rennes, CNRS, Inria, IRISA - UMR 6074, F-35000 Rennes.

generate thrust in one direction only, and three independent torques to control the 6 DoFs of the system. This is very limiting when they physically interact with the environment [7], since they must modify their attitude to generate lateral forces that allow contact with the environment. To overcome these limitations, some works present new platform designs [8], ranging from the use of tilted rotor platforms [6], [9] or tilting rotors with servomotors [10], [11], to the use of flapping bio-inspired robots [12], [13].

All these developments have enabled the use of aerial manipulators in tasks requiring high precision contact. For example, in [6], the authors present a fully-actuated design for the inspection of hard-to-reach bridges' areas, like ceilings, bearings and beams [14]. The range of tasks recently extended to the manipulation of articulated objects like doors [15], and carts [16], as well as the interaction with humans [17].

In view of the recent improvement of manipulation capabilities of aerial robots, the use of aerial manipulators for physical interaction with flexible elements, such as trees, branches or beams, has been investigated. This field of study presents novel applications, such as the use of aerial manipulators for the collection of environmental DNA from tree branches

<footnote n="1">¹Video of the paper: https://youtu.be/qFYHw3_6gqk</footnote>
# Controlled Shaking Of Trees With An Aerial Manipulator

A. Gonzalez-Morgado, *Graduate Student Member, IEEE*, E. Cuniato, *Graduate Student Member, IEEE*, M. Tognon, *Member, IEEE*, G. Heredia, *Member, IEEE*, R. Siegwart and A. Ollero, *Fellow, IEEE*

*Abstract*—In recent years, the fields of application of aerial manipulators have expanded, ranging from infrastructure inspection to physical interaction with flexible elements such as branches and trees. This paper presents the controlled shaking of a tree with an aerial manipulator. Our work aims at contributing to applications like the identification of tree parameters for environmental health monitoring or the collection of samples and fruits by vibration. To this end, we propose a control strategy for controlled shaking of flexible systems. We adopt a self-excited oscillation strategy that induces vibrations at the natural frequency of the system, at which the greatest amplification and therefore the greatest vibrations occur. Likewise, this work presents a simplified 1 DoF model based on the Rayleigh-Ritz method to analyze dynamic interaction between a tree and the aerial manipulator with the controlled shaking strategy. The proposed control strategy is evaluated through indoor experiments, where an aerial manipulator shakes an indoor tree made of bamboo canes. Experimental results show how the proposed model can estimate properly the amplitude of the vibration and the frequency of the vibration, depending on the grasping point and the control gain of the self-excited oscillation strategy<fnref n="1" />.

*Index Terms*—Aerial Robotics, Aerial Physical Interaction, Control for Manipulation, Control of Flexible Elements, Controlled Shaking

<figure><img src="image_2.png" /><figcaption>Fig. 1. Aerial manipulator controlled shakes a tree. Different frames during the oscillation of the system are shown.</figcaption></figure>

## I. INTRODUCTION

AERIAL manipulators bring together the flexibility and maneuverability of aerial platforms with the functionalities and capabilities of robotic manipulators [1]. They offer a chance to diminish the duration, expenses, and potential hazards associated with operations conducted at elevated, dangerous, or hard-to-reach areas. This led to the expansion of application fields for UAVs (Unmanned Aerial Vehicles) which range from support in post-disaster situations [2] to the inspection of oil and gas refineries [3], chemical plants [4], power lines [5], or civil infrastructures [6]. This shows the huge recent development of aerial manipulators.

However, some challenges need to be addressed to develop secure, robust and outdoor capable aerial platforms. Standard coplanar multirotors are underactuated systems, since they can

This work has been supported by the AEROTRAIN Marie Skłodowska-Curie project (MSCA-ITN-2020-953454), funded by the European Commission, and by the MARTIN project (PID2022-143267OB-100) funded by MCIN/AEI/ 10.13039/501100011033 and by ERDF A way of making Europe.

A. Gonzalez-Morgado (corresponding author: mantonio@us.es.), G. Heredia and A. Ollero are with the GRVC Robotics Lab, Escuela Tecnica Superior de Ingenieria, Universidad de Sevilla, 41012 Sevilla, Spain.

E. Cuniato, and R. Siegwart are with the Autonomous Systems Lab (ASL), ETH Zurich.

M. Tognon is with the Univ Rennes, CNRS, Inria, IRISA - UMR 6074, F-35000 Rennes.

generate thrust in one direction only, and three independent torques to control the 6 DoFs of the system. This is very limiting when they physically interact with the environment [7], since they must modify their attitude to generate lateral forces that allow contact with the environment. To overcome these limitations, some works present new platform designs [8], ranging from the use of tilted rotor platforms [6], [9] or tilting rotors with servomotors [10], [11], to the use of flapping bio-inspired robots [12], [13].

All these developments have enabled the use of aerial manipulators in tasks requiring high precision contact. For example, in [6], the authors present a fully-actuated design for the inspection of hard-to-reach bridges' areas, like ceilings, bearings and beams [14]. The range of tasks recently extended to the manipulation of articulated objects like doors [15], and carts [16], as well as the interaction with humans [17].

In view of the recent improvement of manipulation capabilities of aerial robots, the use of aerial manipulators for physical interaction with flexible elements, such as trees, branches or beams, has been investigated. This field of study presents novel applications, such as the use of aerial manipulators for the collection of environmental DNA from tree branches

<footnote n="1">¹Video of the paper: https://youtu.be/qFYHw3_6gqk</footnote>
[18] or for canopy sampling [19], which helps biodiversity monitoring. Several works exploit the unlimited working space of aerial manipulators to install sensors [20] in trees and plantations, allowing data collection and monitoring to detect early plant stress [21]. In addition, innovative designs have emerged, featuring robots capable of moving within tree canopies [22], where UAVs cannot access due to densely vegetated environments filled with leaves and branches. Other studies focus on the new challenges that open up this field of research, like the coupling dynamics that appear between the aerial platform and the flexible system [23].

Furthermore, in the field of agriculture and environmental science, various applications necessitate the vibration of trees and plants. For instance, the shake-and-catch method [24] is one of the most common methods of fruit harvesting. Here the tree is shaken using a shaker, a mechanical device that induce vibrations to the tree [25]. The tree shaking kinetic energy is transmitted to branches and fruits, which generates enough force to break the fruit-stem interface. These fruits are collected afterwards or fall directly onto a net that is previously placed around the tree. Typically, this method uses a large and powerful shaker, mounted on a tractor, that shakes at a frequency independent of the tree. When the shaker is smaller, with a corresponding reduction in power, it can be operated by humans. Nevertheless, due to the power limits of the shaker, it is important to shake the tree at its natural frequency, at which the greatest amplification and therefore the greatest vibrations occur [26]. Therefore, shakers capable of identifying these natural frequencies have been investigated in order to adapt the shaking frequency during operation [27]. In the same way, many methods for identifying dynamic parameters of trees need to induce vibrations [28] into the trees. The identified parameters (natural frequency, damping, stiffness, ...) are used to monitor the condition of trees and they can help in the detection of possible damages by tracking variations in the characteristic frequency or damping [29]. However, the above applications require large shakers operated by people, who must move to the tree to be vibrated for fruit harvesting or identifying its parameters.

## A. Contributions and novelties

In this work we present the first aerial manipulator capable of controlled shaking flexible elements such as trees or branches (see Figure 1), contributing to the field of agriculture and environmental science. This study opens the door to the use of UAVs for fruit harvesting or the assessment of tree parameters, enabling the remote execution of these tasks without the direct involvement of humans or big machines, like shakers, in areas that are inaccessible to human intervention. We propose the following contributions:

* Design of a control scheme for the UAV to shake flexible systems, like trees or branches, in a controlled way.

* Development of a 1 DoF model for analyzing the closed-loop system behavior when the robot is controlled by the proposed method.

To the best of our knowledge, this is the first work that proposes a control scheme for an aerial manipulator shaking

a tree in a controlled way. The proposed control strategy relies on a self-excited oscillation strategy [30]. This approach generates a force command for the UAV based on its velocity, guaranteeing that the tree-UAV system vibrates at its natural frequency, without requiring knowledge of any parameters of the shaken system, such as length, stiffness, damping, or equivalent mass. We have combined this strategy with an hybrid force/pose controller [11], such to generate oscillations in the interaction direction while controlling the position in the others.

In addition, although there are already several tree models developed in biomechanics fields [31], they are typically too complex to allow a theoretical analysis of the tree-UAV interaction during the controlled shaking. The proposed model in this work is a simple 1 DoF model based on the Rayleigh-Ritz method [32], which allows to reduce an infinite number of degrees-of-freedom of a system into a finite number [33]. Compared with other UAV-flexible elements models [19], [23], our model can easily consider the controlled shaking scheme, making analysis possible and easier. Despite its simplicity, it still predicts the behaviour of the complete tree-UAV system considering the controller gain, in particular the oscillation's frequency and amplitude.

The proposed controlled shaking scheme is evaluated through indoor experiments using a fully-actuated UAV and a custom-made tree. While the control strategy can be validated with any flexible system as it does not require knowledge of any parameters of the interacting system, this setup enables better control of the associated parameters. This allows us to precisely compare the real experimental results with the expected ones, which are obtained by combining the proposed tree model with the control scheme. This comparison shows that the system vibrates at its natural frequency.

The rest of the paper is organised as follows. Section II describes and provides models for the fully-actuated UAV and the tree. Section III presents the control of the UAV, the self-excited oscillation strategy and the theoretical analysis of the tree-UAV system using the proposed method. The experimental results are presented in Section IV, where the model and the self-excited oscillation strategy are evaluated with indoor flights. Finally, the conclusions of this work, along with future research directions, are summarized in Section V.

## II. UAV-TREE SYSTEM MODELING

### A. Platform Description

The fully-actuated UAV used in this work is an Omnidirectional Micro Aerial Vehicle (OMAV) as in [11]. It is an aerial robot composed of two main elements: the main body in the middle that contains most of the system's inertia, and six tiltrotor units that can rotate with respect to the main body. Each unit is composed of a servomotor for the arm rotation, and two coaxial propellers to provide the desired thrust. Each of the arms can be independently actuated, allowing the platform to instantaneously exert a full 6 DoFs control wrench at any orientation, thus making the system fully-actuated. We additionally rigidly attach a carbon fiber stick with a hook to the main body of the vehicle such to enable it to interact with
[18] or for canopy sampling [19], which helps biodiversity monitoring. Several works exploit the unlimited working space of aerial manipulators to install sensors [20] in trees and plantations, allowing data collection and monitoring to detect early plant stress [21]. In addition, innovative designs have emerged, featuring robots capable of moving within tree canopies [22], where UAVs cannot access due to densely vegetated environments filled with leaves and branches. Other studies focus on the new challenges that open up this field of research, like the coupling dynamics that appear between the aerial platform and the flexible system [23].

Furthermore, in the field of agriculture and environmental science, various applications necessitate the vibration of trees and plants. For instance, the shake-and-catch method [24] is one of the most common methods of fruit harvesting. Here the tree is shaken using a shaker, a mechanical device that induce vibrations to the tree [25]. The tree shaking kinetic energy is transmitted to branches and fruits, which generates enough force to break the fruit-stem interface. These fruits are collected afterwards or fall directly onto a net that is previously placed around the tree. Typically, this method uses a large and powerful shaker, mounted on a tractor, that shakes at a frequency independent of the tree. When the shaker is smaller, with a corresponding reduction in power, it can be operated by humans. Nevertheless, due to the power limits of the shaker, it is important to shake the tree at its natural frequency, at which the greatest amplification and therefore the greatest vibrations occur [26]. Therefore, shakers capable of identifying these natural frequencies have been investigated in order to adapt the shaking frequency during operation [27]. In the same way, many methods for identifying dynamic parameters of trees need to induce vibrations [28] into the trees. The identified parameters (natural frequency, damping, stiffness, ...) are used to monitor the condition of trees and they can help in the detection of possible damages by tracking variations in the characteristic frequency or damping [29]. However, the above applications require large shakers operated by people, who must move to the tree to be vibrated for fruit harvesting or identifying its parameters.

## A. Contributions and novelties

In this work we present the first aerial manipulator capable of controlled shaking flexible elements such as trees or branches (see Figure 1), contributing to the field of agriculture and environmental science. This study opens the door to the use of UAVs for fruit harvesting or the assessment of tree parameters, enabling the remote execution of these tasks without the direct involvement of humans or big machines, like shakers, in areas that are inaccessible to human intervention. We propose the following contributions:

* Design of a control scheme for the UAV to shake flexible systems, like trees or branches, in a controlled way.

* Development of a 1 DoF model for analyzing the closed-loop system behavior when the robot is controlled by the proposed method.

To the best of our knowledge, this is the first work that proposes a control scheme for an aerial manipulator shaking

a tree in a controlled way. The proposed control strategy relies on a self-excited oscillation strategy [30]. This approach generates a force command for the UAV based on its velocity, guaranteeing that the tree-UAV system vibrates at its natural frequency, without requiring knowledge of any parameters of the shaken system, such as length, stiffness, damping, or equivalent mass. We have combined this strategy with an hybrid force/pose controller [11], such to generate oscillations in the interaction direction while controlling the position in the others.

In addition, although there are already several tree models developed in biomechanics fields [31], they are typically too complex to allow a theoretical analysis of the tree-UAV interaction during the controlled shaking. The proposed model in this work is a simple 1 DoF model based on the Rayleigh-Ritz method [32], which allows to reduce an infinite number of degrees-of-freedom of a system into a finite number [33]. Compared with other UAV-flexible elements models [19], [23], our model can easily consider the controlled shaking scheme, making analysis possible and easier. Despite its simplicity, it still predicts the behaviour of the complete tree-UAV system considering the controller gain, in particular the oscillation's frequency and amplitude.

The proposed controlled shaking scheme is evaluated through indoor experiments using a fully-actuated UAV and a custom-made tree. While the control strategy can be validated with any flexible system as it does not require knowledge of any parameters of the interacting system, this setup enables better control of the associated parameters. This allows us to precisely compare the real experimental results with the expected ones, which are obtained by combining the proposed tree model with the control scheme. This comparison shows that the system vibrates at its natural frequency.

The rest of the paper is organised as follows. Section II describes and provides models for the fully-actuated UAV and the tree. Section III presents the control of the UAV, the self-excited oscillation strategy and the theoretical analysis of the tree-UAV system using the proposed method. The experimental results are presented in Section IV, where the model and the self-excited oscillation strategy are evaluated with indoor flights. Finally, the conclusions of this work, along with future research directions, are summarized in Section V.

## II. UAV-TREE SYSTEM MODELING

### A. Platform Description

The fully-actuated UAV used in this work is an omnidirectional Micro Aerial Vehicle (OMAV) as in [11]. It is an aerial robot composed of two main elements: the main body in the middle that contains most of the system's inertia, and six tiltrotor units that can rotate with respect to the main body. Each unit is composed of a servomotor for the arm rotation, and two coaxial propellers to provide the desired thrust. Each of the arms can be independently actuated, allowing the platform to instantaneously exert a full 6 DoFs control wrench at any orientation, thus making the system fully-actuated. We additionally rigidly attach a carbon fiber stick with a hook to the main body of the vehicle such to enable it to interact with
the tree. The interacting hook is the result of a 3D-printed design and is completely passive.

<figure><img src="image_3.png" /><figcaption>Fig. 2. OMAV platform with hook for the controlled shaking of trees. The body frame of OMAV is indicated with $\mathcal{F}_B$, with the colored propellers mounted on the arm along the body y-axis. A zoom-in of the hook better shows its details.</figcaption></figure>

## B. OMAV Model

We model our system following a Lagrangian approach as in [11]. We define a world and a body frame denoted by $\mathcal{F}_W$ and $\mathcal{F}_B$, respectively. $\mathcal{F}_B$ is located at the OMAV's center of mass. We then define the position $w \in \mathbb{R}^3$ and orientation $w \mathbf{R}_B \in \text{SO}(3)$ of the OMAV's body with respect to $\mathcal{F}_W$. From these we can define the linear $b \dot{\mathbf{p}} \in \mathbb{R}^3$ and angular $b \omega \in \mathbb{R}^3$ velocities, which we stack together to from the system's twist $b \mathbf{v} \in \mathbb{R}^6$. We also assume the platform is always able to generated a full control wrench $\tau_c \in \mathbb{R}^6$ through its propellers and servomotors, following the control allocation procedure in [11]. In the end, we obtain the Lagrangian model for our OMAV as

$$M_B \dot{\mathbf{v}} + \mathbf{c}(b\omega) + \mathbf{g}(w \mathbf{R}_B) = \tau_c + \tau_{\text{ext}}, \quad (1)$$

where $M_B \in \mathbb{R}^{6 \times 6}$ is the constant, body-frame, inertia matrix, $c(b\omega) \in \mathbb{R}^6$ is the vector of Coriolis forces, $g(w \mathbf{R}_B) \in \mathbb{R}^6$ is the gravity wrench and $\tau_{\text{ext}} \in \mathbb{R}^6$ represents the external wrench applied on the system because of its interaction with the environment.

1) **External Wrench Estimator:** To successfully interact with the tree, we require a measure or estimate of the interaction wrench. To this end, we employ the momentum-based estimator proposed in [11], [34] as

$$\hat{\tau}_{\text{ext}} = \mathbf{K}_E \left( M_B \mathbf{v} - \int (\tau_c - \mathbf{c} - \mathbf{g} + \hat{\tau}_{\text{ext}}) dt \right), \quad (2)$$

where $\mathbf{K}_E \in \mathbb{R}^{6 \times 6}$ is a positive definite and diagonal matrix which represents the estimator gain and $\hat{\tau}_{\text{ext}} \in \mathbb{R}^6$ represents the estimated wrench.

## C. Tree Model

In this work we use the Rayleigh-Ritz method to obtain a 1 DoF model of the tree. We choose as our only DoF $x_a(t) \in \mathbb{R}$, the displacement of the contact point between the UAV and the tree $a \in \mathbb{R}$.

**Assumption 1.** Although the tree is a continuous system, it can be represented by a single degree of freedom during the

interaction with the UAV controlled by the proposed strategy. This assumption is based on the fact that the closed-loop system (UAV with the tree) behavior mainly vibrates in its first mode of vibration, so only one harmonic function is needed to represent it, i.e., only one degree of freedom. This assumption is validated by the experimental results presented in Section IV-B.

Without loss of generality, we consider that the oscillation occurs along the $x_W$ direction, since we can redefine the orientation of $\mathcal{F}_W$ to make $x_W$ and $x_B$ collinear. We consider that the world frame is attached to the tree base, with $y_W$ in the direction of the resting tree trunk, as shown in Figure 3. Using the Lagrange's equation, the dynamics of the system reads:

$$m^{\text{eq}} \ddot{x}_a(t) + c^{\text{eq}} \dot{x}_a(t) + k^{\text{eq}} x_a(t) = f_x(t), \quad (3)$$

where $x_a(t)$ is the displacement of the contact point between the UAV and the tree trunk, $f_x(t) \in \mathbb{R}$ is the force applied by the UAV onto the tree trunk in the $x_W$ direction, and $m^{\text{eq}} \in \mathbb{R}$, $c^{\text{eq}} \in \mathbb{R}$ and $k^{\text{eq}} \in \mathbb{R}$ are respectively the equivalent mass, equivalent damping and equivalent stiffness of the tree. So, in order to analyze and estimate the behaviour between the tree and the UAV, we need to compute the equivalent parameters of the model $m^{\text{eq}}$, $c^{\text{eq}}$ and $k^{\text{eq}}$.

**Assumption 2.** We study the model describing the tree-UAV interaction when this happens at the trunk level$^2$. In practice, most applications such as fruit harvesting and parameter monitoring involve shaking the tree from the trunk.

The parameters $m^{\text{eq}}$ and $k^{\text{eq}}$ are obtained from the total kinetic energy $T \in \mathbb{R}$ and potential energy $K \in \mathbb{R}$ of the system. The equivalent damping $c^{\text{eq}}$ is instead obtained experimentally, as its analytical estimation can be very challenging because there are many sources of dissipation during the motion [35]. For obtaining the kinetic energy $T$ and potential energy $K$, we divide the tree into three different components (see Figure 3): the trunk, the $N_B$ branches connected to the trunk and the $N_F$ fruits, where $N_B \in \mathbb{N}$ and $N_F \in \mathbb{N}$ are the number of branches and fruits of the tree, respectively. If a branch has branches, these can be treated as fruit at their point of attachment, following a process similar to that described in Section II-C2. This way, the kinetic $T$ and potential energy $K$ of the complete tree become:

$$T = T_T + \sum_{i=1}^{N_B} T_{B_i} + \sum_{j=1}^{N_F} T_{F_j} \quad \text{and} \quad K = K_T, \quad (4)$$

where $T_{B_i} \in \mathbb{R}$ is the kinetic energy of the $i$-th branch $B_i$, $T_{F_j} \in \mathbb{R}$ is the kinetic energy of the $j$-th fruit $F_j$, and $T_T \in \mathbb{R}$ and $K_T \in \mathbb{R}$ are the kinetic energy and potential energy of the trunk, respectively. In this case, only the elastic potential energy of the trunk is considered since it is the only element that deforms. In addition, we do not consider the gravitational potential energy, as the movement is perpendicular to the

<footnote>$^2$The proposed model could also be applied to interactions with branches if the trunk is highly rigid and can be considered as a fixed base.</footnote>
the tree. The interacting hook is the result of a 3D-printed design and is completely passive.

<figure><img src="image_3.png" /><figcaption>Fig. 2. OMAV platform with hook for the controlled shaking of trees. The body frame of OMAV is indicated with $\mathcal{F}_B$, with the colored propellers mounted on the arm along the body y-axis. A zoom-in of the hook better shows its details.</figcaption></figure>

## B. OMAV Model

We model our system following a Lagrangian approach as in [11]. We define a world and a body frame denoted by $\mathcal{F}_W$ and $\mathcal{F}_B$, respectively. $\mathcal{F}_B$ is located at the OMAV's center of mass. We then define the position $w \in \mathbb{R}^3$ and orientation $w \mathbf{R}_B \in \text{SO}(3)$ of the OMAV's body with respect to $\mathcal{F}_W$. From these we can define the linear $b \dot{\mathbf{p}} \in \mathbb{R}^3$ and angular $b \omega \in \mathbb{R}^3$ velocities, which we stack together to from the system's twist $b \mathbf{v} \in \mathbb{R}^6$. We also assume the platform is always able to generated a full control wrench $\tau_c \in \mathbb{R}^6$ through its propellers and servomotors, following the control allocation procedure in [11]. In the end, we obtain the Lagrangian model for our OMAV as

$$M_B \dot{\mathbf{v}} + \mathbf{c}(b\omega) + \mathbf{g}(w \mathbf{R}_B) = \tau_c + \tau_{\text{ext}}, \quad (1)$$

where $M_B \in \mathbb{R}^{6 \times 6}$ is the constant, body-frame, inertia matrix, $\mathbf{c}(b\omega) \in \mathbb{R}^6$ is the vector of Coriolis forces, $\mathbf{g}(w \mathbf{R}_B) \in \mathbb{R}^6$ is the gravity wrench and $\tau_{\text{ext}} \in \mathbb{R}^6$ represents the external wrench applied on the system because of its interaction with the environment.

1) External Wrench Estimator: To successfully interact with the tree, we require a measure or estimate of the interaction wrench. To this end, we employ the momentum-based estimator proposed in [11], [34] as

$$\hat{\tau}_{\text{ext}} = \mathbf{K}_E \left( \mathbf{M}_B \mathbf{v} - \int (\tau_c - \mathbf{c} - \mathbf{g} + \hat{\tau}_{\text{ext}}) dt \right), \quad (2)$$

where $\mathbf{K}_E \in \mathbb{R}^{6 \times 6}$ is a positive definite and diagonal matrix which represents the estimator gain and $\hat{\tau}_{\text{ext}} \in \mathbb{R}^6$ represents the estimated wrench.

## C. Tree Model

In this work we use the Rayleigh-Ritz method to obtain a 1 DoF model of the tree. We choose as our only DoF $x_a(t) \in \mathbb{R}$, the displacement of the contact point between the UAV and the tree $a \in \mathbb{R}$.

**Assumption 1.** Although the tree is a continuous system, it can be represented by a single degree of freedom during the

interaction with the UAV controlled by the proposed strategy. This assumption is based on the fact that the closed-loop system (UAV with the tree) behavior mainly vibrates in its first mode of vibration, so only one harmonic function is needed to represent it, i.e., only one degree of freedom. This assumption is validated by the experimental results presented in Section IV-B.

Without loss of generality, we consider that the oscillation occurs along the $\mathbf{x}_W$ direction, since we can redefine the orientation of $\mathcal{F}_W$ to make $\mathbf{x}_W$ and $\mathbf{x}_B$ collinear. We consider that the world frame is attached to the tree base, with $y_W$ in the direction of the resting tree trunk, as shown in Figure 3. Using the Lagrange's equation, the dynamics of the system reads:

$$m^{\text{eq}} \ddot{x}_a(t) + c^{\text{eq}} \dot{x}_a(t) + k^{\text{eq}} x_a(t) = f_x(t), \quad (3)$$

where $x_a(t)$ is the displacement of the contact point between the UAV and the tree trunk, $f_x(t) \in \mathbb{R}$ is the force applied by the UAV onto the tree trunk in the $\mathbf{x}_W$ direction, and $m^{\text{eq}} \in \mathbb{R}$, $c^{\text{eq}} \in \mathbb{R}$ and $k^{\text{eq}} \in \mathbb{R}$ are respectively the equivalent mass, equivalent damping and equivalent stiffness of the tree. So, in order to analyze and estimate the behaviour between the tree and the UAV, we need to compute the equivalent parameters of the model $m^{\text{eq}}$, $c^{\text{eq}}$ and $k^{\text{eq}}$.

**Assumption 2.** We study the model describing the tree-UAV interaction when this happens at the trunk level$^2$. In practice, most applications such as fruit harvesting and parameter monitoring involve shaking the tree from the trunk.

The parameters $m^{\text{eq}}$ and $k^{\text{eq}}$ are obtained from the total kinetic energy $T \in \mathbb{R}$ and potential energy $K \in \mathbb{R}$ of the system. The equivalent damping $c^{\text{eq}}$ is instead obtained experimentally, as its analytical estimation can be very challenging because there are many sources of dissipation during the motion [35]. For obtaining the kinetic energy $T$ and potential energy $K$, we divide the tree into three different components (see Figure 3): the trunk, the $N_B$ branches connected to the trunk and the $N_F$ fruits, where $N_B \in \mathbb{N}$ and $N_F \in \mathbb{N}$ are the number of branches and fruits of the tree, respectively. If a branch has branches, these can be treated as fruit at their point of attachment, following a process similar to that described in Section II-C2. This way, the kinetic $T$ and potential energy $K$ of the complete tree become:

$$T = T_T + \sum_{i=1}^{N_B} T_{B_i} + \sum_{j=1}^{N_F} T_{F_j} \quad \text{and} \quad K = K_T, \quad (4)$$

where $T_{B_i} \in \mathbb{R}$ is the kinetic energy of the $i$-th branch $B_i$, $T_{F_j} \in \mathbb{R}$ is the kinetic energy of the $j$-th fruit $F_j$, and $T_T \in \mathbb{R}$ and $K_T \in \mathbb{R}$ are the kinetic energy and potential energy of the trunk, respectively. In this case, only the elastic potential energy of the trunk is considered since it is the only element that deforms. In addition, we do not consider the gravitational potential energy, as the movement is perpendicular to the

<footnote>$^2$The proposed model could also be applied to interactions with branches if the trunk is highly rigid and can be considered as a fixed base.</footnote>
<figure><img src="image_4.png" /><figcaption>Fig. 3. Schematic of the tree, including the trunk, branches and fruits. a) Real trees: 1) trunk with low stiffness, where the model is only applied for UAV–trunk interaction and 2) trunk with high stiffness, where the model is also valid for UAV–branch interaction, as the trunk can be considered as fixed base. b) Simplified model for the tree, composed of a trunk, branches and fruits. c) Deformed trunk with an applied force $f_x$. d) Branch parameters. e) Detailed view of a branch.</figcaption></figure>

gravity. Thus, the equivalent mass $m^{eq}$ and stiffness $k^{eq}$ become:

$$m^{eq} = m_T^{eq} + \sum_{i=1}^{N_B} m_{B_i}^{eq} + \sum_{j=1}^{N_F} m_{F_j}^{eq} \quad \text{and} \quad k^{eq} = k_T^{eq}, \quad (5)$$

where $m_{B_i}^{eq} \in \mathbb{R}$ is the equivalent mass of the i-th branch $B_i$, $m_{F_j}^{eq} \in \mathbb{R}$ is the equivalent mass of the j-th fruit $F_j$, and $m_T^{eq} \in \mathbb{R}$ and $k_T^{eq} \in \mathbb{R}$ are the equivalent mass and stiffness of the trunk, respectively. The equivalent mass and stiffness of each element can be identified from the kinetic and potential energy functions, expressed as functions of the degree of freedom $x_a(t)$. We compute these equivalent parameters in the following.

1) Trunk parameters: the trunk is modeled as a vertical beam of length $l_T \in \mathbb{R}$.

**Assumption 3.** For simplicity, we suppose that the tree is composed of straight branches and a straight trunk. If this is not the case, the curves describing the shape of the branches and trunk have to be parameterized to accurately represent their form. In practice, various state-of-the-art methods [36], [37] can be employed for this purpose. The morphology of the tree would be considered in the calculation of the kinetic and potential energy of each element, that is in the integrals (6) and (7).

The displacement of any point of the trunk is given by $x(\zeta, t) \in \mathbb{R}$, where $\zeta \in \mathbb{R}$ is the spatial variable that goes from the ground to the tip of the trunk $\zeta \in [0, l_T]$. We consider the UAV applying a force $f_x$ at the point a. The kinetic energy $T_T$ and the potential energy $K_T$ of the trunk are:

$$T_T = \frac{1}{2} \int_{0}^{l_T} \rho_T \dot{x}(\zeta, t)^2 d\zeta \quad (6)$$

$$K_T = \frac{1}{2} \int_{0}^{l_T} EI (x''(\zeta, t))^2 d\zeta, \quad (7)$$

where $\rho_T \in \mathbb{R}$ is the linear density of the trunk, $E \in \mathbb{R}$ is Young's modulus, $I \in \mathbb{R}$ is the section inertia, $\dot{x}(\zeta, t) = dx(\zeta, t)/dt$ and $x''(\zeta, t) = \partial^2 x(\zeta, t)/\partial\zeta^2$. For simplicity, we assume that $\rho_T$, $E$ and $I$ are constant along

the length of the trunk. Likewise, we assume that $x(\zeta, t)$ can be expressed in separable variables<fnref n="3" /> as:

$$x(\zeta, t) = x_a(t)\eta(\zeta) \quad \forall \zeta \in [0, l_T], \quad (8)$$

where $x_a(t)$ is the displacement of point a, chosen as the degree of freedom of the system, and $\eta(\zeta) \in \mathbb{R}$ is the shape of the deformed beam. Doing that, we can identify the expressions of the equivalent mass $m_T^{eq}$ and equivalent stiffness $k_T^{eq}$ as:

$$T_T = \frac{1}{2} \dot{x}_a^2 \int_{0}^{l_T} \rho_T (\eta(\zeta))^2 d\zeta, \quad (9)$$

$$K_T = \frac{1}{2} x_a^2 \underbrace{\int_{0}^{l_T} EI (\eta''(\zeta))^2 d\zeta}_{k_T^{eq}}. \quad (10)$$

Moreover, the expression of the deformation $\eta(\zeta)$ must satisfy the Dirichlet boundary conditions, related with the displacement of the beam $\eta(0) = 0$, and Neumann boundary conditions, related with the rotation of the beam $\eta'(0) = 0$ [38]. As we consider a 1 DoF model, a good approximation is to take $\eta(\zeta)$ as the deformed shape of the beam under a static load applied at the point a [39]. Doing this, $\eta(\zeta)$ becomes:

$$\eta(\zeta) = \begin{cases} \frac{3}{2} \left(\frac{\zeta}{a}\right)^2 - \frac{1}{2} \left(\frac{\zeta}{a}\right)^3 & \text{if } 0 \le \zeta \le a \\ \frac{3}{2} \left(\frac{\zeta}{a}\right) - \frac{1}{2} & \text{if } a \le \zeta \le l_T \end{cases}. \quad (11)$$

From the expression of $\eta(\zeta)$ in (11) and the expression of $m_T^{eq}$ in (9) and $k_T^{eq}$ in (10), we obtain:

$$m_T^{eq} = \rho_T \left[ a \frac{3}{140} + \frac{3}{4} \left( \frac{l_T^3 - a^3}{a^2} \right) + \frac{1}{4} (l_T - a) - \frac{3}{4} \left( \frac{l_T^2 - a^2}{a} \right) \right] \quad (12)$$

<footnote n="3">³This assumption is equivalent to decompose $x(\zeta, t)$ as the sum of its different vibration modes and only keep the first mode. We do this as our model is a 1 DoF model, so we consider only one vibration mode.</footnote>
<figure><img src="image_4.png" /><figcaption>Fig. 3. Schematic of the tree, including the trunk, branches and fruits. a) Real trees: 1) trunk with low stiffness, where the model is only applied for UAV–trunk interaction and 2) trunk with high stiffness, where the model is also valid for UAV–branch interaction, as the trunk can be considered as fixed base. b) Simplified model for the tree, composed of a trunk, branches and fruits. c) Deformed trunk with an applied force $f_x$. d) Branch parameters. e) Detailed view of a branch.</figcaption></figure>

gravity. Thus, the equivalent mass $m^{eq}$ and stiffness $k^{eq}$ become:

$$m^{eq} = m_T^{eq} + \sum_{i=1}^{N_B} m_{B_i}^{eq} + \sum_{j=1}^{N_F} m_{F_j}^{eq} \quad \text{and} \quad k^{eq} = k_T^{eq}, \quad (5)$$

where $m_{B_i}^{eq} \in \mathbb{R}$ is the equivalent mass of the i-th branch $B_i$, $m_{F_j}^{eq} \in \mathbb{R}$ is the equivalent mass of the j-th fruit $F_j$, and $m_T^{eq} \in \mathbb{R}$ and $k_T^{eq} \in \mathbb{R}$ are the equivalent mass and stiffness of the trunk, respectively. The equivalent mass and stiffness of each element can be identified from the kinetic and potential energy functions, expressed as functions of the degree of freedom $x_a(t)$. We compute these equivalent parameters in the following.

1) Trunk parameters: the trunk is modeled as a vertical beam of length $l_T \in \mathbb{R}$.

**Assumption 3.** For simplicity, we suppose that the tree is composed of straight branches and a straight trunk. If this is not the case, the curves describing the shape of the branches and trunk have to be parameterized to accurately represent their form. In practice, various state-of-the-art methods [36], [37] can be employed for this purpose. The morphology of the tree would be considered in the calculation of the kinetic and potential energy of each element, that is in the integrals (6) and (7).

The displacement of any point of the trunk is given by $x(\zeta, t) \in \mathbb{R}$, where $\zeta \in \mathbb{R}$ is the spatial variable that goes from the ground to the tip of the trunk $\zeta \in [0, l_T]$. We consider the UAV applying a force $f_x$ at the point a. The kinetic energy $T_T$ and the potential energy $K_T$ of the trunk are:

$$T_T = \frac{1}{2} \int_{0}^{l_T} \rho_T \dot{x}(\zeta, t)^2 d\zeta \quad (6)$$

$$K_T = \frac{1}{2} \int_{0}^{l_T} EI (x''(\zeta, t))^2 d\zeta, \quad (7)$$

where $\rho_T \in \mathbb{R}$ is the linear density of the trunk, $E \in \mathbb{R}$ is Young's modulus, $I \in \mathbb{R}$ is the section inertia, $\dot{x}(\zeta, t) = dx(\zeta, t)/dt$ and $x''(\zeta, t) = \partial^2 x(\zeta, t)/\partial\zeta^2$. For simplicity, we assume that $\rho_T$, $E$ and $I$ are constant along

the length of the trunk. Likewise, we assume that $x(\zeta, t)$ can be expressed in separable variables<fnref n="3" /> as:

$$x(\zeta, t) = x_a(t)\eta(\zeta) \quad \forall \zeta \in [0, l_T], \quad (8)$$

where $x_a(t)$ is the displacement of point a, chosen as the degree of freedom of the system, and $\eta(\zeta) \in \mathbb{R}$ is the shape of the deformed beam. Doing that, we can identify the expressions of the equivalent mass $m_T^{eq}$ and equivalent stiffness $k_T^{eq}$ as:

$$T_T = \frac{1}{2} \dot{x}_a^2 \int_0^{l_T} \rho_T (\eta(\zeta))^2 d\zeta, \quad (9)$$

$$K_T = \frac{1}{2} x_a^2 \underbrace{\int_0^{l_T} EI (\eta''(\zeta))^2 d\zeta}_{k_T^{eq}}. \quad (10)$$

Moreover, the expression of the deformation $\eta(\zeta)$ must satisfy the Dirichlet boundary conditions, related with the displacement of the beam $\eta(0) = 0$, and Neumann boundary conditions, related with the rotation of the beam $\eta'(0) = 0$ [38]. As we consider a 1 DoF model, a good approximation is to take $\eta(\zeta)$ as the deformed shape of the beam under a static load applied at the point a [39]. Doing this, $\eta(\zeta)$ becomes:

$$\eta(\zeta) = \begin{cases} \frac{3}{2} \left(\frac{\zeta}{a}\right)^2 - \frac{1}{2} \left(\frac{\zeta}{a}\right)^3 & \text{if } 0 \le \zeta \le a \\ \frac{3}{2} \left(\frac{\zeta}{a}\right) - \frac{1}{2} & \text{if } a \le \zeta \le l_T \end{cases}. \quad (11)$$

From the expression of $\eta(\zeta)$ in (11) and the expression of $m_T^{eq}$ in (9) and $k_T^{eq}$ in (10), we obtain:

$$m_T^{eq} = \rho_T \left[ a \frac{3}{140} + \frac{3}{4} \left( \frac{l_T^3 - a^3}{a^2} \right) + \frac{1}{4} (l_T - a) - \frac{3}{4} \left( \frac{l_T^2 - a^2}{a} \right) \right] \quad (12)$$

<footnote n="3">³This assumption is equivalent to decompose $x(\zeta, t)$ as the sum of its different vibration modes and only keep the first mode. We do this as our model is a 1 DoF model, so we consider only one vibration mode.</footnote>
$$k_T^{\text{eq}} = \frac{3EI}{a^3}, \quad (13)$$

which show the dependence with the grasping point a.

2) Branch parameters: the branches are modeled as straight beams (see Assumption 3) embedded in the trunk of the tree, which do not deform, but are moved only by the movement of the trunk. As shown in Figures 3d and 3e, the branch $B_i$, with length $l_{B_i} \in \mathbb{R}$, is embedded in the trunk at the trunk point $R_i$, which is a height $r_i \in \mathbb{R}$ from the ground. Likewise, this branch forms an angle $\beta_i \in \mathbb{R}$ with the ground XZ, while its projection in the XZ plane forms an angle $\theta_i \in \mathbb{R}$ with the X-axis, the direction of the applied force $f_x$. The kinetic energy of a branch can be computed as:

$$T_{B_i} = \frac{1}{2} \int_{0}^{l_{B_i}} \rho_{B_i} \mathbf{v}_{\lambda_i}^{\top} \mathbf{v}_{\lambda_i} d\lambda_i, \quad (14)$$

where $\rho_{B_i} \in \mathbb{R}$ is the linear density of the branch $i$, $\lambda_i \in [0, l_{B_i}]$ is the coordinate that parametrizes any point of the branch $i$ and $\mathbf{v}_{\lambda_i} \in \mathbb{R}^3$ is the velocity of the point $\lambda_i$ of the branch. The velocity $\mathbf{v}_{\lambda_i}$ can be expressed as a function of the linear velocity $\mathbf{v}_{R_i} \in \mathbb{R}^3$ and rotational speed of the trunk $\boldsymbol{\omega}_{R_i} \in \mathbb{R}^3$ at the branch-trunk union $R_i$:

$$\mathbf{v}_{\lambda_i} = \mathbf{v}_{R_i} + \boldsymbol{\omega}_{R_i} \times \mathbf{R}_i \boldsymbol{\lambda}_i, \quad (15)$$

where $\mathbf{R}_i \boldsymbol{\lambda}_i \in \mathbb{R}^3$ is the position vector from the trunk-branch union $R_i$ to the position $\boldsymbol{\lambda}_i$ in the branch-$i$. In addition, the linear velocity $\mathbf{v}_{R_i}$ and rotational speed $\boldsymbol{\omega}_{R_i}$ of the trunk can be expressed as a function of the degree of freedom $\dot{x}_a$:

$$\mathbf{v}_{R_i} = \dot{x}_a \eta(r_i) \mathbf{x}_W = \dot{x}_a \eta_{R_i} \mathbf{x}_W, \quad (16)$$

$$\boldsymbol{\omega}_{R_i} = \frac{d}{dt} \left( \frac{\partial x(\zeta, t)}{\partial \zeta} \right) \bigg|_{\zeta=r_i} \mathbf{z}_W = \dot{x}_a \eta'_{R_i} \mathbf{z}_W. \quad (17)$$

Finally, developing (14) with the velocity $\mathbf{v}_{\lambda_i}$ as a function of $\dot{x}_a$, we can identify the equivalent mass of the branch $B_i$. Doing this, we have:

$$m_{B_i}^{\text{eq}} = \rho_{B_i} \left[ \eta_{R_i}^2 l_{B_i} + \eta_{R_i}'^2 \frac{l_{B_i}^3}{3} \left[ \cos^2(\beta_i) \cos^2(\theta_i) + \sin^2(\beta_i) \right] - \eta_{R_i} \eta'_{R_i} l_{B_i}^2 \sin(\beta_i) \right], \quad (18)$$

where the dependence with the grasping point a appears in the expressions of $\eta_{R_i} = \eta(\zeta = r_i)$ and $\eta'_{R_i} = \frac{\partial \eta(\zeta)}{\partial \zeta} \bigg|_{\zeta=r_i}$.

3) Fruit parameters: the fruits are modeled as point masses attached to the branches. The kinetic energy of the j-th fruit $F_j$, can be expressed as:

$$T_{F_j} = \frac{1}{2} m_{F_j} \mathbf{v}_{F_j}^{\top} \mathbf{v}_{F_j}, \quad (19)$$

where $m_{F_j} \in \mathbb{R}$ is the mass of the fruit $F_j$ and $\mathbf{v}_{F_j} \in \mathbb{R}^3$ its velocity. The velocity $\mathbf{v}_{F_j}$ can be computed using Equation (15) with $\boldsymbol{\lambda}_i = l_{F_j}$, where $l_{F_j} \in \mathbb{R}$ is the distance between the fruit and the trunk-branch union point, as shown in Figure 3d. By doing this and taking out the common factor $\dot{x}_a$, the equivalent mass $m_{F_j}^{\text{eq}}$ of the fruit $F_j$ is obtained as:

$$m_{F_j}^{\text{eq}} = m_{F_j} \left[ \eta_{R_i}^2 - 2\eta_{R_i} \eta'_{R_i} l_{F_j} \sin(\beta_i) + l_{F_j}^2 \eta'_{R_i}^2 \left[ \sin^2(\beta_i) + \cos^2(\beta_i) \cos^2(\theta_i) \right] \right]. \quad (20)$$

<figure><img src="image_5.png" /><figcaption>Fig. 4. OMAV control scheme including the self-excited oscillation strategy for controlled shaking of flexible systems using the force $f_x$. We highlighted in red the novel shaking controller, which we integrate in a state-of-the-art Hybrid Force/Pose controller [11].</figcaption></figure>

### III. CONTROL AND ANALYSIS

We present here the different elements of the control pipeline that allow our aerial robot to control the shaking of a tree-like structure. For reference, the whole control system is visible in Figure 4.

#### A. OMAV controller

To induce vibration with the robot, we rely on a hybrid position/force controller that allows to apply forces in certain directions, while stabilizing the position along the others. In particular, we aim at controlling the tree oscillations in the forward direction (in our case the x direction of the $\mathcal{F}_B$ frame). As in [11], we define an hybrid force/pose controller as:

$$\boldsymbol{\tau}_c = \mathbf{\Lambda} \boldsymbol{\tau}_f + (\mathbf{I} - \mathbf{\Lambda}) \boldsymbol{\tau}_p + \mathbf{c}_{(B\omega)} + \mathbf{g}_{(W\mathbf{R}_B)}, \quad (21)$$

where $\mathbf{\Lambda} \in \mathbb{R}^{6 \times 6}$ is a diagonal and binary selection matrix with an element one or zero depending if we enable the force $\boldsymbol{\tau}_f$ or position $\boldsymbol{\tau}_p$ control action on the respective axis. The rest of the command compensates Coriolis and gravity wrenches.

1) Pose control: For the pose controller, we define the geometric pose errors as

$$\mathbf{e}_p = W \mathbf{R}_B^\top (W \mathbf{p} - W \mathbf{p}^{\text{ref}}), \quad (22a)$$

$$\mathbf{e}_R = \frac{1}{2} \text{vec} \left( W \mathbf{R}_B^{\text{ref},T} W \mathbf{R}_B - W \mathbf{R}_B^\top W \mathbf{R}_B^{\text{ref}} \right), \quad (22b)$$

with $W \mathbf{p}^{\text{ref}}$ and $W \mathbf{R}_B^{\text{ref}}$ the reference position and orientation, respectively, and $\text{vec}(\*)$ the operator to extract a vector from a skew-symmetric matrix. We adopt a PD control action as

$$\boldsymbol{\tau}_p = \mathbf{K}_p \tilde{\mathbf{e}} - \mathbf{K}_{DB} \mathbf{v}, \quad (23)$$

with $\mathbf{K}_D, \mathbf{K}_P \in \mathbb{R}^{6 \times 6}$ positive definite gain matrices and $\tilde{\mathbf{e}} = [\mathbf{e}_p^\top \mathbf{e}_R^\top]^\top$ the stacked pose error vector. Also, since we do not provide velocity references to the controller (only poses), the damping term multiplies only the current system's velocity, and not the velocity error.

2) Force control: For the force controller, we define the force error as

$$\mathbf{e}_f = \hat{\boldsymbol{\tau}}_{\text{ext}} - \boldsymbol{\tau}^{\text{ref}}, \quad (24)$$

where $\boldsymbol{\tau}^{\text{ref}} \in \mathbb{R}^6$ is the desired reference wrench. Then we define a PI controller as

$$\boldsymbol{\tau}_f = \mathbf{K}_F \mathbf{e}_f + \mathbf{K}_I \int \mathbf{e}_f dt - \boldsymbol{\tau}^{\text{ref}}, \quad (25)$$
$$k_T^{\text{eq}} = \frac{3EI}{a^3}, \qquad (13)$$

which show the dependence with the grasping point a.

2) Branch parameters: the branches are modeled as straight beams (see Assumption 3) embedded in the trunk of the tree, which do not deform, but are moved only by the movement of the trunk. As shown in Figures 3d and 3e, the branch $B_i$, with length $l_{B_i} \in \mathbb{R}$, is embedded in the trunk at the trunk point $R_i$, which is a height $r_i \in \mathbb{R}$ from the ground. Likewise, this branch forms an angle $\beta_i \in \mathbb{R}$ with the ground XZ, while its projection in the XZ plane forms an angle $\theta_i \in \mathbb{R}$ with the X-axis, the direction of the applied force $f_x$. The kinetic energy of a branch can be computed as:

$$T_{B_i} = \frac{1}{2} \int_{0}^{l_{B_i}} \rho_{B_i} \mathbf{v}_{\lambda_i}^{\top} \mathbf{v}_{\lambda_i} d\lambda_i, \qquad (14)$$

where $\rho_{B_i} \in \mathbb{R}$ is the linear density of the branch $i$, $\lambda_i \in [0, l_{B_i}]$ is the coordinate that parametrizes any point of the branch $i$ and $\mathbf{v}_{\lambda_i} \in \mathbb{R}^3$ is the velocity of the point $\lambda_i$ of the branch. The velocity $\mathbf{v}_{\lambda_i}$ can be expressed as a function of the linear velocity $\mathbf{v}_{R_i} \in \mathbb{R}^3$ and rotational speed of the trunk $\boldsymbol{\omega}_{R_i} \in \mathbb{R}^3$ at the branch-trunk union $R_i$:

$$\mathbf{v}_{\lambda_i} = \mathbf{v}_{R_i} + \boldsymbol{\omega}_{R_i} \times \mathbf{R}_i \boldsymbol{\lambda}_i, \qquad (15)$$

where $\mathbf{R}_i \boldsymbol{\lambda}_i \in \mathbb{R}^3$ is the position vector from the trunk-branch union $R_i$ to the position $\boldsymbol{\lambda}_i$ in the branch-i. In addition, the linear velocity $\mathbf{v}_{R_i}$ and rotational speed $\boldsymbol{\omega}_{R_i}$ of the trunk can be expressed as a function of the degree of freedom $\dot{x}_a$:

$$\mathbf{v}_{R_i} = \dot{x}_a \eta(r_i) \mathbf{x}_W = \dot{x}_a \eta_{R_i} \mathbf{x}_W, \qquad (16)$$

$$\boldsymbol{\omega}_{R_i} = \frac{d}{dt} \left( \frac{\partial x(\zeta, t)}{\partial \zeta} \right) \bigg|_{\zeta=r_i} \mathbf{z}_W = \dot{x}_a \eta'_{R_i} \mathbf{z}_W. \qquad (17)$$

Finally, developing (14) with the velocity $\mathbf{v}_{\lambda_i}$ as a function of $\dot{x}_a$, we can identify the equivalent mass of the branch $B_i$. Doing this, we have:

$$m_{B_i}^{\text{eq}} = \rho_{B_i} \left[ \eta_{R_i}^2 l_{B_i} + \eta_{R_i}^{\prime 2} \frac{l_{B_i}^3}{3} \left[ \cos^2(\beta_i) \cos^2(\theta_i) + \sin^2(\beta_i) \right] - \eta_{R_i} \eta'_{R_i} l_{B_i}^2 \sin(\beta_i) \right], \quad (18)$$

where the dependence with the grasping point a appears in the expressions of $\eta_{R_i} = \eta(\zeta = r_i)$ and $\eta'_{R_i} = \frac{\partial \eta(\zeta)}{\partial \zeta} \bigg|_{\zeta=r_i}$.

3) Fruit parameters: the fruits are modeled as point masses attached to the branches. The kinetic energy of the j-th fruit $F_j$, can be expressed as:

$$T_{F_j} = \frac{1}{2} m_{F_j} \mathbf{v}_{F_j}^{\top} \mathbf{v}_{F_j}, \qquad (19)$$

where $m_{F_j} \in \mathbb{R}$ is the mass of the fruit $F_j$ and $\mathbf{v}_{F_j} \in \mathbb{R}^3$ its velocity. The velocity $\mathbf{v}_{F_j}$ can be computed using Equation (15) with $\boldsymbol{\lambda}_i = l_{F_j}$, where $l_{F_j} \in \mathbb{R}$ is the distance between the fruit and the trunk-branch union point, as shown in Figure 3d. By doing this and taking out the common factor $\dot{x}_a$, the equivalent mass $m_{F_j}^{\text{eq}}$ of the fruit $F_j$ is obtained as:

$$m_{F_j}^{\text{eq}} = m_{F_j} \left[ \eta_{R_i}^2 - 2\eta_{R_i} \eta'_{R_i} l_{F_j} \sin(\beta_i) + l_{F_j}^2 \eta'_{R_i}^2 \left[ \sin^2(\beta_i) + \cos^2(\beta_i) \cos^2(\theta_i) \right] \right]. \quad (20)$$

<figure><img src="image_5.png" /><figcaption>Fig. 4. OMAV control scheme including the self-excited oscillation strategy for controlled shaking of flexible systems using the force $f_x$. We highlighted in red the novel shaking controller, which we integrate in a state-of-the-art Hybrid Force/Pose controller [11].</figcaption></figure>

### III. CONTROL AND ANALYSIS

We present here the different elements of the control pipeline that allow our aerial robot to control the shaking of a tree-like structure. For reference, the whole control system is visible in Figure 4.

#### A. OMAV controller

To induce vibration with the robot, we rely on a hybrid position/force controller that allows to apply forces in certain directions, while stabilizing the position along the others. In particular, we aim at controlling the tree oscillations in the forward direction (in our case the x direction of the $\mathcal{F}_B$ frame). As in [11], we define an hybrid force/pose controller as:

$$\boldsymbol{\tau}_c = \boldsymbol{\Lambda} \boldsymbol{\tau}_f + (\mathbf{I} - \boldsymbol{\Lambda}) \boldsymbol{\tau}_p + \mathbf{c}_{(B\omega)} + \mathbf{g}_{(W\mathbf{R}_B)}, \qquad (21)$$

where $\boldsymbol{\Lambda} \in \mathbb{R}^{6 \times 6}$ is a diagonal and binary selection matrix with an element one or zero depending if we enable the force $\boldsymbol{\tau}_f$ or position $\boldsymbol{\tau}_p$ control action on the respective axis. The rest of the command compensates Coriolis and gravity wrenches.

1) Pose control: For the pose controller, we define the geometric pose errors as

$$\mathbf{e}_p = W \mathbf{R}_B^{\top} (W \mathbf{p} - W \mathbf{p}^{\text{ref}}), \qquad (22a)$$

$$\mathbf{e}_R = \frac{1}{2} \text{vec} \left( W \mathbf{R}_B^{\text{ref},T} W \mathbf{R}_B - W \mathbf{R}_B^{\top} W \mathbf{R}_B^{\text{ref}} \right), \quad (22b)$$

with $W \mathbf{p}^{\text{ref}}$ and $W \mathbf{R}_B^{\text{ref}}$ the reference position and orientation, respectively, and $\text{vec}(\*)$ the operator to extract a vector from a skew-symmetric matrix. We adopt a PD control action as

$$\boldsymbol{\tau}_p = \mathbf{K}_p \tilde{\mathbf{e}} - \mathbf{K}_{DB} \mathbf{v}, \qquad (23)$$

with $\mathbf{K}_D, \mathbf{K}_P \in \mathbb{R}^{6 \times 6}$ positive definite gain matrices and $\tilde{\mathbf{e}} = [\mathbf{e}_p^{\top} \mathbf{e}_R^{\top}]^{\top}$ the stacked pose error vector. Also, since we do not provide velocity references to the controller (only poses), the damping term multiplies only the current system's velocity, and not the velocity error.

2) Force control: For the force controller, we define the force error as

$$\mathbf{e}_f = \hat{\boldsymbol{\tau}}_{\text{ext}} - \boldsymbol{\tau}^{\text{ref}}, \qquad (24)$$

where $\boldsymbol{\tau}^{\text{ref}} \in \mathbb{R}^6$ is the desired reference wrench. Then we define a PI controller as

$$\boldsymbol{\tau}_f = \mathbf{K}_F \mathbf{e}_f + \mathbf{K}_I \int \mathbf{e}_f dt - \boldsymbol{\tau}^{\text{ref}}, \qquad (25)$$
with $K_F, K_I \in \mathbb{R}^{6 \times 6}$ positive definite gain matrices. In practice, as we only need forces in the forward direction to oscillate the tree, the reference force vector will only have the first element non-zero $\tau^{\text{ref}} = [f_x^{\text{ref}} \mathbf{0}_5]^\top$, where $\mathbf{0}_5 \in \mathbb{R}^5$ is a vector of zeros of dimension 5.

## B. Self-Excited Oscillation Controller

For inducing vibrations into the tree, we use a self-excited oscillation strategy [30]. This feedback control force would generate self-excited oscillations at the natural frequency of the system [40]. Self-excited vibrations in a mechanical system are fundamentally attributed to the presence of an inherent state-dependent force, i.e., the forcing function is a function of the displacement, velocity, or acceleration of the system. In our case, this force is the UAV force applied to the tree, and its reference depends on the velocity of the contact point as follows:

$$f_x^{\text{ref}}(t) = k_c \text{sgn}(\dot{x}_a(t)), \quad (26)$$

where $k_c \in \mathbb{R}$ is the self-excited control gain and sgn is the signum function. As a result, a new dynamic equilibrium, called limit cycle, is attained, and the system tree-UAV system continues to vibrate autonomously. As shown in [40], using this control law with $k_c > 0$ into the system (3), the system exhibits stable vibrations with:

$$x_a(t) = R \cos(\omega_T t) \text{ with } R = \frac{4k_c}{\pi c^{\text{eq}} \omega_T}, \quad \omega_T = \sqrt{\frac{k^{\text{eq}}}{m^{\text{eq}}}}, \quad (27)$$

where $\omega_T \in \mathbb{R}$ is the oscillation frequency of the vibration, and it is equal to the natural frequency of the system (3). The amplitude of oscillation $R$ can be controlled by regulating the control gain $k_c$. Therefore, to induce vibrations in the tree at its natural frequency, the UAV will have the reference force presented in Equation (26). Finally, Figure 4 shows the complete control scheme of the OMAV, including the relay-feedback control force, where we consider that the vibration is induced in the forward direction of the aerial platform. This controlled shaking scheme does not require any additional information about what is being manipulated, only the velocity and interaction force, which are easily estimable or measurable. This allows it to be used with any flexible system, such as trees or branches. This control scheme combines the pose controller of (23) for tracking the position and attitude of the platform and the force controller of (25), whose reference is modified based on (26).

## C. Interaction Model

The self-excited vibration strategy presented before would induce oscillations into tree at its natural frequency, which can be obtained using Equation (27), as shown above. In addition, the presented strategy assumes that the force applied to the tree by the UAV $f_x(t)$ can instantaneously change its direction, as the function sgn is used. However, the UAV has non-negligible dynamics, which do not allow the force to change instantaneously. For this reason, the frequency of vibration might be different from the natural frequency of the tree and from the one estimated by (27).

To improve our estimation of the tree-UAV interaction, it is necessary to analyze the translation dynamics of the UAV in the X direction. According to the dynamics in (1), the translation dynamics in X are:

$$m_U \ddot{x}_U = f_x^{\text{Tree}} + f_x^{\text{prop}} \quad \rightarrow \quad f_x = f_x^{\text{prop}} - m_U \ddot{x}_U, \quad (28)$$

where $m_U \in \mathbb{R}$ is the UAV mass, $\ddot{x}_U$ is the X-acceleration of the UAV, $f_x^{\text{Tree}}$ the force applied by the tree into the UAV and $f_x^{\text{prop}}$ the control wrench generated by the propellers in the X-direction. Applying the action-reaction principle $f_x^{\text{Tree}} = -f_x$ and considering that the $\ddot{x}_U = \ddot{x}_a$ as the UAV graps the tree at the contact point a, the dynamic equation of the tree (3) with the UAV remains:

$$(m^{\text{eq}} + m_U) \ddot{x}_a + c^{\text{eq}} \dot{x}_a + k^{\text{eq}} x_a = f_x^{\text{prop}}. \quad (29)$$

For having a more precise interaction model, we assume that the desired interaction force $f_x^{\text{ref}}$ in (26) is the actuated force $f_x^{\text{prop}}$. This is a valid assumption since the actuated force $f_x^{\text{prop}}$ has a much faster dynamics than the translational one. Therefore, Equation (29) in combination with $f_x^{\text{prop}} = k_c \text{sgn}(\dot{x}_a(t))$, conforms our tree-UAV interaction model with the self-excited oscillation strategy. This results in the following motion for the contact point a:

$$x_a(t) = R \cos(\omega_{UT} t) \text{ with} \quad (30a)$$

$$R = \frac{4k_c}{\pi c^{\text{eq}} \omega_{UT}} \quad \text{and} \quad \omega_{UT} = \sqrt{\frac{k^{\text{eq}}}{m^{\text{eq}} + m_U}}, \quad (30b)$$

which now considers the translation dynamics of the UAV through the UAV mass $m_U$.

## D. Performance Analysis

Finally, this section analyzes the controlled shaking of the tree using the tree-UAV interaction model and the self-excited oscillation strategy. Two main elements are observed during the task executions: the tree oscillation frequency obtained and its amplitude at the tree's tip. In particular, their dependence with the grasping point a and with the self-excited control gain $k_c$ are evaluated.

1) Oscillation frequency: The oscillation frequency can be computed developing the expression (30b):

$$\omega_{UT} = \sqrt{\frac{k_T^{\text{eq}}}{m_T^{\text{eq}} + \sum_{i=1}^{N_B} m_{B_i}^{\text{eq}} + \sum_{j=1}^{N_F} m_{F_j}^{\text{eq}} + m_U}}, \quad (31)$$

where the equivalent parameters can be computed using the expressions presented in Section II-C. This expression shows that the oscillation frequency only depends on the grasping point a, and not the control gain $k_c$.

Analyzing the behavior of this expression with respect to the grasping point a, we can establish that as the grasping point moves to higher values, the oscillation frequency decreases. This can be explained analyzing the trunk stiffness depending on the grasping point a. As (13) shows, the trunk has a higher stiffness in the lower part. For this reason, the trunk can counteract the UAV force sooner than in the higher part of the trunk, as it needs a lower displacement. Hence, we obtain higher oscillation frequencies in the lower part of the trunk.
with $K_F, K_I \in \mathbb{R}^{6 \times 6}$ positive definite gain matrices. In practice, as we only need forces in the forward direction to oscillate the tree, the reference force vector will only have the first element non-zero $\tau^{\text{ref}} = [f_x^{\text{ref}} \mathbf{0}_5]^\top$, where $\mathbf{0}_5 \in \mathbb{R}^5$ is a vector of zeros of dimension 5.

## B. Self-Excited Oscillation Controller

For inducing vibrations into the tree, we use a self-excited oscillation strategy [30]. This feedback control force would generate self-excited oscillations at the natural frequency of the system [40]. Self-excited vibrations in a mechanical system are fundamentally attributed to the presence of an inherent state-dependent force, i.e., the forcing function is a function of the displacement, velocity, or acceleration of the system. In our case, this force is the UAV force applied to the tree, and its reference depends on the velocity of the contact point as follows:

$$f_x^{\text{ref}}(t) = k_c \text{sgn}(\dot{x}_a(t)), \quad (26)$$

where $k_c \in \mathbb{R}$ is the self-excited control gain and sgn is the signum function. As a result, a new dynamic equilibrium, called limit cycle, is attained, and the system tree-UAV system continues to vibrate autonomously. As shown in [40], using this control law with $k_c > 0$ into the system (3), the system exhibits stable vibrations with:

$$x_a(t) = R \cos(\omega_T t) \text{ with } R = \frac{4k_c}{\pi c^{\text{eq}} \omega_T}, \quad \omega_T = \sqrt{\frac{k^{\text{eq}}}{m^{\text{eq}}}}, \quad (27)$$

where $\omega_T \in \mathbb{R}$ is the oscillation frequency of the vibration, and it is equal to the natural frequency of the system (3). The amplitude of oscillation $R$ can be controlled by regulating the control gain $k_c$. Therefore, to induce vibrations in the tree at its natural frequency, the UAV will have the reference force presented in Equation (26). Finally, Figure 4 shows the complete control scheme of the OMAV, including the relay-feedback control force, where we consider that the vibration is induced in the forward direction of the aerial platform. This controlled shaking scheme does not require any additional information about what is being manipulated, only the velocity and interaction force, which are easily estimable or measurable. This allows it to be used with any flexible system, such as trees or branches. This control scheme combines the pose controller of (23) for tracking the position and attitude of the platform and the force controller of (25), whose reference is modified based on (26).

## C. Interaction Model

The self-excited vibration strategy presented before would induce oscillations into tree at its natural frequency, which can be obtained using Equation (27), as shown above. In addition, the presented strategy assumes that the force applied to the tree by the UAV $f_x(t)$ can instantaneously change its direction, as the function sgn is used. However, the UAV has non-negligible dynamics, which do not allow the force to change instantaneously. For this reason, the frequency of vibration might be different from the natural frequency of the tree and from the one estimated by (27).

To improve our estimation of the tree-UAV interaction, it is necessary to analyze the translation dynamics of the UAV in the X direction. According to the dynamics in (1), the translation dynamics in X are:

$$m_U \ddot{x}_U = f_x^{\text{Tree}} + f_x^{\text{prop}} \quad \rightarrow \quad f_x = f_x^{\text{prop}} - m_U \ddot{x}_U, \quad (28)$$

where $m_U \in \mathbb{R}$ is the UAV mass, $\ddot{x}_U$ is the X-acceleration of the UAV, $f_x^{\text{Tree}}$ the force applied by the tree into the UAV and $f_x^{\text{prop}}$ the control wrench generated by the propellers in the X-direction. Applying the action-reaction principle $f_x^{\text{Tree}} = -f_x$ and considering that the $\ddot{x}_U = \ddot{x}_a$ as the UAV graps the tree at the contact point $a$, the dynamic equation of the tree (3) with the UAV remains:

$$(m^{\text{eq}} + m_U)\ddot{x}_a + c^{\text{eq}}\dot{x}_a + k^{\text{eq}}x_a = f_x^{\text{prop}}. \quad (29)$$

For having a more precise interaction model, we assume that the desired interaction force $f_x^{\text{ref}}$ in (26) is the actuated force $f_x^{\text{prop}}$. This is a valid assumption since the actuated force $f_x^{\text{prop}}$ has a much faster dynamics than the translational one. Therefore, Equation (29) in combination with $f_x^{\text{prop}} = k_c \text{sgn}(\dot{x}_a(t))$, conforms our tree-UAV interaction model with the self-excited oscillation strategy. This results in the following motion for the contact point $a$:

$$x_a(t) = R \cos(\omega_{UT} t) \text{ with} \quad (30a)$$

$$R = \frac{4k_c}{\pi c^{\text{eq}} \omega_{UT}} \quad \text{and} \quad \omega_{UT} = \sqrt{\frac{k^{\text{eq}}}{m^{\text{eq}} + m_U}}, \quad (30b)$$

which now considers the translation dynamics of the UAV through the UAV mass $m_U$.

## D. Performance Analysis

Finally, this section analyzes the controlled shaking of the tree using the tree-UAV interaction model and the self-excited oscillation strategy. Two main elements are observed during the task executions: the tree oscillation frequency obtained and its amplitude at the tree's tip. In particular, their dependence with the grasping point $a$ and with the self-excited control gain $k_c$ are evaluated.

1) Oscillation frequency: The oscillation frequency can be computed developing the expression (30b):

$$\omega_{UT} = \sqrt{\frac{k_T^{\text{eq}}}{m_T^{\text{eq}} + \sum_{i=1}^{N_B} m_{B_i}^{\text{eq}} + \sum_{j=1}^{N_F} m_{F_j}^{\text{eq}} + m_U}}, \quad (31)$$

where the equivalent parameters can be computed using the expressions presented in Section II-C. This expression shows that the oscillation frequency only depends on the grasping point $a$, and not the control gain $k_c$.

Analyzing the behavior of this expression with respect to the grasping point $a$, we can establish that as the grasping point moves to higher values, the oscillation frequency decreases. This can be explained analyzing the trunk stiffness depending on the grasping point $a$. As (13) shows, the trunk has a higher stiffness in the lower part. For this reason, the trunk can counteract the UAV force sooner than in the higher part of the trunk, as it needs a lower displacement. Hence, we obtain higher oscillation frequencies in the lower part of the trunk.
<figure><img src="image_6.png" /><figcaption>Fig. 5. Bamboo tree built for experiments, consisting of a trunk and two branches. The different grasping points are marked with red circles in both pictures. Configuration 1: the UAV applies the force perpendicular to the plane of the tree. Configuration 2: the UAV applies the force in the same plane of the tree.</figcaption></figure>

2) Oscillation amplitude at the tree tip: As presented in Section III-C, the oscillation of the grasping point is given by (30). The movement of the tree tip $x_{Tip}(t)$ can be computed using (8) with $\zeta = l_T$. Doing that, $x_{Tip}(t) = x_a(t)\eta(l_T)$, where $\eta(l_T)$ can be computed using (11). The expression of $R$ can be computed using the equivalent damping $c^{eq}$ as a function of the damping ratio $\zeta^{eq} \in \mathbb{R}$, with $c^{eq} = 2m^{eq}\zeta^{eq}\omega_T$, where $m^{eq}$ is the equivalent mass of the tree and $\omega_T = \sqrt{k^{eq}/m^{eq}}$ is the natural frequency of the tree without considering the UAV $^4$. Using this relation, the amplitude of the oscillation at the tree tip $X_{Tip}$ is:

$$X_{Tip} = \frac{2\eta(l_T)k_c}{\pi m^{eq} \zeta^{eq} \omega_T \omega_{UT}} = \frac{2\eta(l_T)k_c}{\pi \zeta^{eq} k^{eq}} \sqrt{1 + \frac{m_U}{m^{eq}}}, \quad (32)$$

where the damping ratio $\zeta^{eq}$ can be identified experimentally.

Analysing this expression, we can establish that the oscillation amplitude $X_{Tip}$ increases linearly with the control gain $k_c$, as they have a linear dependence. In addition, the oscillation amplitude $X_{Tip}$ increases for higher grasping points $a$, as it has an inverse relation with $k^{eq}$. This curve trend can be explained analyzing the trunk stiffness. As Equation (13) establishes, the trunk stiffness decreases as the grasping point $a$ moves upward. Hence, for a given control gain $k_c$, the trunk will suffer higher oscillations in the higher part of the trunk than in the lower part.

Finally, using also the conclusions from the oscillation frequency, we can establish that if we are interested in having high frequency movements the UAV has to interact with the tree in the low part of the trunk, whereas if we are interested in achieving high amplitude movements the UAV has to interact with the tree in the high part of the trunk.

## IV. RESULTS

### A. Description of the experiments

The aerial platform presented in Section II-A has been used for the experiments. Also, in order to validate the tree model presented in Section II-C, a tree composed of bamboo canes

has been built. Notice that instead of using a real tree, we construct one ourselves, allowing for better control over its parameters and enabling us to properly validate our model. The tree is composed of a trunk and two branches, both made of bamboo, as shown in Figure 5. The tree also includes two apples as fruits at the branches. To make the union between the branches and the trunk, we used 3D printed pieces. We also use a motion capture system as an external positioning system to obtain ground truth positions of the tree and the aerial robot.

During the experiments, the UAV grasps the tree trunk at different points. In total, four different grasping points have been used, as shown in Figure 5. These points are $a_1 = 0.25$ m, $a_2 = 0.4$ m, $a_3 = 0.75$ m and $a_4 = 1.00$ m, which correspond to the distance $a$ from the base of the tree, as marked in Figure 5. In addition, different values of the gain $k_c$ (see (26)) are tested for each grasping point, to obtain different oscillation amplitudes. For each value of $k_c$, we wait until the steady state is reached, in which the system oscillates continuously. After capturing the data in the steady state, the next $k_c$ value is selected. For each experiment, we measured the the UAV's position and velocity, the position of the tip point of the trunk and the force applied to the tree by the UAV. In addition, we consider two configurations depending on the orientation of the UAV force with respect to the tree. In the Configuration 1, the force is perpendicular to the plane containing the tree, and in the Configuration 2, the force is coplanar to it. Figure 5 shows both configurations. This allows changing the dynamics of the system, by changing the equivalent mass of the branches and fruits, as shown in Section II-C. Each configuration is parameterized with the $\theta_i$ angle associated with each branch. For configuration 1 it was possible to perform the experiments at all the grasping points. However, for configuration 2, it was only possible to perform the experiments at $a_1$ and $a_4$, since at $a_2$ and $a_3$ the UAV does not have much clearance and hits the branch.

Figure 6 shows an example of the experimental data acquired during an experiment. In this case, the aerial robot interact with the tree thought the grasping point $a_4$ and in the configuration 1. This figure shows how the system

<footnote>$^4$We do not consider here the UAV as the damping ratio has been identified experimentally for the tree without UAV</footnote>
<figure><img src="image_6.png" /><figcaption>Fig. 5. Bamboo tree built for experiments, consisting of a trunk and two branches. The different grasping points are marked with red circles in both pictures. Configuration 1: the UAV applies the force perpendicular to the plane of the tree. Configuration 2: the UAV applies the force in the same plane of the tree.</figcaption></figure>

2) Oscillation amplitude at the tree tip: As presented in Section III-C, the oscillation of the grasping point is given by (30). The movement of the tree tip $x_{Tip}(t)$ can be computed using (8) with $\zeta = l_T$. Doing that, $x_{Tip}(t) = x_a(t)\eta(l_T)$, where $\eta(l_T)$ can be computed using (11). The expression of $R$ can be computed using the equivalent damping $c^{eq}$ as a function of the damping ratio $\zeta^{eq} \in \mathbb{R}$, with $c^{eq} = 2m^{eq}\zeta^{eq}\omega_T$, where $m^{eq}$ is the equivalent mass of the tree and $\omega_T = \sqrt{k^{eq}/m^{eq}}$ is the natural frequency of the tree without considering the UAV $^4$. Using this relation, the amplitude of the oscillation at the tree tip $X_{Tip}$ is:

$$X_{Tip} = \frac{2\eta(l_T)k_c}{\pi m^{eq} \zeta^{eq} \omega_T \omega_{UT}} = \frac{2\eta(l_T)k_c}{\pi \zeta^{eq} k^{eq}} \sqrt{1 + \frac{m_U}{m^{eq}}}, \quad (32)$$

where the damping ratio $\zeta^{eq}$ can be identified experimentally.

Analysing this expression, we can establish that the oscillation amplitude $X_{Tip}$ increases linearly with the control gain $k_c$, as they have a linear dependence. In addition, the oscillation amplitude $X_{Tip}$ increases for higher grasping points $a$, as it has an inverse relation with $k^{eq}$. This curve trend can be explained analyzing the trunk stiffness. As Equation (13) establishes, the trunk stiffness decreases as the grasping point $a$ moves upward. Hence, for a given control gain $k_c$, the trunk will suffer higher oscillations in the higher part of the trunk than in the lower part.

Finally, using also the conclusions from the oscillation frequency, we can establish that if we are interested in having high frequency movements the UAV has to interact with the tree in the low part of the trunk, whereas if we are interested in achieving high amplitude movements the UAV has to interact with the tree in the high part of the trunk.

## IV. RESULTS

### A. Description of the experiments

The aerial platform presented in Section II-A has been used for the experiments. Also, in order to validate the tree model presented in Section II-C, a tree composed of bamboo canes

has been built. Notice that instead of using a real tree, we construct one ourselves, allowing for better control over its parameters and enabling us to properly validate our model. The tree is composed of a trunk and two branches, both made of bamboo, as shown in Figure 5. The tree also includes two apples as fruits at the branches. To make the union between the branches and the trunk, we used 3D printed pieces. We also use a motion capture system as an external positioning system to obtain ground truth positions of the tree and the aerial robot.

During the experiments, the UAV grasps the tree trunk at different points. In total, four different grasping points have been used, as shown in Figure 5. These points are $a_1 = 0.25$ m, $a_2 = 0.4$ m, $a_3 = 0.75$ m and $a_4 = 1.00$ m, which correspond to the distance $a$ from the base of the tree, as marked in Figure 5. In addition, different values of the gain $k_c$ (see (26)) are tested for each grasping point, to obtain different oscillation amplitudes. For each value of $k_c$, we wait until the steady state is reached, in which the system oscillates continuously. After capturing the data in the steady state, the next $k_c$ value is selected. For each experiment, we measured the the UAV’s position and velocity, the position of the tip point of the trunk and the force applied to the tree by the UAV. In addition, we consider two configurations depending on the orientation of the UAV force with respect to the tree. In the Configuration 1, the force is perpendicular to the plane containing the tree, and in the Configuration 2, the force is coplanar to it. Figure 5 shows both configurations. This allows changing the dynamics of the system, by changing the equivalent mass of the branches and fruits, as shown in Section II-C. Each configuration is parameterized with the $\theta_i$ angle associated with each branch. For configuration 1 it was possible to perform the experiments at all the grasping points. However, for configuration 2, it was only possible to perform the experiments at $a_1$ and $a_4$, since at $a_2$ and $a_3$ the UAV does not have much clearance and hits the branch.

Figure 6 shows an example of the experimental data acquired during an experiment. In this case, the aerial robot interact with the tree thought the grasping point $a_4$ and in the configuration 1. This figure shows how the system

<footnote>$^4$We do not consider here the UAV as the damping ratio has been identified experimentally for the tree without UAV</footnote>
reaches, after a transient, a new equilibrium when changing from $k_c = 0.2$ N to $k_c = 0.3$ N. It also shows as the interaction force $f_x$ has a sinusoidal behaviour, while the wrench commanded to the propellers $f_x^{prop}$ changes faster due to the force reference $f_x^{ref}$ changes.

<figure><img src="image_7.png" /><figcaption>Fig. 6. Results of the experiment in the grasping point $a_4$ and with configuration 1, during a change of the control gain $k_c$. The upper figure represents the force reference $f_x^{ref}$, the estimated interaction force $f_x$ and the wrench commanded to the propellers in the interaction direction $f_x^{prop}$. The lower figure represents the displacement of the UAV $x_a$ and its velocity $\dot{x}_a$ in the interaction direction X.</figcaption></figure>

Furthermore, Figure 7 shows amplitude spectrum of $x_a(t)$ during the experiment at the rasping point $a_4$, in configuration 1 and with control gain $k_c = 0.2$ N. The spectrum shows that the system is oscillating with only one mode of vibration, as most of the amplitude is concentrated in a specific frequency ($\omega = 4.714 \text{ rad s}^{-1}$) which validates Assumption 1.

<figure><img src="image_8.png" /><figcaption>Fig. 7. Amplitude spectrum of $x_a(t)$ during the experiment at the rasping point $a_4$, in configuration 1 and with control gain $k_c = 0.2$ N.</figcaption></figure>

To precisely estimate the behavior of the closed-loop system with the proposed controller, we should know the true value of the parameters in the proposed model. Table I shows the experimental values of the parameters needed for the proposed model, where the damping ratio $\zeta^{eq}$ has been identified experimentally with the logarithmic decrement method [41]. In the case that no parameter is known, it can be estimated using any of the existing techniques: morphology identification of trees [36], [37], position estimation of apples [42], [43] or mechanical properties database [44], [45]. However, in this work we suppose that all the parameters are known and the parameter estimation is out of scope.

<table><caption>TABLE I
PARAMETERS NEEDED FOR THE MODEL</caption><thead><tr><th>Parameter</th><th>Value</th></tr></thead><tbody><tr><td>UAV mass, m<sub>U</sub></td><td>2.5 kg</td></tr><tr><td>Grasping points, a</td><td>[0.25, 0.4, 0.75, 1.00] m</td></tr><tr><td>Trunk length, l<sub>T</sub></td><td>2 m</td></tr><tr><td>Trunk density, ρ<sub>T</sub></td><td>0.13 kg m<sup>-1</sup></td></tr><tr><td>Branch densities, ρ<sub>B1</sub>, ρ<sub>B2</sub></td><td>ρ<sub>B1</sub> = ρ<sub>B2</sub> = 0.11 kg m<sup>-1</sup></td></tr><tr><td>Trunk flexural rigidity, EI</td><td>31.58 N m<sup>2</sup></td></tr><tr><td>Branch lengths, l<sub>B1</sub>, l<sub>B2</sub></td><td>l<sub>B1</sub> = l<sub>B2</sub> = 2 m</td></tr><tr><td>Branch-Trunk unions, r<sub>1</sub>, r<sub>2</sub></td><td>r<sub>1</sub> = r<sub>2</sub> = 0.5 m</td></tr><tr><td>Branch angles β<sub>i</sub>, β<sub>1</sub>, β<sub>2</sub></td><td>β<sub>1</sub> = 10°, β<sub>2</sub> = 170°</td></tr><tr><td>Branch angles in Conf. 1, θ<sub>1</sub><sup>C1</sup>, θ<sub>2</sub><sup>C1</sup></td><td>θ<sub>1</sub><sup>C1</sup> = 90°, θ<sub>2</sub><sup>C1</sup> = -90°</td></tr><tr><td>Branch angles in Conf. 2, θ<sub>1</sub><sup>C2</sup>, θ<sub>2</sub><sup>C2</sup></td><td>θ<sub>1</sub><sup>C2</sup> = 0°, θ<sub>2</sub><sup>C2</sup> = 180°</td></tr><tr><td>Fruit masses, m<sub>F1</sub>, m<sub>F2</sub></td><td>m<sub>F1</sub> = 0.13 kg, m<sub>F2</sub> = 0.12 kg</td></tr><tr><td>Fruit position, l<sub>F1</sub>, l<sub>F2</sub></td><td>l<sub>F1</sub> = 0.45 m, l<sub>F2</sub> = 1.3 m</td></tr><tr><td>Damping ratio, ζ<sup>eq</sup></td><td>0.098</td></tr></tbody></table>

## B. Experimental Validation

This section compares the experimental results with the predictions obtained using the proposed model and the self-excited oscillation strategy, presented in Section III-D.

1) **Oscillation frequency:** The oscillation frequencies for each control gain $k_c$ and for the different configurations and grasping points are represented in Figure 8 with markers of different colors. As we can see, the achieved oscillation frequency is independent of the control gain for a given grasping point and configuration, confirming the theoretical analysis of Section III-D. In addition, we represented in Figure 8 the mean value of the obtained frequencies. These mean values are compared with the theoretical values obtained with Equation (31), as shown in Figure 9. This curve shows that the proposed model properly estimates the oscillation frequency for a given grasping point and an orientation with respect to the tree, represented in the branch angles $\theta_i$. In addition, the curve trend shows as the oscillation frequency decreases when the UAV interacts with the tree at higher grasping points. This result validates the theoretical analysis done in Section III-D. Therefore, in order to achieve high-frequency oscillations, the UAV has to interact with the tree in the lower part of the trunk.

<table><caption>TABLE II
ESTIMATION ERROR FOR THE DIFFERENT CONFIGURATIONS.</caption><thead><tr><th>Mag.</th><th>Conf.</th><th>Error</th><th>a<sub>1</sub></th><th>a<sub>2</sub></th><th>a<sub>3</sub></th><th>a<sub>4</sub></th></tr></thead><tbody><tr><td rowspan="4">ω</td><td rowspan="2">1</td><td>e<sub>Abs</sub> [rad s<sup>-1</sup>]</td><td>0.66</td><td>1.15</td><td>0.86</td><td>0.49</td></tr><tr><td>e<sub>Rel</sub> [%]</td><td>3.61</td><td>8.02</td><td>12.87</td><td>9.94</td></tr><tr><td rowspan="2">2</td><td>e<sub>Abs</sub> [rad s<sup>-1</sup>]</td><td>1.10</td><td>-</td><td>-</td><td>0.22</td></tr><tr><td>e<sub>Rel</sub> [%]</td><td>8.26</td><td>-</td><td>-</td><td>5.01</td></tr><tr><td rowspan="4">X<sub>Tip</sub>/k<sub>c</sub></td><td rowspan="2">1</td><td>e<sub>Abs</sub> [mm N<sup>-1</sup>]</td><td>2.00</td><td>0.92</td><td>3.55</td><td>13.70</td></tr><tr><td>e<sub>Rel</sub> [%]</td><td>24.29</td><td>5.69</td><td>4.27</td><td>8.39</td></tr><tr><td rowspan="2">2</td><td>e<sub>Abs</sub> [mm N<sup>-1</sup>]</td><td>1.25</td><td>-</td><td>-</td><td>4.40</td></tr><tr><td>e<sub>Rel</sub> [%]</td><td>26.34</td><td>-</td><td>-</td><td>3.82</td></tr></tbody></table>

Finally, Table II summarizes the absolute estimation error.
reaches, after a transient, a new equilibrium when changing from $k_c = 0.2$ N to $k_c = 0.3$ N. It also shows as the interaction force $f_x$ has a sinusoidal behaviour, while the wrench commanded to the propellers $f_x^{prop}$ changes faster due to the force reference $f_x^{ref}$ changes.

<figure><img src="image_7.png" /><figcaption>Fig. 6. Results of the experiment in the grasping point $a_4$ and with configuration 1, during a change of the control gain $k_c$. The upper figure represents the force reference $f_x^{ref}$, the estimated interaction force $f_x$ and the wrench commanded to the propellers in the interaction direction $f_x^{prop}$. The lower figure represents the displacement of the UAV $x_a$ and its velocity $\dot{x}_a$ in the interaction direction X.</figcaption></figure>

Furthermore, Figure 7 shows amplitude spectrum of $x_a(t)$ during the experiment at the rasping point $a_4$, in configuration 1 and with control gain $k_c = 0.2$ N. The spectrum shows that the system is oscillating with only one mode of vibration, as most of the amplitude is concentrated in a specific frequency ($\omega = 4.714$ rad s$^{-1}$) which validates Assumption 1.

<figure><img src="image_8.png" /><figcaption>Fig. 7. Amplitude spectrum of $x_a(t)$ during the experiment at the rasping point $a_4$, in configuration 1 and with control gain $k_c = 0.2$ N.</figcaption></figure>

To precisely estimate the behavior of the closed-loop system with the proposed controller, we should know the true value of the parameters in the proposed model. Table I shows the experimental values of the parameters needed for the proposed model, where the damping ratio $\zeta^{eq}$ has been identified experimentally with the logarithmic decrement method [41]. In the case that no parameter is known, it can be estimated using any of the existing techniques: morphology identification of trees [36], [37], position estimation of apples [42], [43] or mechanical properties database [44], [45]. However, in this work we suppose that all the parameters are known and the parameter estimation is out of scope.

<table><caption>TABLE I
PARAMETERS NEEDED FOR THE MODEL</caption><thead><tr><th>Parameter</th><th>Value</th></tr></thead><tbody><tr><td>UAV mass, m<sub>U</sub></td><td>2.5 kg</td></tr><tr><td>Grasping points, a</td><td>[0.25, 0.4, 0.75, 1.00] m</td></tr><tr><td>Trunk length, l<sub>T</sub></td><td>2 m</td></tr><tr><td>Trunk density, ρ<sub>T</sub></td><td>0.13 kg m<sup>-1</sup></td></tr><tr><td>Branch densities, ρ<sub>B1</sub>, ρ<sub>B2</sub></td><td>ρ<sub>B1</sub> = ρ<sub>B2</sub> = 0.11 kg m<sup>-1</sup></td></tr><tr><td>Trunk flexural rigidity, EI</td><td>31.58 N m<sup>2</sup></td></tr><tr><td>Branch lengths, l<sub>B1</sub>, l<sub>B2</sub></td><td>l<sub>B1</sub> = l<sub>B2</sub> = 2 m</td></tr><tr><td>Branch-Trunk unions, r<sub>1</sub>, r<sub>2</sub></td><td>r<sub>1</sub> = r<sub>2</sub> = 0.5 m</td></tr><tr><td>Branch angles β<sub>i</sub>, β<sub>1</sub>, β<sub>2</sub></td><td>β<sub>1</sub> = 10°, β<sub>2</sub> = 170°</td></tr><tr><td>Branch angles in Conf. 1, θ<sub>1</sub><sup>C1</sup>, θ<sub>2</sub><sup>C1</sup></td><td>θ<sub>1</sub><sup>C1</sup> = 90°, θ<sub>2</sub><sup>C1</sup> = -90°</td></tr><tr><td>Branch angles in Conf. 2, θ<sub>1</sub><sup>C2</sup>, θ<sub>2</sub><sup>C2</sup></td><td>θ<sub>1</sub><sup>C2</sup> = 0°, θ<sub>2</sub><sup>C2</sup> = 180°</td></tr><tr><td>Fruit masses, m<sub>F1</sub>, m<sub>F2</sub></td><td>m<sub>F1</sub> = 0.13 kg, m<sub>F2</sub> = 0.12 kg</td></tr><tr><td>Fruit position, l<sub>F1</sub>, l<sub>F2</sub></td><td>l<sub>F1</sub> = 0.45 m, l<sub>F2</sub> = 1.3 m</td></tr><tr><td>Damping ratio, ζ<sup>eq</sup></td><td>0.098</td></tr></tbody></table>

## B. Experimental Validation

This section compares the experimental results with the predictions obtained using the proposed model and the self-excited oscillation strategy, presented in Section III-D.

1) **Oscillation frequency:** The oscillation frequencies for each control gain $k_c$ and for the different configurations and grasping points are represented in Figure 8 with markers of different colors. As we can see, the achieved oscillation frequency is independent of the control gain for a given grasping point and configuration, confirming the theoretical analysis of Section III-D. In addition, we represented in Figure 8 the mean value of the obtained frequencies. These mean values are compared with the theoretical values obtained with Equation (31), as shown in Figure 9. This curve shows that the proposed model properly estimates the oscillation frequency for a given grasping point and an orientation with respect to the tree, represented in the branch angles $\theta_i$. In addition, the curve trend shows as the oscillation frequency decreases when the UAV interacts with the tree at higher grasping points. This result validates the theoretical analysis done in Section III-D. Therefore, in order to achieve high-frequency oscillations, the UAV has to interact with the tree in the lower part of the trunk.

<table><caption>TABLE II
ESTIMATION ERROR FOR THE DIFFERENT CONFIGURATIONS.</caption><thead><tr><th>Mag.</th><th>Conf.</th><th>Error</th><th>a<sub>1</sub></th><th>a<sub>2</sub></th><th>a<sub>3</sub></th><th>a<sub>4</sub></th></tr></thead><tbody><tr><td rowspan="4">ω</td><td rowspan="2">1</td><td>e<sub>Abs</sub> [rad s<sup>-1</sup>]</td><td>0.66</td><td>1.15</td><td>0.86</td><td>0.49</td></tr><tr><td>e<sub>Rel</sub> [%]</td><td>3.61</td><td>8.02</td><td>12.87</td><td>9.94</td></tr><tr><td rowspan="2">2</td><td>e<sub>Abs</sub> [rad s<sup>-1</sup>]</td><td>1.10</td><td>-</td><td>-</td><td>0.22</td></tr><tr><td>e<sub>Rel</sub> [%]</td><td>8.26</td><td>-</td><td>-</td><td>5.01</td></tr><tr><td rowspan="4">X<sub>Tip</sub>/k<sub>c</sub></td><td rowspan="2">1</td><td>e<sub>Abs</sub> [mm N<sup>-1</sup>]</td><td>2.00</td><td>0.92</td><td>3.55</td><td>13.70</td></tr><tr><td>e<sub>Rel</sub> [%]</td><td>24.29</td><td>5.69</td><td>4.27</td><td>8.39</td></tr><tr><td rowspan="2">2</td><td>e<sub>Abs</sub> [mm N<sup>-1</sup>]</td><td>1.25</td><td>-</td><td>-</td><td>4.40</td></tr><tr><td>e<sub>Rel</sub> [%]</td><td>26.34</td><td>-</td><td>-</td><td>3.82</td></tr></tbody></table>

Finally, Table II summarizes the absolute estimation error.
<figure><img src="image_9.png" /><figcaption>Fig. 8. Oscillation frequency of the tree ω obtained with different control gains $k_c$, for the four grasping points a and both configurations. The markers represent the experimental results, while the dashed lines represent the mean value of the experimental data.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 10. Oscillation amplitude at the tree’s tip $X_{Tip}$ obtained with different control gains $k_c$, for the four grasping points a and both configurations. The markers represent the experimental results, while the dashed lines represent the interpolation line of the experimental data.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 9. Oscillation frequency ω for different grasping points a. The proposed model (see Equation (31)) are represented with a blue line for configuration 1 and with a red line for configuration 2, while the markers represent the mean value of the experimental data, represented in Figure 8</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 11. Ratio $X_{Tip}/k_c$ for different grasping points a. The proposed model (see Equation (32)) is represented in a blue line, while the markers represent the slope of the linear interpolation of the experimental data, represented in Figure 10.</figcaption></figure>

error $e_{Abs}$ and the relative estimation error $e_{Rel}$ obtained during the experiments for the different grasping points. For a given magnitude $f$ in a specific grasping point $a$, these errors are computed as $e_{Abs}(f) = |f_{Est} - \overline{f_{Exp}}|$ and $e_{Rel}(f) = 100|f_{Est} - \overline{f_{Exp}}|/f_{Exp}$, where $f_{Est}$ is the estimated value of $f$ for the grasping point $a$ and $\overline{f_{Exp}} = \sum_{i=0}^{N} f_{Exp,k_i}/N$ is the average value of $f$ for the different control gains $k_i \in \{k_1, k_2, ..., k_N\}$ in the experiments. Note that the estimation of $f$ does not depend on the control gain $k_c$, as it is the case with the oscillation frequency ω and the ratio $X_{Tip}/k_c$. This table shows that the proposed model can estimate the oscillation frequency properly, with errors as low as 3.61%. The higher errors in the estimation can be produced by parameter errors (damping, lengths, densities, etc.), by neglected higher vibration modes that could be also slightly excited during the experiments and/or the difference between the force reference $f_x^{ref}$ and the actuated one $f_x^{prop}$.

2) Oscillation amplitude at the tree tip: The oscillation amplitudes at the tree tip for each control gain $k_c$ and for the different configurations and grasping points are represented in

Figure 10 with markers of different colors. As we can see, in all cases the oscillation amplitude increases with the control gain, as the model predicts. In addition, our model establishes a linear relation between the control gain and the oscillation amplitude (see Equation (32)). Therefore we use a linear regressor to interpolate the data sets. The linear regressions are represented with dashed lines in Figure 10.

In addition, we compare the experimental results with the estimations of our model. In order to do that, we represent the linear regression slope obtained from the experimental results with the estimations of the oscillation amplitude for the different grasping points, as presented in Figure 11. This ratio $X_{Tip}/k_c$, expressed in units of m N⁻¹, can be understood as the mechanical admittance of the tree for a given grasping point $a$, as it represents the amplitude of deformation given a force (control gain). The curve trend shows that the oscillation amplitude grows not only with the control gain $k_c$, but also with the grasping point $a$, validating the theoretical analysis of Section III-D.

In the same way as before, Table II includes the absolute
<figure><img src="image_9.png" /><figcaption>Fig. 8. Oscillation frequency of the tree ω obtained with different control gains $k_c$, for the four grasping points a and both configurations. The markers represent the experimental results, while the dashed lines represent the mean value of the experimental data.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 10. Oscillation amplitude at the tree’s tip $X_{Tip}$ obtained with different control gains $k_c$, for the four grasping points a and both configurations. The markers represent the experimental results, while the dashed lines represent the interpolation line of the experimental data.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>Fig. 9. Oscillation frequency ω for different grasping points a. The proposed model (see Equation (31)) are represented with a blue line for configuration 1 and with a red line for configuration 2, while the markers represent the mean value of the experimental data, represented in Figure 8</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Fig. 11. Ratio $X_{Tip}/k_c$ for different grasping points a. The proposed model (see Equation (32)) is represented in a blue line, while the markers represent the slope of the linear interpolation of the experimental data, represented in Figure 10.</figcaption></figure>

error $e_{Abs}$ and the relative estimation error $e_{Rel}$ obtained during the experiments for the different grasping points. For a given magnitude $f$ in a specific grasping point $a$, these errors are computed as $e_{Abs}(f) = |f_{Est} - \overline{f_{Exp}}|$ and $e_{Rel}(f) = 100|f_{Est} - \overline{f_{Exp}}|/f_{Exp}$, where $f_{Est}$ is the estimated value of $f$ for the grasping point $a$ and $\overline{f_{Exp}} = \sum_{i=0}^{N} f_{Exp,k_i}/N$ is the average value of $f$ for the different control gains $k_i \in \{k_1, k_2, ..., k_N\}$ in the experiments. Note that the estimation of $f$ does not depend on the control gain $k_c$, as it is the case with the oscillation frequency $\omega$ and the ratio $X_{Tip}/k_c$. This table shows that the proposed model can estimate the oscillation frequency properly, with errors as low as 3.61%. The higher errors in the estimation can be produced by parameter errors (damping, lengths, densities, etc.), by neglected higher vibration modes that could be also slightly excited during the experiments and/or the difference between the force reference $f_x^{ref}$ and the actuated one $f_x^{prop}$.

2) Oscillation amplitude at the tree tip: The oscillation amplitudes at the tree tip for each control gain $k_c$ and for the different configurations and grasping points are represented in

Figure 10 with markers of different colors. As we can see, in all cases the oscillation amplitude increases with the control gain, as the model predicts. In addition, our model establishes a linear relation between the control gain and the oscillation amplitude (see Equation (32)). Therefore we use a linear regressor to interpolate the data sets. The linear regressions are represented with dashed lines in Figure 10.

In addition, we compare the experimental results with the estimations of our model. In order to do that, we represent the linear regression slope obtained from the experimental results with the estimations of the oscillation amplitude for the different grasping points, as presented in Figure 11. This ratio $X_{Tip}/k_c$, expressed in units of m N⁻¹, can be understood as the mechanical admittance of the tree for a given grasping point $a$, as it represents the amplitude of deformation given a force (control gain). The curve trend shows that the oscillation amplitude grows not only with the control gain $k_c$, but also with the grasping point $a$, validating the theoretical analysis of Section III-D.

In the same way as before, Table II includes the absolute
error and relative error obtained with the proposed model for each grasping point. The table shows that our model can predict the oscillation amplitude with errors as low as 3.82%. Similar to the frequency estimations, the higher error can be produced by parameter errors, higher vibration modes slightly excited or difference between $f_x^{ref}$ and $f_x^{prop}$.

Finally, the experimental results confirm the theoretical analysis performed in Section III-D, which established that high frequencies are achieved interacting with the lower part of the tree trunk while high oscillation amplitudes are obtaining interacting with the higher part of the tree.

## V. CONCLUSION

This paper presents the use of an aerial manipulator to controlled shake of trees, which can contribute in the future as a new method of fruit harvesting or the assessment of tree parameter for health monitoring. First, the paper has presented a 1 DoF dynamic model of a tree composed of a trunk, branches and fruits. The proposed model is based on the Rayleigh-Ritz method, which allows to make predictions on the tree-UAV behavior during the interaction. Then, this work has incorporated a self-induced oscillation strategy into the UAV controller, which guarantees the oscillation of the tree-UAV system at its natural frequency. Combining both, the tree model with the self-induced oscillation strategy, the resulting interaction model can estimate the oscillation frequency and the oscillation amplitude during the tree-UAV interaction and depending on the grasping point. This theoretical study helps in deriving useful information on how to interact with a tree according to the goal (high movements of the tree, high frequency movements, ...) of the task to be performed. The controlled shaking scheme and the proposed model have been evaluated with indoor flights, where the UAV controlled shakes a built tree of bamboo. The results show how the proposed model can properly estimate the oscillation frequency and oscillation amplitude of the tree tip with errors as low as 3.82% and 3.61% respectively.

Future works include the implementation of the proposed controlled shaking scheme in outdoors scenarios and with real trees. We will combine this with various estimation strategies to identify the parameters required for model estimations.

## REFERENCES

[1] A. Ollero, M. Tognon, A. Suarez, D. Lee, and A. Franchi, "Past, present, and future of aerial robotic manipulators," *IEEE Transactions on Robotics*, vol. 38, no. 1, pp. 626-645, 2021.

[2] N. Nikhil, S. Shreyas, G. Vyshnavi, and S. Yadav, "Unmanned aerial vehicles (UAV) in disaster management applications," in *ICSSIT*. IEEE, 2020, pp. 140-148.

[3] M. A. Trujillo, J. R. Martínez-de Dios, C. Martín, A. Viguria, and A. Ollero, "Novel Aerial Manipulator for Accurate and Robust Industrial NDT Contact Inspection: A New Tool for the Oil and Gas Inspection Industry," *Sensors*, vol. 19, no. 6, 2019.

[4] A. Ollero, G. Heredia, A. Franchi, G. Antonelli, K. Kondak, A. Sanfeliu, A. Viguria, J. R. Martinez-de Dios, F. Pierri, J. Cortes, A. Santamaría-Navarro, M. A. Trujillo Soto, R. Balachandran, J. Andrade-Cetto, and A. Rodriguez, "The AEROARMS Project: Aerial Robots with Advanced Manipulation Capabilities for Inspection and Maintenance," *IEEE Robotics Automation Magazine*, vol. 25, no. 4, pp. 12-23, 2018.

[5] A. Suarez, H. Romero, R. Salmoral, J. A. Acosta, J. Zambrano, and A. Ollero, "Experimental Evaluation of Aerial Manipulation Robot for the Installation of Clip Type Bird Diverters: Outdoor Flight Tests," in *AIRPHARO*, 2021, pp. 1-7.

[6] P. J. Sanchez-Cuevas, A. Gonzalez-Morgado, N. Cortes, D. B. Gayango, A. E. Jimenez-Cano, A. Ollero, and G. Heredia, "Fully-Actuated Aerial Manipulator for Infrastructure Contact Inspection: Design, Modeling, Localization, and Control," *Sensors*, vol. 20, no. 17, 2020.

[7] P. Lassen and M. Fumagalli, "Can your drone touch? exploring the boundaries of consumer-grade multirotors for physical interaction," in *2022 International Conference on Robotics and Automation (ICRA)*, 2022, pp. 1-7.

[8] M. Hamandi, F. Usai, Q. Sablé, N. Staub, M. Tognon, and A. Franchi, "Design of multitoro aerial vehicles: A taxonomy based on input allocation," *The International Journal of Robotics Research*, vol. 40, no. 8-9, pp. 1015-1044, 2021.

[9] R. Rashad, J. Goerres, R. Aarts, J. B. Engelen, and S. Stramigioli, "Fully actuated multirotor UAVs: A literature review," *IEEE RAM*, vol. 27, no. 3, pp. 97-107, 2020.

[10] C. Ding and L. Lu, "A tilting-rotor unmanned aerial vehicle for enhanced aerial locomotion and manipulation capabilities: Design, control, and applications," *IEEE/ASME Transactions on Mechatronics*, vol. 26, no. 4, pp. 2237-2248, 2021.

[11] K. Bodie, M. Brunner, M. Pantic, S. Walser, P. Pfndler, U. Angst, R. Siegwart, and J. Nieto, "An omnidirectional aerial manipulation platform for contact-based inspection." Robotics: Science and Systems Foundation, Jun 2019.

[12] R. Zufferey, J. Tormo-Barbero, M. M. Guzmán, F. J. Maldonado, E. Sanchez-Laulhe, P. Grau, M. Pérez, J. Á. Acosta, and A. Ollero, "Design of the high-payload flapping wing robot e-flap," *IEEE RA-L*, vol. 6, no. 2, pp. 3097-3104, 2021.

[13] H. Huang, W. He, J. Wang, L. Zhang, and Q. Fu, "An all servo-driven bird-like flapping-wing aerial robot capable of autonomous flight," *IEEE/ASME Transactions on Mechatronics*, vol. 27, no. 6, pp. 5484-5494, 2022.

[14] A. González-Morgado, C. Álvarez Cía, G. Heredia, and A. Ollero, "Fully-actuated, corner contact aerial robot for inspection of hard-to-reach bridge areas," in *2023 International Conference on Unmanned Aircraft Systems (ICUAS)*, 2023, pp. 1191-1198.

[15] M. Brunner, G. Rizzi, M. Studiger, R. Siegwart, and M. Tognon, "A planning-and-control framework for aerial manipulation of articulated objects," *IEEE Robotics and Automation Letters*, vol. 7, no. 4, pp. 10689-10696, 2022.

[16] F. Benzi, M. Brunner, M. Tognon, C. Secchi, and R. Siegwart, "Adaptive tank-based control for aerial physical interaction with uncertain dynamic environments using energy-task estimation," *IEEE Robotics and Automation Letters*, vol. 7, no. 4, pp. 9129-9136, 2022.

[17] M. Allenspach, Y. Vyas, M. Rubio, R. Siegwart, and M. Tognon, "Human-state-aware controller for a tethered aerial robot guiding a human by physical interaction," *IEEE Robotics and Automation Letters*, vol. 7, no. 2, pp. 2827-2834, 2022.

[18] E. Aucone, S. Kirchgeorg, A. Valentini, L. Pellissier, K. Deiner, and S. Mintchev, "Drone-assisted collection of environmental dna from tree branches for biodiversity monitoring," *Science Robotics*, vol. 8, no. 74, p. eadd5762, 2023.

[19] J. R. Kutia, K. A. Stol, and W. Xu, "Aerial manipulator interactions with trees for canopy sampling," *IEEE/ASME Transactions on Mechatronics*, vol. 23, no. 4, pp. 1740-1749, 2018.

[20] C. Geckeler, B. A. Pizzani, and S. Mintchev, "Biodegradable origami gripper actuated with gelatin hydrogel for aerial sensor attachment to tree branches," in *2023 IEEE International Conference on Robotics and Automation (ICRA)*, 2023, pp. 5324-5330.

[21] C. Geckeler, S. E. Ramos, M. C. Schuman, and S. Mintchev, "Robotic volatile sampling for early detection of plant stress: Precision agriculture beyond visual remote sensing," *IEEE Robotics & Automation Magazine*, pp. 2-12, 2023.

[22] S. Kirchgeorg and S. Mintchev, "Multimodal aerial-tethered robot for tree canopy exploration," in *2022 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)*, 2022, pp. 6080-6086.

[23] P. H. M. Souza and K. Stol, "Constrained dynamics of an aerial manipulator interacting with flexible cantilever beams," *IEEE/ASME Transactions on Mechatronics*, vol. 28, no. 2, pp. 967-975, 2023.

[24] Z. Zhang, C. Igathinathane, J. Li, H. Cen, Y. Lu, and P. Flores, "Technology progress in mechanical harvest of fresh market apples," *Computers and Electronics in Agriculture* vol. 175, p. 105606, 2020.

[25] C. Ortiz, A. Torregrosa, and S. Castro-García, "Comparison of a lightweight experimental shaker and an orchard tractor mounted trunk shaker for fresh market citrus harvesting," *Agriculture*, vol. 11, no. 11, 2021.

[26] J. Zhou, L. He, Q. Zhang, X. Du, D. Chen, and M. Karkee, "Evaluation of the influence of shaking frequency and duration in mechanical
error and relative error obtained with the proposed model for each grasping point. The table shows that our model can predict the oscillation amplitude with errors as low as 3.82%. Similar to the frequency estimations, the higher error can be produced by parameter errors, higher vibration modes slightly excited or difference between $f_x^{ref}$ and $f_x^{prop}$.

Finally, the experimental results confirm the theoretical analysis performed in Section III-D, which established that high frequencies are achieved interacting with the lower part of the tree trunk while high oscillation amplitudes are obtaining interacting with the higher part of the tree.

## V. CONCLUSION

This paper presents the use of an aerial manipulator to controlled shake of trees, which can contribute in the future as a new method of fruit harvesting or the assessment of tree parameter for health monitoring. First, the paper has presented a 1 DoF dynamic model of a tree composed of a trunk, branches and fruits. The proposed model is based on the Rayleigh-Ritz method, which allows to make predictions on the tree-UAV behavior during the interaction. Then, this work has incorporated a self-induced oscillation strategy into the UAV controller, which guarantees the oscillation of the tree-UAV system at its natural frequency. Combining both, the tree model with the self-induced oscillation strategy, the resulting interaction model can estimate the oscillation frequency and the oscillation amplitude during the tree-UAV interaction and depending on the grasping point. This theoretical study helps in deriving useful information on how to interact with a tree according to the goal (high movements of the tree, high frequency movements, ...) of the task to be performed. The controlled shaking scheme and the proposed model have been evaluated with indoor flights, where the UAV controlled shakes a built tree of bamboo. The results show how the proposed model can properly estimate the oscillation frequency and oscillation amplitude of the tree tip with errors as low as 3.82% and 3.61% respectively.

Future works include the implementation of the proposed controlled shaking scheme in outdoors scenarios and with real trees. We will combine this with various estimation strategies to identify the parameters required for model estimations.

## REFERENCES

[1] A. Ollero, M. Tognon, A. Suarez, D. Lee, and A. Franchi, "Past, present, and future of aerial robotic manipulators," *IEEE Transactions on Robotics*, vol. 38, no. 1, pp. 626-645, 2021.

[2] N. Nikhil, S. Shreyas, G. Vyshnavi, and S. Yadav, "Unmanned aerial vehicles (UAV) in disaster management applications," in *ICSSIT*. IEEE, 2020, pp. 140-148.

[3] M. A. Trujillo, J. R. Martínez-de Dios, C. Martín, A. Viguria, and A. Ollero, "Novel Aerial Manipulator for Accurate and Robust Industrial NDT Contact Inspection: A New Tool for the Oil and Gas Inspection Industry," *Sensors*, vol. 19, no. 6, 2019.

[4] A. Ollero, G. Heredia, A. Franchi, G. Antonelli, K. Kondak, A. Sanfeliu, A. Viguria, J. R. Martinez-de Dios, F. Pierri, J. Cortes, A. Santamaría-Navarro, M. A. Trujillo Soto, R. Balachandran, J. Andrade-Cetto, and A. Rodriguez, "The AEROARMS Project: Aerial Robots with Advanced Manipulation Capabilities for Inspection and Maintenance," *IEEE Robotics Automation Magazine*, vol. 25, no. 4, pp. 12-23, 2018.

[5] A. Suarez, H. Romero, R. Salmoral, J. A. Acosta, J. Zambrano, and A. Ollero, "Experimental Evaluation of Aerial Manipulation Robot for the Installation of Clip Type Bird Diverters: Outdoor Flight Tests," in *AIRPHARO*, 2021, pp. 1-7.

[6] P. J. Sanchez-Cuevas, A. Gonzalez-Morgado, N. Cortes, D. B. Gayango, A. E. Jimenez-Cano, A. Ollero, and G. Heredia, "Fully-Actuated Aerial Manipulator for Infrastructure Contact Inspection: Design, Modeling, Localization, and Control," *Sensors*, vol. 20, no. 17, 2020.

[7] P. Lassen and M. Fumagalli, "Can your drone touch? exploring the boundaries of consumer-grade multirotors for physical interaction," in *2022 International Conference on Robotics and Automation (ICRA)*, 2022, pp. 1-7.

[8] M. Hamandi, F. Usai, Q. Sablé, N. Staub, M. Tognon, and A. Franchi, "Design of multitoro aerial vehicles: A taxonomy based on input allocation," *The International Journal of Robotics Research*, vol. 40, no. 8-9, pp. 1015-1044, 2021.

[9] R. Rashad, J. Goerres, R. Aarts, J. B. Engelen, and S. Stramigioli, "Fully actuated multirotor UAVs: A literature review," *IEEE RAM*, vol. 27, no. 3, pp. 97-107, 2020.

[10] C. Ding and L. Lu, "A tilting-rotor unmanned aerial vehicle for enhanced aerial locomotion and manipulation capabilities: Design, control, and applications," *IEEE/ASME Transactions on Mechatronics*, vol. 26, no. 4, pp. 2237-2248, 2021.

[11] K. Bodie, M. Brunner, M. Pantic, S. Walser, P. Pfndler, U. Angst, R. Siegwart, and J. Nieto, "An omnidirectional aerial manipulation platform for contact-based inspection." Robotics: Science and Systems Foundation, Jun 2019.

[12] R. Zufferey, J. Tormo-Barbero, M. M. Guzmán, F. J. Maldonado, E. Sanchez-Laulhe, P. Grau, M. Pérez, J. Á. Acosta, and A. Ollero, "Design of the high-payload flapping wing robot e-flap," *IEEE RA-L*, vol. 6, no. 2, pp. 3097-3104, 2021.

[13] H. Huang, W. He, J. Wang, L. Zhang, and Q. Fu, "An all servo-driven bird-like flapping-wing aerial robot capable of autonomous flight," *IEEE/ASME Transactions on Mechatronics*, vol. 27, no. 6, pp. 5484-5494, 2022.

[14] A. González-Morgado, C. Álvarez Cía, G. Heredia, and A. Ollero, "Fully-actuated, corner contact aerial robot for inspection of hard-to-reach bridge areas," in *2023 International Conference on Unmanned Aircraft Systems (ICUAS)*, 2023, pp. 1191-1198.

[15] M. Brunner, G. Rizzi, M. Studiger, R. Siegwart, and M. Tognon, "A planning-and-control framework for aerial manipulation of articulated objects," *IEEE Robotics and Automation Letters*, vol. 7, no. 4, pp. 10689-10696, 2022.

[16] F. Benzi, M. Brunner, M. Tognon, C. Secchi, and R. Siegwart, "Adaptive tank-based control for aerial physical interaction with uncertain dynamic environments using energy-task estimation," *IEEE Robotics and Automation Letters*, vol. 7, no. 4, pp. 9129-9136, 2022.

[17] M. Allenspach, Y. Vyas, M. Rubio, R. Siegwart, and M. Tognon, "Human-state-aware controller for a tethered aerial robot guiding a human by physical interaction," *IEEE Robotics and Automation Letters*, vol. 7, no. 2, pp. 2827-2834, 2022.

[18] E. Aucone, S. Kirchgeorg, A. Valentini, L. Pellissier, K. Deiner, and S. Mintchev, "Drone-assisted collection of environmental dna from tree branches for biodiversity monitoring," *Science Robotics*, vol. 8, no. 74, p. eadd5762, 2023.

[19] J. R. Kutia, K. A. Stol, and W. Xu, "Aerial manipulator interactions with trees for canopy sampling," *IEEE/ASME Transactions on Mechatronics*, vol. 23, no. 4, pp. 1740-1749, 2018.

[20] C. Geckeler, B. A. Pizzani, and S. Mintchev, "Biodegradable origami gripper actuated with gelatin hydrogel for aerial sensor attachment to tree branches," in *2023 IEEE International Conference on Robotics and Automation (ICRA)*, 2023, pp. 5324-5330.

[21] C. Geckeler, S. E. Ramos, M. C. Schuman, and S. Mintchev, "Robotic volatile sampling for early detection of plant stress: Precision agriculture beyond visual remote sensing," *IEEE Robotics & Automation Magazine*, pp. 2-12, 2023.

[22] S. Kirchgeorg and S. Mintchev, "Multimodal aerial-tethered robot for tree canopy exploration," in *2022 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)*, 2022, pp. 6080-6086.

[23] P. H. M. Souza and K. Stol, "Constrained dynamics of an aerial manipulator interacting with flexible cantilever beams," *IEEE/ASME Transactions on Mechatronics*, vol. 28, no. 2, pp. 967-975, 2023.

[24] Z. Zhang, C. Igathinathane, J. Li, H. Cen, Y. Lu, and P. Flores, "Technology progress in mechanical harvest of fresh market apples," *Computers and Electronics in Agriculture* vol. 175, p. 105606, 2020.

[25] C. Ortiz, A. Torregrosa, and S. Castro-García, "Comparison of a lightweight experimental shaker and an orchard tractor mounted trunk shaker for fresh market citrus harvesting," *Agriculture*, vol. 11, no. 11, 2021.

[26] J. Zhou, L. He, Q. Zhang, X. Du, D. Chen, and M. Karkee, "Evaluation of the influence of shaking frequency and duration in mechanical
harvesting of sweet cherry,” *Applied Engineering in Agriculture*, vol. 29, no. 5, pp. 607–612, 2013.

[27] L. He, X. Liu, X. Du, and C. Wu, "In-situ identification of shaking frequency for adaptive vibratory fruit harvesting," *Computers and Electronics in Agriculture*, vol. 170, p. 105245, 2020.

[28] I. Kovacic, D. Radomirovic, and M. Zukovic, "Tree vibrations: Determining oscillatory properties by using infra-red marker-tracking system," *Urban Forestry & Urban Greening*, vol. 34, pp. 114–120, 2018.

[29] E. Grande, E. Giordano, and F. Clementi, "Evaluation of dynamic properties of trees subjected to induced vibrations," *Applied Sciences*, vol. 13, no. 12, 2023.

[30] W. Ding, *Self-Excited Vibration: Theory, Paradigms, and Research Methods*. Springer Berlin Heidelberg, 2010.

[31] K. R. James, G. A. Dahle, J. Grabosky, B. Kane, A. Detter et al., "Tree biomechanics literature review: Dynamics," *Arboric. Urban For.*, vol. 40, pp. 1–15, 2014.

[32] Y. Kumar, "The rayleigh-ritz method for linear dynamic, static and buckling behavior of beams, shells and plates: A literature review," *Journal of Vibration and Control*, vol. 24, no. 7, pp. 1205–1227, 2018.

[33] M. Sunar, "Theory of vibration — substructuring," in *Encyclopedia of Vibration*, S. Braun, Ed. Oxford: Elsevier, 2001, pp. 1332–1335.

[34] F. Ruggiero, J. Cacace, H. Sadeghian, and V. Lippiello, "Impedance control of vtol uavs with a momentum-based external generalized forces estimator," in *2014 IEEE international conference on robotics and automation (ICRA)*. IEEE, 2014, pp. 2093–2099.

[35] J. Moore and D. Maguire, "Natural sway frequencies and damping ratios of trees: Concepts, review and synthesis of previous studies," *Trees - Structure and Function*, vol. 18, pp. 195–203, 01 2004.

[36] C. Geckeler, E. Aucone, Y. Schnider, A. Simeon, J.-P. von Bassewitz, Y. Zhu, and S. Mintchev, "Learning occluded branch depth maps in forest environments using rgb-d images," *IEEE Robotics and Automation Letters*, vol. 9, no. 3, pp. 2439–2446, 2024.

[37] C. H. Kim and G. Kantor, "Occlusion reasoning for skeleton extraction of self-occluded tree canopies," in *2023 IEEE International Conference on Robotics and Automation (ICRA)*, 2023, pp. 9580–9586.

[38] S. Mazumder, "Chapter 2 - the finite difference method," in *Numerical Methods for Partial Differential Equations*, S. Mazumder, Ed. Academic Press, 2016, pp. 51–101.

[39] A. M. Kabe and B. H. Sako, "Chapter 6 - multi-degree-of-freedom systems," in *Structural Dynamics Fundamentals and Advanced Applications*, A. M. Kabe and B. H. Sako, Eds. Academic Press, 2020, pp. 333–435.

[40] A. Malas, and S. Chatterjee, "Generating self-excited oscillation in a class of mechanical systems by relay-feedback," *Nonlinear Dynamics*, vol. 76, no. 2, p. 1253–1269, 2014.

[41] L. B. Magalas and T. Malinowski, "Measurement techniques of the logarithmic decrement," *Solid State Phenomena*, vol. 89, pp. 247–260, 2003.

[42] A. Gongal, S. Amatya, M. Karkee, Q. Zhang, and K. Lewis, "Sensors and systems for fruit detection and localization: A review," *Computers and Electronics in Agriculture*, vol. 116, pp. 8–19, 2015.

[43] A. Silwal, M. Karkee, and Q. Zhang, "A hierarchical approach to apple identification for robotic harvesting," *Transactions of the ASABE*, vol. 59, no. 5, pp. 1079–1086, September 2016.

[44] B. Bartolucci, A. Rosa, C. Bertolin, F. Berto, F. Penta, A. Siani, C. Bartolucci, D. Rosa, A. Bertolini, C. Berto, and F. Siani, "Mechanical properties of the most common european woods: a literature review," *Frattura ed Integrità Strutturale*, vol. 14, pp. 249–274, 08 2020.

[45] F. Arriaga, X. Wang, G. Íñiguez González, D. F. Llana, M. Esteban, and P. Niemz, "Mechanical properties of wood: A review," *Forests*, vol. 14, no. 6, 2023.

<figure><img src="image_13.png" /></figure>

**Antonio González-Morgado** received B.Eng. and M.Eng. degrees in Industrial Engineering from University of Seville in 2018 and 2021, respectively. He has also received the Diplôme d'Ingénieur Généraliste and the Master of Aerospace Engineering in 2021 from École Centrale Lyon. He has been intern in Velan, France and CEA, France, in 2020 and 2021, respectively. Currently, he is a Ph.D. student supported the Marie Curie AERO-TRAIN project. He has been visiting researcher at the Autonomous System Lab, ETH Zurich, under

the supervision of Marco Tognon and Roland Siegwart. His main research interests are in design and control of aerial manipulators.

<figure><img src="image_14.png" /></figure>

**Eugenio Cuniato** (Graduate Student Member, IEEE) is currently pursuing a PhD degree at the Autonomous Systems Lab, ETH Zurich, Switzerland with a Marie Skłodowska-Curie (MSCA) scholarship. He is part of the AERO-TRAIN European Project whose goal is to develop innovative aerial manipulation technologies for industrial inspection and maintenance. He received his Master's Degree in Automation Engineering from the University of Naples Federico II (2021). His research focus is on control of aerial manipulators for high-performance aerial physical interaction. He received the Best Paper Award on Safety, Security and Rescue Robotics at the IEEE IROS 2022 conference.

<figure><img src="image_15.png" /></figure>

**Marco Tognon** (Member, IEEE) received the M.Sc. degree in automation engineering from the University of Padua, Italy, in 2014, with a master thesis carried out at the Max Planck Institute for Biological Cybernetics, Germany, and the Ph.D. degree in robotics from the INSA Toulouse, France, in 2018, developing his thesis at LAAS-CNRS, Toulouse. Since November 2022, he has been a Tenured Researcher with the Rainbow Team, INRIA Rennes, Rennes, France. Before that, he was a Postdoctoral Researcher with the Autonomous System Lab, ETH Zürich, Switzerland. From 2018 to 2020, he was a Postdoctoral Researcher with LAAS-CNRS. His current research interests include robotics, aerial physical interaction, multirobot systems, and human-robot physical interaction.

<figure><img src="image_16.png" /></figure>

**Guillermo Heredia** (Member, IEEE) is Full Professor at University of Seville (Spain). He was a Visiting Researcher at the Field Robotics Centre, Carnegie Mellon University (USA), and worked for an international automobile manufacturer (General Motors). He participated as Senior Researcher in 65 R&D projects, leading or co-leading 12 of them, including the FP7 ECSAFEMOBIL and the H2020 AEROBI and HYFLIERS projects. He is author or co-author of more than 100 publications on aerial robotics and aerial manipulation.

<figure><img src="image_17.png" /></figure>

**Roland Siegwart** is a Professor of autonomous mobile robots with ETH Zürich; the Founding Co-Director of the Technology Transfer Center, Wyss Zurich; and a board member of multiple high-tech companies. From 1996 to 2006, he was a Professor with EPFL Lausanne. He has held a visiting position with Stanford University and NASA Ames. He was the Vice President of ETH Zürich, from 2010 to 2014. His research interests include the design, control, and navigation of flying and wheeled and walking robots operating in complex and highly dynamical environments. He received the IEEE RAS Pioneer Award and the IEEE RAS Inaba Technical Award. He is among the most cited scientist in robots world-wide, a co-founder of more than half a dozen spin-off companies, and a strong promoter of innovation and entrepreneurship in Switzerland.

<figure><img src="image_18.png" /></figure>

**Aníbal Ollero** (Fellow, IEEE) is currently a Full Professor and the Head of the GRVC Robotics Laboratory, and the Scientific Advisor of the Center for Aerospace Technologies (CATEC), Seville. He has been a Full Professor at the University of Santiago and the University of Malaga, Spain, and a Researcher with the Robotics Institute, Carnegie Mellon University, Pittsburgh, USA, and LAAS-CNRS, Toulouse, France. He has authored more than 800 publications, including 9 books and more than 235 articles in journals, and has been an editor of 15 books. He led more than 160 research projects, participating in more than 41 projects of the European research programs, being coordinator of 7. He has been recognized with 28 awards, including the Spanish National Research Award in Engineering, the Rei Jaume I in New Technologies, the Overall Information and Communication Technologies Innovation Radar Prize 2017 of the European Commission, and several best paper awards in conferences.
harvesting of sweet cherry,” *Applied Engineering in Agriculture*, vol. 29, no. 5, pp. 607–612, 2013.

[27] L. He, X. Liu, X. Du, and C. Wu, "In-situ identification of shaking frequency for adaptive vibratory fruit harvesting," *Computers and Electronics in Agriculture*, vol. 170, p. 105245, 2020.

[28] I. Kovacic, D. Radomirovic, and M. Zukovic, "Tree vibrations: Determining oscillatory properties by using infra-red marker-tracking system," *Urban Forestry & Urban Greening*, vol. 34, pp. 114–120, 2018.

[29] E. Grande, E. Giordano, and F. Clementi, "Evaluation of dynamic properties of trees subjected to induced vibrations," *Applied Sciences*, vol. 13, no. 12, 2023.

[30] W. Ding, *Self-Excited Vibration: Theory, Paradigms, and Research Methods*. Springer Berlin Heidelberg, 2010.

[31] K. R. James, G. A. Dahle, J. Grabosky, B. Kane, A. Detter et al., "Tree biomechanics literature review: Dynamics," *Arboric. Urban For.*, vol. 40, pp. 1–15, 2014.

[32] Y. Kumar, "The rayleigh-ritz method for linear dynamic, static and buckling behavior of beams, shells and plates: A literature review," *Journal of Vibration and Control*, vol. 24, no. 7, pp. 1205–1227, 2018.

[33] M. Sunar, "Theory of vibration — substructuring," in *Encyclopedia of Vibration*, S. Braun, Ed. Oxford: Elsevier, 2001, pp. 1332–1335.

[34] F. Ruggiero, J. Cacace, H. Sadeghian, and V. Lippiello, "Impedance control of vtol uavs with a momentum-based external generalized forces estimator," in *2014 IEEE international conference on robotics and automation (ICRA)*. IEEE, 2014, pp. 2093–2099.

[35] J. Moore and D. Maguire, "Natural sway frequencies and damping ratios of trees: Concepts, review and synthesis of previous studies," *Trees - Structure and Function*, vol. 18, pp. 195–203, 01 2004.

[36] C. Geckeler, E. Aucone, Y. Schnider, A. Simeon, J.-P. von Bassewitz, Y. Zhu, and S. Mintchev, "Learning occluded branch depth maps in forest environments using rgb-d images," *IEEE Robotics and Automation Letters*, vol. 9, no. 3, pp. 2439–2446, 2024.

[37] C. H. Kim and G. Kantor, "Occlusion reasoning for skeleton extraction of self-occluded tree canopies," in *2023 IEEE International Conference on Robotics and Automation (ICRA)*, 2023, pp. 9580–9586.

[38] S. Mazumder, "Chapter 2 - the finite difference method," in *Numerical Methods for Partial Differential Equations*, S. Mazumder, Ed. Academic Press, 2016, pp. 51–101.

[39] A. M. Kabe and B. H. Sako, "Chapter 6 - multi-degree-of-freedom systems," in *Structural Dynamics Fundamentals and Advanced Applications*, A. M. Kabe and B. H. Sako, Eds. Academic Press, 2020, pp. 333–435.

[40] A. Malas, and S. Chatterjee, "Generating self-excited oscillation in a class of mechanical systems by relay-feedback," *Nonlinear Dynamics*, vol. 76, no. 2, p. 1253–1269, 2014.

[41] L. B. Magalas and T. Malinowski, "Measurement techniques of the logarithmic decrement," *Solid State Phenomena*, vol. 89, pp. 247–260, 2003.

[42] A. Gongal, S. Amatya, M. Karkee, Q. Zhang, and K. Lewis, "Sensors and systems for fruit detection and localization: A review," *Computers and Electronics in Agriculture*, vol. 116, pp. 8–19, 2015.

[43] A. Silwal, M. Karkee, and Q. Zhang, "A hierarchical approach to apple identification for robotic harvesting," *Transactions of the ASABE*, vol. 59, no. 5, pp. 1079–1086, September 2016.

[44] B. Bartolucci, A. Rosa, C. Bertolin, F. Berto, F. Penta, A. Siani, C. Bartolucci, D. Rosa, A. Bertolini, C. Berto, and F. Siani, "Mechanical properties of the most common european woods: a literature review," *Frattura ed Integrità Strutturale*, vol. 14, pp. 249–274, 08 2020.

[45] F. Arriaga, X. Wang, G. Íñiguez González, D. F. Llana, M. Esteban, and P. Niemz, "Mechanical properties of wood: A review," *Forests*, vol. 14, no. 6, 2023.

<figure><img src="image_13.png" /></figure>

**Eugenio Cuniato** (Graduate Student Member, IEEE) is currently pursuing a PhD degree at the Autonomous Systems Lab, ETH Zurich, Switzerland with a Marie Skłodowska-Curie (MSCA) scholarship. He is part of the AERO-TRAIN European Project whose goal is to develop innovative aerial manipulation technologies for industrial inspection and maintenance. He received his Master's Degree in Automation Engineering from the University of Naples Federico II (2021). His research focus is on control of aerial manipulators for high-performance aerial physical interaction. He received the Best Paper Award on Safety, Security and Rescue Robotics at the IEEE IROS 2022 conference.

<figure><img src="image_14.png" /></figure>

**Marco Tognon** (Member, IEEE) received the M.Sc. degree in automation engineering from the University of Padua, Italy, in 2014, with a master thesis carried out at the Max Planck Institute for Biological Cybernetics, Germany, and the Ph.D. degree in robotics from the INSA Toulouse, France, in 2018, developing his thesis at LAAS-CNRS, Toulouse. Since November 2022, he has been a Tenured Researcher with the Rainbow Team, INRIA Rennes, Rennes, France. Before that, he was a Postdoctoral Researcher with the Autonomous System Lab, ETH Zürich, Switzerland. From 2018 to 2020, he was a Postdoctoral Researcher with LAAS-CNRS. His current research interests include robotics, aerial physical interaction, multirobot systems, and human-robot physical interaction.

<figure><img src="image_15.png" /></figure>

**Guillermo Heredia** (Member, IEEE) is Full Professor at University of Seville (Spain). He was a Visiting Researcher at the Field Robotics Centre, Carnegie Mellon University (USA), and worked for an international automobile manufacturer (General Motors). He participated as Senior Researcher in 65 R&D projects, leading or co-leading 12 of them, including the FP7 ECSAFEMOBIL and the H2020 AEROBI and HYFLIERS projects. He is author or co-author of more than 100 publications on aerial robotics and aerial manipulation.

<figure><img src="image_16.png" /></figure>

**Roland Siegwart** is a Professor of autonomous mobile robots with ETH Zürich; the Founding Co-Director of the Technology Transfer Center, Wyss Zurich; and a board member of multiple high-tech companies. From 1996 to 2006, he was a Professor with EPFL Lausanne. He has held a visiting position with Stanford University and NASA Ames. He was the Vice President of ETH Zürich, from 2010 to 2014. His research interests include the design, control, and navigation of flying and wheeled and walking robots operating in complex and highly dynamical environments. He received the IEEE RAS Pioneer Award and the IEEE RAS Inaba Technical Award. He is among the most cited scientist in robots world-wide, a co-founder of more than half a dozen spin-off companies, and a strong promoter of innovation and entrepreneurship in Switzerland.

<figure><img src="image_17.png" /></figure>

**Aníbal Ollero** (Fellow, IEEE) is currently a Full Professor and the Head of the GRVC Robotics Laboratory, and the Scientific Advisor of the Center for Aerospace Technologies (CATEC), Seville. He has been a Full Professor at the University of Santiago and the University of Malaga, Spain, and a Researcher with the Robotics Institute, Carnegie Mellon University, Pittsburgh, USA, and LAAS-CNRS, Toulouse, France. He has authored more than 800 publications, including 9 books and more than 235 articles in journals, and has been an editor of 15 books. He led more than 160 research projects, participating in more than 41 projects of the European research programs, being coordinator of 7. He has been recognized with 28 awards, including the Spanish National Research Award in Engineering, the Rei Jaume I in New Technologies, the Overall Information and Communication Technologies Innovation Radar Prize 2017 of the European Commission, and several best paper awards in conferences.

<figure><img src="image_18.png" /></figure>

**Antonio González-Morgado** received B.Eng. and M.Eng. degrees in Industrial Engineering from University of Seville in 2018 and 2021, respectively. He has also received the Diplôme d'Ingénieur Généraliste and the Master of Aerospace Engineering in 2021 from École Centrale Lyon. He has been intern in Velan, France and CEA, France, in 2020 and 2021, respectively. Currently, he is a Ph.D. student supported the Marie Curie AERO-TRAIN project. He has been visiting researcher at the Autonomous System Lab, ETH Zurich, under

the supervision of Marco Tognon and Roland Siegwart. His main research interests are in design and control of aerial manipulators.