Energy Efficient Path Planning for 3D Aerial Inspections

The use of Unmanned Aerial Vehicles (UAVs) in different inspection tasks is increasing. This technology reduces inspection costs and collects high quality data of distinct structures, including areas that are not easily accessible by human operators. However, the reduced energy available on the UAVs limits their flight endurance. To increase the autonomy of a single flight, it is important to optimize the path to be performed by the UAV, in terms of energy loss. Therefore, this work presents a novel formulation of the Travelling Salesman Problem (TSP) and a path planning algorithm that uses a UAV energy model to solve this optimization problem. The novel TSP formulation is defined as Asymmetric Travelling Salesman Problem with Precedence Loss (ATSP-PL), where the cost of moving the UAV depends on the previous position. The energy model relates each UAV movement with its energy consumption, while the path planning algorithm is focused on minimizing the energy loss of the UAV, ensuring that the structure is fully covered. The developed algorithm was tested in both simulated and real scenarios. The simulated experiments were performed with realistic models of wind turbines and a UAV, whereas the real experiments were performed with a real UAV and an illumination tower. The inspection paths generated presented improvements over 24% and 8%, when compared with other methods, for the simulated and real experiments, respectively, optimizing the energy consumption of the UAV.


I. INTRODUCTION
Traditional methods for inspecting and monitoring diverse types of structures are often expensive, repetitive, time consuming and potentially dangerous for human operators [1]. To overcome this issue, several inspection tasks are currently being carried out by Unmanned Aerial Vehicles (UAVs) [2], [3], which represents a reduction of the inspection costs and of the total time consumed [4]. UAVs are autonomous aerial robots that can collect high quality data of structures [5], such as images (visual, thermal, among others) and point clouds. These vehicles can also easily access the entire sur-The associate editor coordinating the review of this manuscript and approving it for publication was Mohammad Alshabi . face of inspection. However, one of the challenges related with the use of UAVs in autonomous inspections is related with their limited battery capacity [5]. To allow a UAV to perform a wider inspection in a single flight, it is important to increase its endurance by selecting trajectories that optimize the energy consumption during an inspection procedure. A path that minimizes the total displacement performed by a UAV does not guarantee an optimized energy consumption. Therefore, a preferable approach consists in selecting a path that directly minimizes the energy consumption of a UAV (that may not be the shortest path). This can be achieved by smoothing the trajectory of the UAV, by reducing its accelerations during the performance of the path Therefore, this work presents an algorithm for 3D path planning for inspection of structures with UAVs, minimizing the energy loss. As a result, the contributions of this work are the following: 1) A novel Travelling Salesman Problem (TSP) formulation, named ATSP-PL (Asymmetric Travelling Salesman Problem, with Precedence Loss), of which the Precedence Loss designation is related to the objective function, where the loss of moving from a current position to the next one depends on the previous one. 2) An energy model of a UAV that relates each movement that the vehicle performs with its correspondent energy loss, and a resulting algorithm for inspection path planning that optimizes the energy consumption of a UAV, enabling it to autonomously perform a full coverage of the structure being inspected; 3) Validation of the developed algorithm in both simulated and real scenarios, using realistic models of wind turbines and of a UAV in the first case, and a real UAV along with an illumination tower in the second one. This article is organized as follows: Section II presents a brief review of the related work. Afterwards, section III presents the algorithm developed to obtain the inspection path and the respective UAV energy model. Experimental evaluation is demonstrated in sections IV and V, which includes experiments in both simulation scenarios and in a real environment. Finally, section VI clarifies the most important conclusions of this research.

II. RELATED WORK
The use of UAVs in different types of inspection tasks is increasing significantly. The most common sensors used for inspections with UAVs are both visual and thermal cameras [6], [7]. These sensors allow to execute non-destructive tests, capable of detecting the most significant in-service defects [7]. To improve the quality of the data collected by a camera sensor during a inspection, a LiDAR sensor can also be used [8], [9]. A 3D LiDAR represents a relevant method to obtain the characteristics of the inspected structure, given its capability of directly obtaining depth, resulting in accurate 3D point clouds of the surrounding environment [10]. These 3D sensors are used in robotics to improve the accuracy in distinct tasks, such as in object detection [11], [12], scene reconstruction and navigation [13], [14] Although a UAV equipped with the required sensors is able to collect the necessary information for an inspection task, its limited flight endurance represents a considerable drawback. To increase the autonomy of the UAVs, it is important to optimize the path that the vehicles will perform as a function of its energy consumption, while assuring that the entire surface of inspection is completely verified. The energy model of a UAV represents the energy loss relative to each individual movement of the vehicle, such as the cost of displacement and the cost of directional changes [15]. Each UAV has its own particular characteristics, such as weight, dimensions, propellers and motors type. For that reason, it is not easy to derive a general parametric model that represents the energy consumption of each different UAV. Therefore, one possible way consists in performing tests in each particular UAV. Franco and Buttazzo [16], [17] performed a complete analysis on the energy loss based on real measurements of the UAV IRIS quadrotor, powered by a 3S LiPo battery with 5.5 Ah of capacity. Each movement, including moving forward, ascending/descending, and acceleration/deceleration, was analysed to create a model of the energy consumed relative to each UAV action. Besides its movement, the influence of environmental conditions also takes an important role on the autonomy of the UAV, such as the wind speed and its direction, and the atmospheric flow conditions [5]. These factors interfere with the stability of the UAV, resulting in energy losses required to compensate those external forces.
To achieve the least expensive sequence of inspection, it is possible to apply several path planning algorithms. The most adequate approach can depend on the shape of the structure to be inspected. A geometric method may result in appropriate paths for regular structures (like a cube or a cylinder) with a reduced computation time, whereas an heuristic approach may result in optimized paths for more complex structures. Given that the inspection path can be resampled in a set of viewpoints of the structure [18], [19], this problem can be formulated as the well-known Travelling Salesman Problem [15], [20]. Examples of used optimization techniques are genetic algorithms [15], ant colony [21], and particle swarm optimization (PSO) based algorithms [20], [22], among other heuristics and metaheuristics. A possible approach to reduce the complexity of the problem consists in dividing the structure in multiple sub-spaces of inspection. First it is computed the order that each sub-space will be inspected, and secondly it is obtained the path within each sub-space [23].

III. ENERGY EFFICIENT INSPECTION PATH PLANNING
Several steps need to be performed to achieve an energetically optimized path to inspect complex structures with a UAV, which are depicted in Fig. 1, that takes a wind turbine model as an example. Given a structure's model, it is possible to discretize its surface and define a map of the inspection poses, named viewpoints. A UAV's energy model can be defined by relating each movement with the correspondent energy loss. With these two formulation steps completed, an algorithm for path planning was developed, ensuring that the entire surface of interest is covered during the inspection, optimizing the energy consumption. By analysing the data collected during the inspection it is possible to reconstruct the structure's model, as well as to draw a confidence map that allows to pinpoint the confidence on the quality of the data in each location.

A. DEFINING THE MAP OF THE INSPECTION POSES
The map of the inspection poses consists in defining a set of viewpoints that ensure that the entire surface of interest of the structure is completely covered. The definition of these viewpoints depends upon the characteristics of the sensors VOLUME 11, 2023 FIGURE 1. Scheme of the developed inspection path components, including the mapping of the relevant viewpoints, the UAV energy model, the generated inspection path and the correspondent model reconstruction and confidence map. used in the inspection, such as cameras and 3D LiDARs, and also on the structure's shape. To guarantee the coverage of the entire surface of the inspected structure, these sensors constrain the definition of the inspection map, through the following characteristics [16]: • α -Field of View (FoV), represents the angular aperture of the data captured by the sensor; • (I x , I y ) -data resolution, expresses the size of the captured data in the smallest addressable element; • ρ = I x I y -aspect ratio, represents the relation between the width and height of the data captured by the sensor; • f S -sampling frequency, represents the time expenditure between two consecutive data collections. When the UAV is moving at a given working distance (d w ) from the target, a sensor acquires data that covers a projected area (L x × L y ) that depends on its field of view (α) and aspect ratio (ρ), as can be seen in Fig. 2. This relation is expressed in (1) [16]: To ensure a minimum spacial resolution (R = max( I x L x , I y L y ), the maximum working distance (d max w ) must be as follows (2) [16]: (2) To discretize the entire surface of the inspected model, ensuring that the whole structure is covered, the discretization rate (D r ) must be lower than or equal to min(L x , L y ). After the discretization is performed, each point represents an inspection point of the structure that the UAV must meet. The last constraint that the inspection sensors influence is the maximum speed at which the UAV can navigate. This limitation is related to the sampling frequency (f S ) of the sensor, given that the UAV cannot move faster than the ability of the sensor to gather the information in each inspection point. Given that, the speed of the UAV must be limited according to (3) [16]: Each inspection point has a correspondent viewpoint that consists of the pose where the UAV must be located and oriented to collect the required information of that inspection point. Each viewpoint is obtained according to the working distance (d w ) and to the orientation of the normal (⃗ n, with only the yaw component) relative to each inspection point. Fig. 3 exhibits the method used to calculate each viewpoint. The total inspection path, defined by the complete set of viewpoints, is obtained by applying this method to the whole 3D surface of the inspected structure.

B. UAV ENERGY MODEL
The energy model relates the energy loss of each movement of the UAV, which can be calculated according to (4): where the total energy loss (E total ) corresponds to the sum of the following components: • E d -displacement cost: -E d xy -linear displacement on the xOy plane; -E d z -linear displacement on the z axis; -E d yaw -angular displacement on the z axis.
• E a -acceleration cost: -E a xy -linear acceleration on the xOy plane; -E a z -linear acceleration on the z axis; -E a yaw -angular acceleration on the z axis. Regarding the displacement cost, the linear movement according to the xOy plane and z axis was analysed, as well as the angular motion along the z axis, that corresponds to yaw in Euler angles. The change of direction of the UAV is related with its acceleration and deceleration, and corresponds to another significant factor of energy consumption. The acceleration costs that were analysed depend on the linear acceleration on both xOy plane and z axis, and also on angular acceleration on the z axis (yaw). This cost model enables to objectively quantify the cost of moving from a current pose to another pose. Note that, for calculation purposes, the speed of each movement of the UAV is considered constant, and that each acceleration/deceleration is considered instantaneous.

1) DISPLACEMENT COST
The displacement cost is related with linear and angular movements performed by the UAV, from a given pose p i to another pose p i+1 . Fig. 4 exhibits the displacement diagram of all considered movements that a UAV can perform. d xy corresponds to the displacement on the xOy plane, d z + and d z − corresponds to the ascending and descending displacement on the z axis, respectively, and θ yaw corresponds to the rotation on the z axis. The relation of these movements with its displacement cost (E d ) can be expressed according to (5): where E d xy corresponds to the linear displacement on the xOy plane, E d z to the linear displacement on the z axis, and E d yaw to the angular displacement on the z axis. The energy cost of each linear displacement component can be calculated according to (6) [16]: The energy loss caused by the movement according to the xOy plane (E d xy ) is related to its linear velocity (v xy ), its displacement ( d xy ), and the power consumption as a function of the linear velocity (P d xy (v xy )).
The energy loss caused by the movement according to the z axis (E d z ) depends on its direction, resulting in a higher energy cost when the UAV is ascending (E d z + ) than when it is descending (E d z − ). This cost is related to the ascending/descending velocity (v z + and v z − ), its displacement, ( d z + and d z − ), and the power consumption as a function of the linear velocity (P d z + (v z + ) and P d z − (v z − )). Note that the total energy cost of the movement along the z axis (E d z ) corresponds to (7): The energy cost of the angular displacement component can be calculated according to (8) [16]: The energy loss caused by the rotation according to the z axis (E d yaw ) is related to its angular velocity (ω yaw ), its displacement ( θ yaw ), and the power consumption as a function of the angular velocity (P d yaw (ω yaw )).

2) ACCELERATION COST
The acceleration cost is related with the change of direction of the UAV during its path. Given a current pose p i , a previous pose p i−1 , and the next pose p i+1 , it is possible to detect a change of direction and the correspondent acceleration or deceleration by calculating the distance between the displace- Fig. 5 exhibits an acceleration diagram for the considered changes of direction, and the corresponding accelerations and decelerations. v l and a l correspond to the velocity and acceleration components on axis l = {x, y}, respectively. Considering that the UAV has a given speed v x > 0 and v y = 0 at the current pose, moving to neighbour pose 2 means keeping its current movement, not being necessary to accelerate or decelerate in the direction of any axis (a x = a y = 0), which does not imply any acceleration cost. On the other hand, moving to neighbour pose 1 implies a change of direction, meaning that the UAV has to decelerate in the x axis (a x < 0) and accelerate in the y axis (a y > 0), which causes an energy loss. The relation of these movements with its acceleration cost (E a ) can be expressed according to (9): where E a xy corresponds to the cost of the linear acceleration on the xOy plane, E a z to the cost of the linear acceleration on the z axis, and E a yaw to the cost of the angular acceleration on the z axis. φ represents the difference on the polar angle between two consecutive displacement vectors. Its value allows to detect change of direction between the xOy plane and the z axis. As cos(φ) 2 + sin(φ) 2 = 1, (9) represents the relation between the acceleration along both directions.
The energy cost of each linear acceleration component (E a τ ) can be calculated according to (10) [16]: where K τ is a multiplication factor that depends on the change of direction, P a τ (v τ ) is the power consumption as a function of the velocity, and t a τ is the time that it takes to accelerate/decelerate until the desired velocity is reached. When the UAV is moving in the xOy plane, it is possible to detect changes of direction by calculating the difference of the direction between two consecutive displacement vectors.
|θ| ∈ [0, π] represents that change on the direction of movement in the xOy plane. The energy loss caused by the linear acceleration/deceleration according to the xOy plane (E a xy ) is related to its power consumption as a function of the linear velocity (P a xy (v xy )) and to the time that it takes to accelerate/decelerate until the desired velocity is reached ( t a xy ). K xy ∈ [0, 1] is a multiplication factor that assumes the values presented in (11): (11) Figure 6 helps to explain the values that K xy takes, where θ = 0 is represented with the dark red arrow and corresponds to no change of direction, 0 < |θ| < π 2 is represented in light blue and corresponds to a smooth change of direction, and π 2 ≤ |θ| ≤ π is represented in purple and corresponds to a sharp change of direction. When θ = 0 there is no acceleration in the xOy plane in the UAV movement (a xy = 0), so K xy takes value 0. When π 2 ≤ |θ| ≤ π the UAV has to decelerate and then accelerate (a xy ̸ = 0), as it has to completely change its direction of movement, resulting in K xy FIGURE 6. Different directions of movement in the xOy plane, according to θ . θ = 0 is represented with the dark red arrow and corresponds to no change of direction, 0 < |θ| < π 2 is represented in light blue and corresponds to a smooth change of direction, while π 2 ≤ |θ| ≤ π is represented in purple and corresponds to a sharp change of direction. equal to 1. When 0 < |θ| < π 2 , the UAV's acceleration is an increasing function of θ (a xy ̸ = 0). K xy changes according to the ellipsoidal function presented in (12), that depends on θ: The acceleration/deceleration cost according to the z axis (E a z ) depends if the UAV is ascending (E a z + ) or descending (E a z − ). It is possible to detect changes of direction by calculating the difference of the direction between two consecutive displacement vectors. |φ| ∈ [0,π] represents that change on the direction of movement in the z axis. The acceleration cost is related to its power consumption as a function of the linear velocity (P a z + (v z + ) and P a z − (v z − )) and to the time that it takes to accelerate/decelerate until the desired velocity is reached ( t a z + and t a z − ).
is a multiplication factor that assumes the values presented in (13): When φ = 0 there is no acceleration in the z axis in the UAV movement (a z = 0), so K ζ takes value 0. When φ ̸ = 0 the UAV has to accelerate/decelerate on the z axis (a z ̸ = 0), resulting in K ζ equal to 1. Note that the total energy cost given by the acceleration along the z axis (E a z ) corresponds to (14): The energy cost of the angular acceleration component can be calculated according to (15) [16]: It is possible to detect changes of angular directions by calculating the difference of the angular direction between two consecutive displacement vectors. |ψ| ∈ [0,π] represents that change on the angular direction of movement in the z axis. The energy loss caused by the angular acceleration according to the z axis (E a yaw ) is related to its power consumption as a function of the angular velocity (P a yaw (ω yaw )) and to the time that it takes to accelerate/decelerate until the desired angular velocity is reached ( t a yaw ). K yaw ∈ [0, 1], is a multiplication factor that assumes the values presented in (16): When ψ = 0 there is no angular acceleration in the UAV movement according to the z axis (a yaw = 0), so K yaw takes value 0. When ψ ̸ = 0 the UAV has to accelerate/decelerate on the z axis (a yaw ̸ = 0), resulting in K yaw equal to 1.

C. INSPECTION PATH PLANNING ALGORITHM
Given a set of viewpoints and the UAV energy model that relates the energy loss of moving from a viewpoint to another, it is possible to develop an algorithm that generates a path for the UAV to inspect a structure, optimizing its energy loss. This problem was formulated as a new variant of the Travelling Salesman Problem that is defined as ATSP-PL (Asymmetric Travelling Salesman Problem, with Precedence Loss). The Asymmetric classification is due to the fact that the cost of moving from viewpoint p i to viewpoint p i+1 is not necessarily the same of moving from viewpoint p i+1 to viewpoint p i . The Precedence Loss designation is related with the fact that the cost of moving from the current viewpoint (p i ) to the next viewpoint (p i+1 ) depends on the previous viewpoint (p i−1 ), given the significant difference that the acceleration cost implies in the calculation of the total energy loss.
The algorithm developed in this research, named EEPPA (Energy Efficient Path Planning Algorithm), proposes a method to solve this novel TSP formulation (the ATSP-PL), and is divided into two main procedures: • Spiral Procedure -a greedy and time efficient method, specially suited to cover regular shaped structures, such as cylinders or cubes; • ATSP Solver followed by a Local Search Procedurethe original ATSP-PL problem is relaxed and solved as a traditional ATSP, followed by a Local Search Procedure that adapts the obtained solution to take into account the Precedence Loss property.

1) ATSP-PL FORMULATION
As in the classical TSP, the solution of the ATSP-PL consists of a sequence of viewpoints organized by its order of inspection. The first viewpoint is the only viewpoint whose position in the sequence is previously known, as it is user defined. Each solution sequence has an associated total energy cost and total linear displacement. While performing the path, the UAV has to meet all viewpoints exactly once, with the objective of minimizing the energy loss during the execution of the complete path. With this information it is possible to modulate this problem as follows [24]: x ij = 1, if the path goes from viewpoint i to j, 0, otherwise, i∈Q j̸ =i, j∈Q x ij ≤ |Q| − 1, ∀Q ⊆ {1, · · · , N }, |Q| ≥ 2, (20) where N is the total number of viewpoints. The constraint presented in (20) ensures that no sub-sequences (Q) are formed, resulting in a single sequence of viewpoints [24]. VOLUME 11, 2023 The objective function is the following: where c kij is the energy cost of moving from viewpoint i to j, given a previous viewpoint k. Note that the constraints presented in (17) to (20) are equivalent to the standard ATSP formulation, while only the objective function, presented in (21), has been altered to fit the ATSP-PL formulation. This new formulation minimizes the cost of moving from the current viewpoint to the next viewpoint, given a previous viewpoint, justifying the Precedence Loss designation.

2) EEPPA -ENERGY EFFICIENT PATH PLANNING ALGORITHM
The developed EEPPA algorithm, whose flowchart is presented in Fig. 7, takes as input the map of viewpoints of the inspected structure and the energy model of the UAV, and outputs the inspection path as a sequence of viewpoints. The map of viewpoints comprises all inspection poses that the UAV must meet, while the energy model enables to calculate the energy consumption of moving from one pose to the following one, which not only iteratively assists on selecting the subsequent viewpoint, but also allows to determine the total energy loss of the path produced. The inspected structure can be divided into sub-structures to allow the application of different methods that may fit better depending on the shape of the sub-structure to be inspected. This division is performed in the first stage of the algorithm (represented with the blue diamond in Fig. 7) recurring to model fitting and segmentation algorithms [25], capable of detecting objects with regular shapes, such as cylinders or cubes, allowing to divide the complete structure into regular shaped and complex shaped sub-structures. The path planning for the regular shaped sub-structures is generated with the Spiral Procedure method (represented in blue on the left side of Fig. 7), due to its ability to achieve a feasible solution in a short computation time. This greedy method is led by a simple dispatching rule with a considerable cost-benefit ratio between the quality of the solution and the time spent to achieve it. The flowchart of Spiral Procedure method is presented in Fig. 8. The Spiral Procedure method is similar to a cheapest neighbour approach. A solution is generated by iteratively selecting the least expensive neighbour viewpoint to move to, according to the energy loss. Each viewpoint neighbourhood is defined by the euclidean distance within a subset of viewpoints with similar height. When all the viewpoints of that subset are selected, the last viewpoint's neighbourhood is expanded to the next subset of viewpoints, and this process is repeated until all viewpoints have been covered.
A different approach is used for the inspection of a sub-structure with a complex shape. By relaxing the Precedence Loss property of the ATSP-PL, it is possible to formulate and solve this problem as a traditional ATSP (represented with the green background on the right side of Fig. 7). As the ATSP is a well-known problem, using one of the several optimized open source solvers is the most suitable option. The selected open source ATSP Solver was the OR-Tools from Google Developers [26], in particular the Routing algorithm. This solver allows to select the greedy algorithm that generates the initial solution, and the metaheuristic algorithm capable of improving it. Its input is a cost matrix with a size of N × N (N represents the number of viewpoints in the inspection map), in which each cell has the energy cost of moving from one viewpoint to another. For the purpose of this problem, the Savings algorithm [27] was selected for the initial solution algorithm, and the Guided Local Search for the metaheuristic. This metaheuristic is generally the most efficient metaheuristic for routing problems [26], as it enables the solver to escape from a local minimum, proceeding with the search for better solutions. To adapt and improve the solution obtained with the ATSP Solver for the ATSP-PL formulation, a Local Search method is applied (represented in red on the right side of Fig. 7). The flowchart of this method is presented in Fig. 9. The local search is a well-known approach that is able to enhance a solution by searching for an improving one within its neighbourhood. In the developed Local Search method, two different approaches were followed. One is based on the Relocate algorithm [28], and consists of swapping pairs of viewpoints of a sequence within a defined neighbourhood (represented on the left side of Fig. 9). The other is based on the 2-Opt algorithm [29], [30], and consists in flipping the order of a sub-sequence between two viewpoints within a defined neighbourhood (represented on the right side of Fig. 9). Fig. 10 exhibits an example of a flipped sub-sequence (3,6). This Local Search method will be continuously seeking for an enhanced solution, as long as each of its approaches is able to find an improving solution at each iteration.

IV. SIMULATION SCENARIO RESULTS
The EEPPA algorithm was applied in the Gazebo simulator, which allows to replicate real environment conditions and to represent the UAV and its sensors with realistic models, as well as a structure to be inspected. This section displays the generated map of viewpoints and the obtained path of inspection, presenting the correspondent energy  consumption, total displacement and time spent during the inspection. By processing the data collected during the execution of the inspection path, it is possible to reconstruct the structure's model and to map the confidence on the quality of the data at each location. The results obtained in the conducted simulation were performed on a computer with the following specifications: • CPU -Intel ® Core ™ i7 7700HQ Quad-Core, 2.80 GHz; • GPU -Nvidia GeForce GTX 1050 4GB GDDR5; • RAM -16GB DDR4 SDRAM. The UAV model used in the simulator corresponds to the IRIS quadrotor controlled by an ArduPilot autopilot, powered by a 3S LiPo battery (5.5 Ah) [16]. For the inspection map VOLUME 11, 2023 definition, the considered sensors were a monocular camera and a 3D LiDAR, with the following characteristics: For the purpose of this demonstration, the selected structure to be inspected consists of an offshore wind turbine, whose model is depicted in Fig. 11. This model consists of the WindFloat design [31], that has a column height of 95 m and a blade length of 82 m. Different tests were performed with two different blades configurations, shown in Fig. 11b) and Fig. 11c). The top figure corresponds to the horizontal configuration, and the bottom one to the vertical configuration. Given the previously defined parameters and this wind turbine model, it is possible to generate the map of viewpoints of the structure's surface, for the UAV inspection. Fig. 12 exhibits the sequence of operations to obtain the map of viewpoints for the horizontal configuration. Fig. 12a) shows a section of the wind turbine model, while Fig. 12b) depicts the correspondent inspection points of its surface (colored dots),  with a discretization rate (D r ) of 2 m. Fig. 12c) shows the viewpoints obtained (red arrows) for each inspection point, for a working distance (d w ) of 5 m. Fig. 12d) depicts the same viewpoints map and the correspondent wind turbine model. The number of viewpoints generated to inspect this wind turbine model is close to 1300. This procedure is analogous for both vertical and horizontal blade configuration, and for any other structure to be inspected.

A. PATH FOR THE INSPECTION OF THE STRUCTURE
Given the map of viewpoints of the structure and the energy model of the UAV, it is possible to apply the developed path planning algorithm. The results obtained by the EEPPA algorithm are compared with the results obtained with two distinct algorithms, OR-Tools and Cheapest Neighbour (CN). The OR-Tools corresponds to the previously mentioned ATSP solver [26], that obtains a solution based on the relaxation of the original problem (does not consider the acceleration cost of the UAV, that corresponds to the Precedence Loss characteristic of the ATSP-PL problem formulation). The Cheapest Neighbour algorithm is a greedy approach that generates a solution by iteratively selecting the least expensive neighbour viewpoint to move to, according to the energy loss. Each viewpoint neighbourhood is defined by the euclidean distance. Fig. 13 depicts the paths obtained for both wind turbine configurations, with the EEPPA, CN and OR-Tools algorithms. The paths generated are represented by the green line and each viewpoint is represented by a red arrow. By analysing this figure it is noticeable that the paths generated by EEPPA are more sequenced and with smoother movements, while the paths generated by CN and OR-Tools present rougher movements and a less sequential path. The inspection of both wind turbine configurations was performed by the UAV in the simulation environment, where the paths generated by EEPPA algorithm were followed. 1 Tables 2 and 3 present a summary of the values obtained 1 Videos of both inspection paths available at https://youtu.be/ZK6-UtJwNns (horizontal configuration) and https://youtu.be/on-gqbnurrM (vertical configuration), displayed with ten times of the real speed. The wind turbine model is represented in grey and the UAV in yellow. On the left it is visible the images captures by the camera. The point cloud obtained by the LiDAR are represented by the white dots. The red arrows represent the pose (position and orientation) of the UAV at each instance of the inspection. in the simulation scenario for the paths with the vertical and horizontal configurations, respectively, regarding energy loss, displacement and time spent during the inspection.
For the vertical configuration, the EEPPA algorithm resulted in a improvement of 27.4% regarding the CN, and of 26.6% regarding the OR-Tools. For the horizontal configuration, the improvement regarding the CN and the OR-Tools are of 28.3% and 24.8%, respectively. The difference in the energy loss between the paths obtained with the two wind turbine configurations can be explained by the orientation of the blades. In the vertical configuration one of the blades is positioned perpendicularly to the ground, whereas in the horizontal configuration one of the blades is positioned in parallel. The horizontal blade allows the UAV to travel in a straight line more often, not having the necessity to change its direction of movement multiple times, resulting in a smoother path. The vertical blade forces the UAV to successively accelerate and decelerate during its inspection, increasing its energy loss. The paths generated by EEPPA not only allowed a higher energy efficiency, which is the mains goal, but also to achieve faster and shorter paths. This can be explained by the fact that the EEPPA generates smoother and more sequenced paths, avoiding longer distances between two consecutive viewpoints. Fig. 14 shows an example of the improvement generated by the Local Search method in the EEPPA algorithm. Fig. 14a) shows a segment of the path generated by the OR-Tools ATSP Solver, that is optimized to only reduce the displacement cost. Fig. 14b) shows the path obtained for the same path segment after applying the Local Search method, that also takes into account the acceleration cost. Although this change increases the total displacement, the loop generated by the Local Search decreases the number of changes in direction, which smooths the movement of the UAV and decreases the energy loss that this segment of the path would produce.

B. 3D MODEL RECONSTRUCTION
During the inspection of the wind turbine in the simulation scenario, the UAV was equipped with a LiDAR sensor which is able to collect a 3D point cloud of the surrounding environment. With a point cloud registration method, it is possible to concatenate successive point clouds collected at each instance, reconstructing the entire structure of inspection, which results in a complete representation of the wind turbine surface. The used method consists in pre-aligning each point cloud based on the relative position between the UAV and the turbine (given by a RTK GNSS sensor), and then fine tuning the alignment with the ICP algorithm [32], [33], [34], [35]. This 3D reconstruction allows to identify deviations on the blades or tower, or even detect if any part is missing or damaged. Fig. 15 depicts the point cloud of the reconstructed model of the inspected wind turbine, for the vertical blade configuration. 2 A confidence parameter is assigned to each inspection point, allowing to create a confidence map of the inspection data. This confidence parameter is related with the fitting ratio returned by the ICP algorithm, that depends on the pre-alignment given by the relative position between the UAV and the wind turbine obtained by the RTK sensor. This confidence map identifies the certainty on the location of the collected data for the reconstructed model, mapping the areas with less reliable information that may need to be double checked. Figure 16 displays the confidence map of the inspection of the wind turbine for the vertical configuration. The green dots represent the poses with the highest confidence in the quality of the collected data, whereas red dots represent the poses with the lowest confidence in the data that was collected.

V. REAL SCENARIO RESULTS
The developed algorithm for inspection path planning was also tested in a real scenario. The used aerial robot is INESC TEC's CROW (Copter Robot for Offshore Wind-farms) UAV. This UAV is based on a quadcopter frame (with a wingspan of 0.7 m, and a maximum payload of 2 Kg) with a PixHawk controller running the ArduPilot software. Its navigation system is composed of a GNSS/RTK sensor, that has a positioning accuracy of 2.5 m and of 0.025 m when paired with a base station. It has an onboard processing unit (Raspberry Pi 4B 8Gb). CROW is supplied by a 6S LiPo battery of 14 Ah, which enables a total flight time of at least 25 minutes (depending on the payload). The values defined for the linear and angular speeds of the UAV during the inspection are the same used in the simulation experiments. For the inspection map definition, the considered sensors were a monocular camera and a 3D LiDAR, with the following characteristics:  Figure 17 depicts CROW with the described sensors attached.
The structure inspected in this experiment was a 25 m high illumination tower of the football field of the Faculty of Sport of the University of Porto (FADEUP), depicted in Fig. 18. Given that a 3D model of the illumination tower was not available, a coarse sketch of its surface was obtained by collecting several point clouds at distinct points of the structure, recurring to the LiDAR mentioned above. To approximate this experience to the problem of inspecting a larger structure (as the wind turbine model used in the simulation experiments) the discretization rate was considerably reduced to increase the number of inspection points, also resulting in the increase of the total number of viewpoints. In particular,  the discretization rate (D r ) was set to 0.2 m horizontally and 0.5 m vertically, which resulted in about 800 viewpoints. The discretized inspection surface of the illumination tower is depicted in Fig. 19a). The working distance (d w ) was set to 3 m, which allowed to obtain the map of viewpoints exhibited in Fig. 19b). Note that the irregularity and roughness of the obtained viewpoints map is due to the coarse sketch used as model of the illumination tower. Similarly to the simulation scenario, the results obtained by the EEPPA algorithm were compared with the results obtained with the OR-Tools and the CN algorithms. Due to the regular shape of the structure used in these tests, the results obtained when applying only the Spiral Procedure algorithm were also analysed. The path obtained with the EEPPA algorithm for the inspection of the   illumination tower is displayed in Fig. 19f), whereas the paths obtained with the CN, OR-Tools and Spiral Procedure are depicted in Fig. 19c), Fig. 19d) and Fig. 19e), respectively. Table 4 presents a summary of the values obtained in this real scenario, for the paths generated for the inspection of the illumination tower, regarding energy loss, displacement and time spent during the inspection.
The EEPPA algorithm resulted in an improvement of 20.7%, of 8.7% and of 4.6%, when compared with the OR-Tools, CN and Spiral Procedure algorithms, respectively. Although the total displacement generated by the path of EEPPA is considerable higher than the displacement of CN and OR-Tools, the energy loss is considerably lower. This can be explained once again by the fact that EEPPA is an algorithm focused on minimizing the energy consumption during the inspection, instead of the total displacement. The path generated by EEPPA is smoother and more sequential, preventing the UAV from successively accelerating and decelerating its movement. Due to the regular shape of the illumination tower, the results obtained with EEPPA are similar to the ones obtained with only the Spiral Procedure. The main explanation for these results is related with the fact that EEPPA algorithm applies the Spiral Procedure to produce the path for the inspection of the column, only differing in the method used to generate the path to inspect the lamps on the top of the column. Nevertheless, even this small difference in the path generated for the lamps resulted in a more energy efficient path. Fig 20 depicts CROW during an inspection of the illumination tower, and a sample of the collected data. By registering the point cloud collected during the inspection it was possible to reconstruct the illumination tower module that is depicted in Fig. 21.

VI. CONCLUSION
This work presented an energy oriented path planning algorithm that minimizes the energy loss of a UAV during the inspection of a structure. The surface of a given model to be inspected was discretized, allowing to obtain an inspection map. This map is composed of a set of viewpoints that depend on the structure surface shape, and also on the sensors' characteristics. To optimize the path to be performed in terms of energy loss, a complete energy model of the UAV was defined. This model considers the energy consumption of distinct components, such as displacement and acceleration, related with the changes of direction of the movement of the UAV. With the list of viewpoints and the energy model defined, it was possible to solve the path planning problem with the developed algorithm, named EEPPA. The experiments performed in both simulation and real scenarios allowed to validate the developed algorithm, which obtained improvements over 24% and 8%, respectively, when compared with other methods. The simulated experiments were performed with realistic models of wind turbines and a UAV, while the real experiments were performed with a real UAV and an illumination tower. The collected inspection data, in particular the point cloud, allowed to perform model reconstructions of the inspected structures, as well as creating a confidence map that relates the confidence on the data collected in each viewpoint with its localization.