Distributed approaches for coverage missions with multiple heterogeneous UAVs for coastal areas.
Abstract
This Thesis focuses on a high-level framework proposal for heterogeneous aerial, fixed wing teams of robots, which operate in complex coastal areas. Recent advances in the computational capabilities of modern processors along with the decrement of small scale aerial platform manufacturing costs, have given researchers the opportunity to propose efficient and low-cost solutions to a wide variety of problems. Regarding marine sciences and more generally coastal or sea operations, the use of aerial robots brings forth a number of advantages, including information redundancy and operator safety. This Thesis initially deals with complex coastal decomposition in relation with a vehicles’ on-board sensor. This decomposition decreases the computational complexity of planning a flight path, while respecting various aerial or ground restrictions. The sensor-based area decomposition also facilitates a team-wide heterogeneous solution for any team of aerial vehicles. Then, it proposes a novel algorithmic approach of partitioning any given complex area, for an arbitrary number of Unmanned Aerial Vehicles (UAV). This partitioning schema, respects the relative flight autonomy capabilities of the robots, providing them a corresponding region of interest. In addition, a set of algorithms is proposed for obtaining coverage waypoint plans for those areas. These algorithms are designed to afford the non-holonomic nature of fixed-wing vehicles and the restrictions their dynamics impose. Moreover, this Thesis also proposes a variation of a well-known path tracking algorithm, in order to further reduce the flight error of waypoint following, by introducing intermediate waypoints and providing an autopilot parametrisation. Finally, a marine studies test case of buoy information extraction is presented, demonstrating in that manner the flexibility and modular nature of the proposed framework.
Full text
Proyecto Fin de Carrera Ingeniería de Telecomunicación Formato de Publicación de la Escuela Técnica Superior de Ingeniería Autor: F. Javier Payán Somet Tutor: Juan José Murillo Fuentes Dep. Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2013 Tesis Doctoral Ingeniería de Automática, Electrónica y de Telecomunicación Distributed approaches for coverage missions with multiple heterogeneous UAVs for coastal areas Autor: Fotios Balampanis Directores: Jesús Ivan Maza Alcañiz, Aníbal Ollero Baturone Departamento de Ingeniería de Sistemas y Automática Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, Mayo 2018
Tesis Doctoral Ingeniería de Automática, Electrónica y de Telecomunicación Distributed approaches for coverage missions with multiple heterogeneous UAVs for coastal areas Autor: Fotios Balampanis Director: Jesús Ivan Maza Alcañiz Profesor Titular de Universidad Director: Aníbal Ollero Baturone Catedrático de Universidad Departamento de Ingeniería de Sistemas y Automática Escuela Técnica Superior de Ingeniería Universidad de Sevilla Mayo 2018
Tesis Doctoral: Distributed approaches for coverage missions with multiple heterogeneous UAVs for coastal areas Autor: Fotios Balampanis Director: Jesús Ivan Maza Alcañiz Director: Aníbal Ollero Baturone El tribunal nombrado para juzgar la Tesis arriba indicada, compuesto por los siguientes doctores: Presidente: Vocales: Secretario: acuerdan otorgarle la calificación de: El Secretario del Tribunal Fecha:
A mi hermano
Acknowledgments Wr iting an acknowledgements section to a doctoral dissertation may seem, at a personal level, more important than the Thesis itself. Perhaps it is the only chance to thank all the people who joined you and are an integral part of it. What we are, in essence, is a set of experiences and personalities we met on our way. This is also an element of science: it does not belong to anyone specific because it belongs to everyone. Everyone has in their own way contributed to the evolution of human intelligence. I wish I was laconic on my acknowledgements, but I will probably not be. How could I omit to thank 1 The Queen, The Who, Pavlos (x2), Richard Bach, Matisse, Nikolas, Tom, Darwin, Aristophanes and Robin? Initially, I would like to thank my professors, Ivan Maza and Anibal Ollero, for their confidence and guidance. Their experience as well as their clear vision always showed me the path, while they managed to understand me very early on. I must also thank them for their patience. Then, I would like to thank the colleagues at the University of Seville, the group and the laboratory, for their smile and mood to help me made me feel at home from the very first moment. I would also like to thank the colleagues and travelling companions of the MarineUAS program, a group of extremely smart people. They made the trip extremely pleasant, and they managed to open my mind a little more, they took it a little further. A special thank you to Margarita Costa Faria whose friendship and polite nature embraced the moments of my nagging and reminded me of trusting people. In particular, I would like to thank Leopoldo Rodriguez, the parallel life co-traveller of these three years. Without him this dissertation would not have been realised. In all the critical moments he made me remember that things are simpler and there is always room for a little more chili. I would also like to thank Juan Braga and Rafael Salmoral, partners in crime, for their help and patience with regard to my Spanish level in this journey. Finally, I would like to thank all the ETSI staff and in particular the cafeteria staff who have given me countless quantities of coffee for three years. I want to apologize to them because I never learned 1 The cover depicts the port of Agios Kirikos and the statue of Icarus in Icaria island. The photo has been shot by Chrisa Karafti III
Abstract Th is Thesis focuses on a high-level framework proposal for heterogeneous aerial, fixed wing teams of robots, which operate in complex coastal areas. Recent advances in the computational capabilities of modern processors along with the decrement of small scale aerial platform manufacturing costs, have given researchers the opportunity to propose efficient and low-cost solutions to a wide variety of problems. Regarding marine sciences and more generally coastal or sea operations, the use of aerial robots brings forth a number of advantages, including information redundancy and operator safety. This Thesis initially deals with complex coastal decomposition in relation with a vehicles’ on-board sensor. This decomposition decreases the computational complexity of planning a flight path, while respecting various aerial or ground restrictions. The sensor-based area decomposition also facilitates a team-wide heterogeneous solution for any team of aerial vehicles. Then, it proposes a novel algorithmic approach of partitioning any given complex area, for an arbitrary number of Unmanned Aerial Vehicles (UAV). This partitioning schema, respects the relative flight autonomy capabilities of the robots, providing them a corresponding region of interest. In addition, a set of algorithms is proposed for obtaining coverage waypoint plans for those areas. These algorithms are designed to afford the non-holonomic nature of fixed-wing vehicles and the restrictions their dynamics impose. Moreover, this Thesis also proposes a variation of a well-known path tracking algorithm, in order to further reduce the flight error of waypoint following, by introducing intermediate waypoints and providing an autopilot parametrisation. Finally, a marine studies test case of buoy information extraction is presented, demonstrating in that manner the flexibility and modular nature of the proposed framework. XI
Resumen Es ta tesis se centra en la propuesta de un marco de alto nivel para equipos heterogéneos de robots de ala fija que operan en áreas costeras complejas. Los avances recientes en las capacidades computacionales de los procesadores modernos, junto con la disminución de los costes de fabricación de plataformas aéreas a pequeña escala, han brindado a los investigadores la oportunidad de proponer soluciones eficientes y de bajo coste para enfrentar un amplio abanico de cuestiones. Con respecto a las ciencias marinas y, en términos más generales, a las operaciones costeras o marítimas, el uso de robots aéreos conlleva una serie de ventajas, incluidas la redundancia de la información y la seguridad del operador. Esta tesis trata inicialmente con la descomposición de áreas costeras complejas en relación con el sensor a bordo de un vehículo. Esta descomposición disminuye la complejidad computacional de la planificación de una trayectoria de vuelo, al tiempo que respeta varias restricciones aéreas o terrestres. La descomposición del área basada en sensores también facilita una solución heterogénea para todo el equipo para cualquier equipo de vehículos aéreos. Luego, propone un novedoso enfoque algorítmico de partición de cualquier área compleja dada, para un número arbitrario de vehículos aéreos no tripulados (UAV). Este esquema de partición respeta las capacidades relativas de autonomía de vuelo de los robots, proporcionándoles una región de interés correspondiente. Además, se propone un conjunto de algoritmos para obtener planes de puntos de cobertura para esas áreas. Estos algoritmos están diseñados teniendo en cuenta la naturaleza no holonómica de los vehículos de ala fija y las restricciones que impone su dinámica. En ese sentido, esta Tesis también ofrece una variación de un algoritmo de seguimiento de rutas bien conocido, con el fin de reducir aún más el error de vuelo del siguiente punto de recorrido, introduciendo puntos intermedios y proporcionando una parametrización del piloto automático. Finalmente, se presenta un caso de prueba de estudios marinos de extracción de información de boyas, que demuestra de esa manera la flexibilidad y el carácter modular del marco propuesto. XIII
Σύνοψη Ηπαρούσα διατριβή επικεντρώνεται στην πρόταση ενός δομικού πλαισίου υψηλού επιπέδου για ετερογενή, σταθερής πτέρυγας εναέρια ρομπότ, που λειτουργούν σε πολύπλοκες παράκτιες περιοχές. Οι πρόσφατες εξελίξεις στις υπολογιστικές δυνατότητες των σύγχρονων επεξεργαστών, σε συνδυασμό με το μειωμένο κόστος κατασκευής των αεροσκαφών μικρής κλίμακας, δίνουν στους ερευνητές την ευκαιρία να προτείνουν αποτελεσματικές και χαμηλού κόστους λύσεις σε μια ευρεία ποικιλία προβλημάτων. ΄Οσον αφορά τις θάλασσιες επιστήμες αλλά και γενικότερα τις παράκτιες ή θαλάσσιες δραστηριότητες, η χρήση εναέριων ρομπότ παρουσιάζει ορισμένα πλεονεκτήματα, μεταξύ των οποίων είναι και η πλεονάζουσα πληροφορία και η ασφάλεια των χειριστών. Αρχικά, αυτή η διατριβή ασχολείται με την αποσύνθεση μιας πολύπλοκης παράκτιας περιοχής, σε σχέση με τον αισθητήρα της ρομποτικής πλατφόρμας. Αυτή η αποσύνθεση έχει σαν αποτέλεσμα τη μείωση της υπολογιστικής πολυπλοκότητας του σχεδιασμού ενός σχεδίου πτήσης καθώς και την ασφαλή αντιμετώπιση διάφορων εναέριων ή επίγειων περιορισμών. Η αποσύνθεση της περιοχής που βασίζεται στον αισθητήρα του οχήματος, ενισχύει επίσης τη συνοχή μιας καθολικής λύσης για ετερογενείς ομάδες εναέριων οχημάτων. Στη συνέχεια, προτείνει μια νέα αλγοριθμική προσέγγιση για τη διαίρεση μιας περιοχής, για οποιοδήποτε αριθμό μη επανδρωμένων εναέριων οχημάτων. Αυτή η πρόταση διαίρεσης και διαμοιρασμού της περιοχής, σέβεται τις σχετικές ικανότητες αυτονομίας της ομάδας των οχημάτων, παρέχοντας σε κάθε ένα μια περιοχή ενδιαφέροντος. Επιπλέον, προτείνεται ένα σύνολο αλγόριθμων για την απόκτηση σχεδίων πτήσης κάλυψης μιας ολόκληρης υπο-περιοχής. Αυτοί οι αλγόριθμοι έχουν σχεδιαστεί λαμβάνοντας υπ΄όψιν τον μη ολονομικό χαρακτήρα των οχημάτων σταθερής πτέρυγας και τους περιορισμούς που επιβάλουν τα δυναμικά χαρακτηριστικά τους. Επίσης, αυτή η διατριβή προτείνει μια παραλλαγή ενός γνωστού path tracking αλγόριθμου, προκειμένου να μειωθεί περαιτέρω το σφάλμα πτήσης, εισάγοντας ενδιάμεσα σημεία αναφοράς στο σχέδιο πτήσης και παρέχοντας μια παραμετροποίηση του αυτόματου πιλότου. XV
XVI Σύνοψη Τέλος, παρουσιάζεται μια μελέτη σε ένα σενάριο θαλάσσιων επιστημών, η οποία αποσκοπεί στην εξαγωγή πληροφορίας από θαλάσσιους σημαντήρες (buoys), παρουσιάζοντας έτσι την ευελιξία και την αρθρωτή φύση του προτεινόμενου πλαισίου.
Contents Abstract XI Resumen XIII Σύνοψη XV 1 Introduction 1 1.1 Motivation 1 1.2 Objectives 3 1.3 Outline and main contributions 4 1.4 Thesis framework 6 1.5 Scientific output 8 Journal publications 8 Conference publications 8 University stay 9 1.6 Participation in training and other dissemination activities 9 Workshops 9 Poster presentations 10 2 Decomposition of complex coastal areas 11 2.1 Introduction 11 2.2 Related work 12 2.3 Complex coastal areas 16 2.3.1 A safe separation strategy 17 2.3.2 A Decomposition proposal 19 2.4 A Constrained Delaunay Triangulation method 19 2.4.1 FoV based decomposition 22 Simple sensor scenario 22 More complex scenarios 24 2.4.2 Sensor specific - orientation agnostic segregation cells 25 Considered model for the on-board sensors 26 XVII
XVIII Contents 2.4.3 Lloyd optimisation grid enhancement 27 2.5 Results 29 Comparison with a square grid method 30 2.6 Conclusions 32 3 Area partitioning for multiple vehicles 33 3.1 Related work 34 3.2 Partitioning 35 3.2.1 Baseline Area Partitioning Algorithm 37 3.2.2 Reverse Watershed Schema 39 3.3 Deadlock handling 39 3.4 Results 44 3.5 Conclusions 55 4 Coverage plans 57 4.1 Related work 58 4.2 Waypoint prioritization 59 4.2.1 Spiral coverage plans 60 4.2.2 Spiral plans and area characteristics 61 4.3 Results 65 4.3.1 Simulation framework 66 4.3.2 Comparison with a grid decomposition and simulations 69 Inlet coefficient method 70 Sample area 71 Actual area scenario 74 Configuration comparison 76 4.3.3 Field experiments 76 4.4 Conclusions 78 5 Path tracking 81 5.1 Related work 81 5.2 Path tracking and following 83 5.3 An offline-online hybrid approach 85 5.4 Results 88 5.4.1 Simulations 88 Lookahead table extraction 88 Simulated test flight 90 5.4.2 Field experiments 92 A tuning control strategy 93 Pure pursuit flight 98 5.5 Conclusions 99 6 A decentralized framework for coastal missions 101
Contents XIX 6.1 A market based decentralized architecture 102 6.2 Test case: a multi-layer framework for buoy coordinated data acquisition 106 6.2.1 Data acquisition layer and buoy specifications 107 6.2.2 Communication layer as a feedback communication control model 108 6.2.3 Single UAV dynamic role reconfiguration 109 Dynamic team reconfiguration 111 6.3 Experimental setup 112 6.3.1 Buoy 113 6.3.2 Dynamic model for the relay UAVs 113 6.3.3 Energy estimation 115 6.3.4 Throughput estimation 116 6.4 Results 117 6.4.1 QoS Simulations 117 Enhancing link connectivity 118 Saving total propulsion energy 119 6.4.2 Role reassignment and team reconfiguration experiments 120 Team reconfiguration 120 6.5 Conclusions 124 7 Conclusions and future development 125 7.1 Conclusions 125 7.1.1 Decomposition of complex coastal areas 126 7.1.2 Area partitioning for multiple vehicles 127 7.1.3 Coverage Plans 127 7.1.4 Path tracking 128 7.1.5 A distributed framework 128 7.2 Future work 129 List of Figures 131 List of Tables 135 Bibliography 137 Index 145 Glossary 147
6Chapter 1. Introduction holonomic properties of fixed-wing aerial robots and the constraints those impose. Thus, the waypoints lists try to minimise the total path length by not visiting a waypoint or a sub-area twice. At the same time, they try to maximise the overall turn angle of each waypoint-to-waypoint flight plan in order for the vehicles to perform smoother trajectories and reduce the flight over restricted zones. Chapter 5 proposes a variant of a known path tracking algorithm in order to reduce the flight error between waypoints in the flight plan. Limitations of kinematics and dynamics of fixed-wing vehicles dictate the need of methods to follow a predefined path in an optimal manner. While literature provides a wide variety of path tracking algorithms, practise with Commercial-Off-The-Shelf (COTS) platforms shows that uncertainties in motion along with the need of fine tuning of controllers, allow further investigation for practices which can reduce the error. In that Chapter, an off-line method of preflight fine tuning is presented; then, intermediate waypoints at a specific rate are sent during flight. Chapter 6 introduces the aforementioned methods in a framework which can be part of a market-based architecture. In order to validate the flexibility of the methods proposed as well as the capability of being part of a variety of applications, a test-case is presented where a team of UAVs is used in a multi-layered architecture for data acquisition from sea buoys. Finally, Chapter 7 summarises and concludes the Thesis, while proposing some guidelines for further research and experimentation. 1.4 Thesis framework This Thesis has been conducted as part of MarineUAS 1 , a European Union’s Horizon 2020 research and innovation programme, under the Marie Sklodowska-Curie grant agreement No 642153. The purpose of this program is to strengthen the research and training on UAS, having in mind the needs of marine and coastal monitoring of the European shores and remote oceanic and polar regions. The program’s structure included a wide range of academic and industrial partners from several countries, in order to provide the fellows and researchers a high-level training regarding scientific training and skills acquisition. The main motivation of this program is the need to protect and manage the natural environment and marine resources in a sustainable manner, as has been manifested in European legislation such as the European Strategy for Marine and Maritime Research. This programme is characterised as an "Innovative Training Network". As such, its goal is to recruit doctoral fellows and provide them intense training by the following means: 1http://marineuas.eu
1.4 Thesis framework 7 •Training-by-research. • Complementary, interdisciplinary and transferable skills training, e.g. presentation skills training. •Secondment and academic visits. •Hands-on UAS operator training. •Network-wide training events that cover UAS technology. •Rules and regulations. •Operations in non-segregated airspace. •Air traffic management. •Marine and coastal monitoring and science. •Integration of the air, surface and underwater segments. The programme had consortium partners from five universities across Europe as well as seven institutes, organisations and companies. Each of the enlisted early stage researcher fellows was responsible for a project; these projects were mainly focused in control theory and aeronautics. Iceberg detection, performance optimisation of vehicles, cooperative control and cooperative motion planning are some indicating keywords of the projects’ spectrum. In order to fulfil the main objectives of the programme and since the majority of the consortium projects acquired a low-level engineering approach, a complementary organizational level for this Thesis was chosen. As such, this Thesis proposes a highlevel, decentralised framework for planning coverage tasks, having in mind the modularity and flexibility that heterogeneous teams of vehicles should have. The research mainly focuses on non-holonomic, COTS fixed-wing platforms, in order to follow the experimental framework used by the consortium partners and other fellows. Moreover, the choice of focusing mainly on fixed-wing vehicles was taken due to the fact that they present larger autonomy and durability, a crucial consideration in long range missions and harsh environments. This implies that the proposed methods could easily be adopted in other research efforts in the project as well as in real world applications. Furthermore, a critical consideration is adopted in this Thesis: the need for successfully integrating the use of UAS in the future Air Traffic Management (ATM) context. The increased use of these systems in a variety of scientific and commercial applications, implies an increased risk for collisions between them as well as with current vehicles utilising the airspace. Thus, the safe separation between vehicles for collision mitigation and the
8Chapter 1. Introduction appropriate mechanisms for monitoring each of them, implies a proper decomposition and area partition dynamic schema. The proper treatment of complex coastal areas as well as the reduction of flight error are the main considerations for the developed framework. In the next section, a list of the published and submitted peer-reviewed publications is presented, along with a description of the MarineUAS university stay. 1.5 Scientific output The following peer-reviewed publications have been presented in the context of this Thesis: Journal publications • Balampanis, F., Maza, I., & Ollero, A. (2017). Area Partition for Coastal Regions with Multiple UAS. Journal of Intelligent & Robotic Systems, 88 (2-4), 751–766. • Balampanis, F., Maza, I., & Ollero, A. (2017). Coastal Areas Division and Coverage with Multiple UAVs for Remote Sensing. Sensors, 17(4), 808. Conference publications • Balampanis, F., Maza, I., & Ollero, A. (2016, June). Area decomposition, partition and coverage with multiple remotely piloted aircraft systems operating in coastal regions. In Unmanned Aircraft Systems (ICUAS), 2016 International Conference on (pp. 275-283). IEEE. • Rodriguez, L., Balampanis, F., Cobano, J. A., Maza, I., & Ollero, A. (2017). Wind efficient path planning and reconfiguration of UAS in future ATM. In Twelfth USA/Europe Air Traffic Management Research and Development Seminar (ATM 2017) (p. 140). • Balampanis, F., Maza, I., & Ollero, A. (2017, June). Spiral-like coverage path planning for multiple heterogeneous UAS operating in coastal regions. In Unmanned Aircraft Systems (ICUAS), 2017 International Conference on (pp. 617-624). IEEE. • Braga, J., Balampanis, F., Aguiar, P., Sousa, J., Maza, I., & Ollero, A. (2017). Coordinated Efficient Buoys Data Collection in Large Complex Coastal Environments using UAVs. In OCEANS MTS/IEEE, Anchorage, AL, 2017 (OCEANS’17 MTS/IEEE).
1.6 Participation in training and other dissemination activities 9 • Rodríguez, L., Balampanis, F., Cobano, J.A., Maza, I., & Ollero, A. (2017, September). Energy-efficient trajectory generation with spline curves considering environmental and dynamic constraints for small UAS. In 2017 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), Vancouver, BC, 2017, pp. 1739-1745. • Balampanis, F., Aguiar, A. P., Maza, I., & Ollero, A. (2017, October). Path tracking for waypoint lists based on a pure pursuit method for fixed wing UAS. In Workshop on Research, Education and Development of Unmanned Aerial Systems (RED-UAS 2017), Linkoping, pp. 55-59. • Braga, J., Balampanis, F., R. Praveen Jain, Aguiar, P., Sousa, J., Maza, I., & Ollero, A. (2018, May). Dynamic Reconfiguration and Load Management for a UAV Based Buoy Data Extraction Architecture. In OCEANS MTS/IEEE, Kobe (OCEANS’18 MTS/IEEE). • Balampanis, F., Rodríguez, L., Cobano, J.A., Maza, I., & Ollero, A. (2018). A path enhancement experimental framework for heterogeneous non holonomic flying vehicles. Submitted - Under Review. Moreover, in the context of the objectives of the MarineUAS project, the following academic stay and workshops have been attended: University stay An academic visit was made at the Underwater Systems and Technology Laboratory (LSTS) of the Faculty of Engineering in the University of Porto in Portugal, from April of 2017 until July of the same year. During this visit, several path following methods and algorithms have been researched and a path tracking method has been developed in order to reduce the error in flight for non-holonomic vehicles following straight line trajectories. This method has been tested in simulations during the visit, while actual experiments have been conducted after its completion. This study has been published in the Workshop on Research, Education and Development of Unmanned Aerial Systems (RED-UAS) of the same year, as noted in the publications list. In the next section, a list of the participation in training events of the MarineUAS project is presented. 1.6 Participation in training and other dissemination activities Workshops
10 Chapter 1. Introduction • UAS Technology Workshop (2016, February), Andalucian Foundation for Aerospace Technology - Advanced Centre for Aerospace Technologies (FADA-CATEC), Villacarillo, Spain. • UAS Operator Training (2016, April), Northern Research Institute (NORUT), Tromsø, Norway. • Marine and Coastal Surveillance Workshop (2016, April), Maritime Robotics, Trondheim, Norway • Entrepreneurship workshop (2016, April), Maritime Robotics, Trondheim, Norway • Ethics workshop (2016, April), Norwegian University of Science and Technology (NTNU), Trondheim, Norway • Marine and Coastal Science Workshop (2016, June), Instituto Superior Tecnico (IST, Lisboa) and Instituto do Mar (IMAR, Azores), Azores, Portugal • Workshop on UAS in non-segregated airspace (2017, January), Honeywell, Brno, Czech Republic • MarineUAS First summer school (2017, June), Instituto Superior Tecnico (IST, Lisboa), Lisbon, Portugal • Presentation skills workshop (2017, June), Instituto Superior Tecnico (IST, Lisboa), Lisbon, Portugal Poster presentations • Balampanis, F., Maza, I., & Ollero, A. (2016, April). Area Decomposition and Path Planning in Coastal Regions, Remote Controlled and Autonomous Measurement Platforms Flagship (ReCAMP) Workshop, Northern Research Institute Tromsø AS (NORUT), Tromsø, Norway. • Balampanis, F., Maza, I., & Ollero, A. (2017, January). Distributed approaches for coverage and tracking missions with multiple heterogeneous UAVs for coastal areas, MarineUAS Mid-Term Internal Research Review, Honeywell, Brno, Czech Republic. • Rodríguez, L., Balampanis, F., Cobano, J.A., Maza, I., & Ollero, A. (2016, May) MarineUAS Innovative Training Network on Autonomous Unmanned Aerial Systems for Marine and Costal Monitoring, Spanish RPAS Association Conference UNVEX 2016, Madrid, Spain.
2 Decomposition of complex coastal areas Let me see you stripped down to the bone Depeche Mode Tw o main steps are required for a heterogeneous team of robots to be able to perform a coverage task in a given area: (i) a fair partitioning of the area according to the robotic capabilities and (ii) a feasible path planning method. In this Chapter, we present a decomposition method of the area in order to reduce the complexity of those steps and facilitate an online implementation of the solution. In particular, the identified properties of coastal areas for aerial vehicles show that an exact decomposition is needed to ensure complete coverage. Furthermore, a proper decomposition schema is a crucial consideration in order to mitigate the risks of flying unmanned robots near populated areas. Finally, by decreasing the complexity of the area by the means of a CDT, the partition and path planning problems can be reduced to graph search problems. 2.1 Introduction Area decomposition in robotics for coverage, is a widely used method which is performed by a discretisation of a previously known area. Decomposing an area into a collection of sensor-sized cells and treating it like a complete, sensor specific graph, some common tasks 11
12 Chapter 2. Decomposition of complex coastal areas are: obtaining a complete coverage path, increase the information gain and decrease mission time and energy needs. Moreover this decomposition reduces the overall computational complexity by simplifying the area and thus , it can be used to assign the task to a number of robotic vehicles. Then, by formulating a path formation strategy in the decomposed parts, various algorithms are used for the aforementioned tasks. In literature, several approaches can be identified for area decomposition, with Grid based methods and Cellular Decomposition be the most common. As will be described in the following sections, several of these approaches do not account for complex shapes or they simplify the problem, by not considering the heterogeneity of the robotic vehicles and their on-board sensors. Regarding coastal areas, the complex geographical attributes, the vast populated coasts and the increasing interest for activities in even remote off-shore locations, open questions regarding marine environment protection and sustainable management. The attributes of these areas, as well as the increased complexity of airspace management raise challenges in the integration of autonomous or remotely piloted aerial vehicles. More specifically, the integration of UASs into future airspace is one of the greatest challenges in ATM of the Single European Sky ATM Research (SESAR) 1 and the National Aeronautics and Space Administration (NASA) [ 1 ]. The use of UAS to cover wide areas implies the consideration of airspace restrictions and static and dynamic obstacle avoidance. This results in complex shapes that need to be decomposed adequately to ensure coverage by using the variety of heterogeneous UAVs in use today. Having in mind the aforementioned properties of coastal areas, the method used in this Thesis is a CDT [ 2 ]. This method manages to decompose complex areas by having nonconvex polygons with holes as an input and performs an exact decomposition regardless of the area complexity. In that manner an area is decomposed in a grid-like manner, having each of the triangular cells as an on-board sensor Field of View (FoV) sized triangle. Treating the resulting grid as a graph, the coverage problem can be reduced in a graph search problem, where every edge must be visited at least once. These algorithms are described in detail in the next Chapters. 2.2 Related work Most current decomposition approaches reduce or simplify the complexity of the area to be visited. In that manner, the vehicle has discrete states that can be reached by choosing different actions. Different decisions are required for tasks such as coverage or tracking, 1http://www.sesarju.eu/
2.2 Related work 13 (a) (b) Figure 2.1 A problematic found in grid decomposition methods. In (a) , Fig.6 from [ 5 ] is presented, where in order to achieve complete coverage, the path passes over cells which are partially not in the region of interest, like cells 4, 22, 32 etc. In (b) , a case as adopted from [ 4 ]. Big regions are considered as part of the obstacles, even though they are not. and cost-based path planning algorithms can be applied to reach each state [ 3 ]. Numerous studies have extended these algorithms to multiple robots, swarms or distributed systems. Regarding the decomposition process, the following approaches can be identified in literature: Grid based methods, Cellular decomposition and Boustrophedon (or Morse) decomposition. In grid based methods, which are also referred as “resolution complete methods” [ 4 ], the area partitioning is performed by applying a grid overlay on top of the area. In that way, there is a discrete simplicial complex where if all the cells are visited, then complete coverage is assumed. The decomposition problem is somehow bypassed, as for areas that are partially covered by the cells, the vehicle might not cover them at all in order to avoid collisions, or will go over them, ignoring the fact that they are obstacles or no-fly zones [5, 6] (see Fig. 2.1). In [ 7 ] even though triangles are used for the decomposition of the area, the same problem is still present. By changing the shape of the cells, while arranging them regardless of the area to be decomposed, the problem of having cells which are partially outside the domain of interest still remains. In general, all grid based methods tend to “pixelate” the area, a strategy that does not properly solve real world scenarios; areas do not have rectangular borders, neither do their obstacles and no-fly zones. Galceran et al. states in [ 4 ] that this problem can be solved by increasing, locally or globally, the resolution of the grid, creating smaller cells. This actually does not solve the problem in our case, since the dynamics of a fixed wing UAV have to be considered. Increasing the resolution does not guarantee a successful passage over the smaller cells since this operation implies sharper turns in a smaller area. Moreover, we have not identified in literature an algorithm that manages to
14 Chapter 2. Decomposition of complex coastal areas visit all cells of varying resolution and to avoid repeated coverage. This Thesis considers the size of the projected FoV of the UAV sensor for the decomposition method, where the resulting graph is also part of the path planning algorithms, as will be described in the following sections and chapters. As noted before, an exact decomposition in cells which the vehicle can visit, increases computational efficiency and ensures full coverage. If the cells are to be smaller, there is no guarantee that a feasible path for complete coverage can be produced. In cellular decomposition, a non-convex polygonal area is decomposed into a collection of convex polygons. Even though it manages to decompose complex areas, it does not always consider obstacles [ 8 ] within the area or it assumes that the resulting cells will be covered by simple and coherent back and forth lawnmower motions. In Boustrophedon decomposition [ 9 ] the problem of only producing many convex polygons without holes has already been addressed. In this strategy and its derivatives, critical points in the area are found and used in order to divide the area in cells. Even though this method is widely used, it still does not address the problem of covering the whole area by respecting the boundaries, as simple boustrophedon movements of the vehicles might not be able to cover complex spaces [ 10 ], or the resulting path planning algorithm produces repeated coverage paths [ 11 ]. Furthermore, some studies consider movements outside the region of interest in order to take turns and thus provide a full coverage pattern (see Fig.2 in [ 8 ] and Fig.2, Fig.5 in [ 12 ]). In [ 13 ], the authors use the aforementioned methods in order to perform complete coverage of an area. By using a variation of Morse functions, the area is segregated by performing a sweeping motion. In that manner, a set of trapezoidal cells is created and treated as a flow network. Nevertheless the results show that more complex areas and obstacles might not be treated properly by the coverage paths produced. In almost all of the above cases, it can be noticed that the surrounding area is rectangular and fairly simple. Moreover, most consider a coverage by lawnmower or boustrophedon paths, assuming that the vehicles will always be able to compensate uncertainties in motion or external factors like wind drifts for aerial vehicles, or surface anomalies for ground vehicles. As it can be seen in Fig. 2.2, coastal areas are usually far from being simple, convex polygons. In [ 14 ], the authors provide an optimal decomposition and path planning solution using an Integer Linear Programming (ILP) solver, by taking into account camera-sized grid decomposition. Their solution manages to obtain the desired optical samples, although they do not provide the computational time needed for the ILP to find a solution. In the same context, the authors in [ 8 ] use an enhanced exact cellular decomposition method for an area and provide a coverage path consistent with the on-board camera of the UAV. Although
2.2 Related work 15 Figure 2.2 A coastal area complexity example, Salamina island and Piraeus port in Greece. Numerous residential and commercial areas in a 12Km2 area, define narrow passages and complex shore shapes. their solution manages to produce smooth paths with minimal turns, their algorithms are tested only over convex polygon areas. The work presented in [ 15 ] deals with the same problem of area coverage for photogrammetric sensing. The authors include energy, speed and image resolution constraints in their proposed algorithms, such as an energy fail-safe mechanism for the safe return to the landing point. However, the provided solution and experiments do not account for complex, non-convex polygonal areas. Literature also provides a number of bio-inspired approaches that manage to apply strategies that have been observed in nature. In these studies like [ 16 ], the behaviour of insects leaving pheromones on their path is used to draw the boundaries of the covered area. In [ 17 ] the authors provide a particle filter strategy for a team of robots using bio-inspired swarm tactics in order to communicate the shape of the region of interest or anomaly as they name it in their study. In both cases, the areas are simple environments and a complete coverage is not guaranteed. The work presented in [ 11 ] decomposes an area by using convex decomposition and produces parallel lines in the decomposed parts, which are used as straight line paths. The algorithm presented tries to minimize the total amount of turns and provides a complete coverage plan. This strategy produces promising results, but complete coverage is not always achieved. Moreover, in some cases, repeated coverage is performed in order for the vehicle to visit the initial position of the next decomposed
22 Chapter 2. Decomposition of complex coastal areas Figure 2.6 The coordinate frames of the system. In presence of a gimbal which compensates for the movements of the UAS, the sensor always points down, parallel to the ground. 2.4.1 FoV based decomposition Let us consider that function P(t) calculates the area projection of the on-board sensor Field of View (FoV) to the ground at a specific moment t during the flight. This function is dependent on the relative rotation matrices among the UAS frame F U , sensor frame FS and ground FGcoordinate frames, as it can be seen in Fig. 2.6. The usual convention in aviation is used for the UAV movements; roll ( γ ), pitch ( β ) and yaw ( α ). The UAV coordinate frame has its x -axis pointing forwards in relation with the movement, the y -axis is given by the right-hand rule, while the z -axis points downwards. There are also cases where a gimbal is used in order for the UAV attitude to become irrelevant and the gimbal to compensate for the UAV movements, or its specific orientation to translate the sensing angle. Simple sensor scenario In a simple scenario where a UAV is equipped with a camera, mounted in a gimbal in order to compensate for the movements of the platform and facing downwards in a zero angle 0olike in Fig. 2.6, the FoV projection will be similar to the one shown in Fig. 2.7.
2.4 A Constrained Delaunay Triangulation method 23 As such, the FoV angles of θx and θy , draw the projected rectangle P(xp,yp) in relation with the flight altitude h: xp=2htanθx 2(2.2) and yp=2htanθy 2(2.3) Then, the largest inscribed circle CP of P is the one that has a radius r equal to the smallest side of P:Pmin . Hence, the largest triangle W that can be inscribed in this circle is the equilateral triangle WPand its side Wdis Wd=√3Pmin 2(2.4) Thus, as it will be described in the next section, the desired triangulation for the mesh generation is a CDT with maximum side equal to Figure 2.7 Calculating the CDT side constraint in an simple scenario. In the existence of a gimbal, the FoV projection is in a zero angle 0o , parallel to the ground. Thus, the maximum triangle that can be inscribed, is the inscribed equilateral triangle of the trapezoid’s inscribed circle.
24 Chapter 2. Decomposition of complex coastal areas Wd= 2√3htanθmin 2 2(2.5) If the camera takes a picture in a d ω interval, which depends on its frame rate, in the area constrained by the previously mentioned triangle, then by the time Ω it will have finished or stopped its task and the covered area will be A=ZΩ 0 P(ω)dω=ZΩ 0"√3 4W2 d#dω=ZΩ 0 2√3htanθmin 2 2 dω(2.6) in relation with the FoV smaller angle and altitude. This coverage metric only measures the area covered by the triangles and not by the UAS or the complete FoV projection. More complex scenarios In case the sensor is tilted in one or two axes, with the yaw axis z=0 and the gimbal compensates for the platform movements in order to keep the sensor stable, then the resulting FoV projection shape is either a trapezoid or an oblique projection. In the first case, where the sensor is only tilted in the YS axis in relation with the ground as seen in Fig. 2.8, a pitch angle β of the platform along with the bisector of the FoV angle θ, compose the front angle ψ, as it has been adopted by [8]. The ground FoV projection forms a trapezoid T=a,b,c,d , having the bottom side a equal to a=2htan(θy 2) sin(ψ−β+θx 2)(2.7) and upper side bequal to b=2htan(θy 2) sin(ψ+β+θx 2)(2.8) where θy and θx are the FoV angles of the camera on the YS and XS axe respectively, of Fig. 2.6. Then, the largest triangle W that can be inscribed in the trapezoid circle is the equilateral triangle WPand its side Wdis
2.4 A Constrained Delaunay Triangulation method 25 Figure 2.8 Calculating the CDT side constraint in case the sensor is tilted in one or more axis [ 8 ]. The angle θ is the FoV angle as in the previous examples. Pitch angle β along with the bisector of θ compose the front mounted angle ψ . Then, the ground projection is a trapezoid and the maximum triangle that can be inscribed, is the inscribed equilateral triangle of the trapezoid’s inscribed circle. Wd=√3ab 2=r32htan(θy 2) sin(ψ−β+θx 2) 2htan(θy 2) sin(ψ+β+θx 2) 2(2.9) Once again, as it will be described in the next section, the desired triangulation for mesh generation is a CDT with maximum side equal to 2.9. The calculation of the coverage is done in the same way as in Eq.( 2.6), replacing the size of Wd. Even more complicated scenarios are tackled by the authors in [ 29 ] and [ 10 ], while [ 8 ] provides a detailed calculation method. In this Thesis, we use a downward facing, zero angle, parallel to the ground scenario, where a gimbal compensates for the UAV movements. 2.4.2 Sensor specific - orientation agnostic segregation cells The aforementioned calculations manage to provide complete coverage formalisation in the simplicial complex. Having specific information about the sensor, namely its FoV and its sampling rate, the problem can be reduced to a discrete sampling scenario, by creating a grid where every cell is a sample. By performing a discretisation of the area according to the sensing capabilities and characteristics of the sensor, the CDT is applied by using a maximum triangle side constraint which can be calculated from the FoV of the
26 Chapter 2. Decomposition of complex coastal areas sensor as described in the previous section. Since the calculated triangle is inscribed in the projection’s inscribed circle, if every centre of each circle is considered as a waypoint, the configuration manages to produce segregation cells that are always covered by the sensor FoV if all waypoints are visited, no matter the flight direction; thus, this strategy is identified as orientation agnostic. Then, the integral function in Eq.( 2.6) for coverage calculation in continuous space can be simplified as the sum of m cells, for every segregated cell vk of the triangulation; since we take a sample in every discretised cell, the sum of all cells provides a complete coverage metric ∑m k=1P(vk), where the area of each cell vkcan be calculated as vk=ps(s−va)(s−vb)(s−vc)(2.10) having s as the semiperimeter of the triangle and a , b , c as the sides; hence the decomposition cells are platform and sensor specific. Please note that this calculation rule also applies in case the resulting FoV projection is not a square or a trapezoid, but a tangential quadrilateral. It is always possible to obtain the inscribed circle and the respective inscribed triangle as a function of the FoV sides. In case the FoV projection is a square, the inscribed triangle calculation is trivial. Considered model for the on-board sensors Regarding the use of on-board sensors, the literature mostly refers to cameras [ 30 ]. In those cases, by knowing the length, width and focal length of the camera, as well as the altitude from the sea, the shape of the projection of the FoV of the camera can be calculated based on the attitude of the UAV. This is not the case for point or side scan beam sensors, which have a wide width, but a really narrow length scanning profile. As explained before, in this Thesis the term FoV is used in order to refer to the projection of the FoV of a generic camera on the sea. As such, the on-board sensor along with a gimbal, point downwards at a zero angle and parallel to the ground. Thus, pitch and roll angles of the UAV with respect to the horizontal plane are not relevant to the FoV projection, since the gimbal compensates for these angles, as in Fig. 2.6. This is also true in other scenarios which involve communication sensors. Interesting oceanographic data is continuously collected and stored by buoys. Since a continuous operation of the communication devices of the buoys would waste considerable amounts of energy, a more conservative method is applied. Typically, data traffic in such scenarios is sparse and it is only needed for data acquisition or extraction. Hence, an interrogation method is applied where the communication circuitry is enabled either in intervals or operating in low energy mode until it is waken-up by a required signal. The latter should
2.4 A Constrained Delaunay Triangulation method 27 (a) (b) Figure 2.9 In ( a ), the data extraction structure scenario where in ( b ), the communication range cone area RA of the buoy on a side (left) and a top (right) view. The top view of this range at a specific altitude which makes the transmission feasible, creates a circle in which an inscribed triangle can be used as the maximum resolution size and constraint of a triangulated mesh. Then, each centroid of this collection of triangles, serves as a waypoint for the coverage task. have an appropriate frequency, while the communication range of the wake-up sensor is limited. Figure 2.9a shows the on-range area RA to enable transmission that is determined by the transceiver specifications. Since the location of the buoys is not always known due to sea drifts, the proposed decomposition method will perform the aforementioned CDT and each triangular cell will have the size of the expected buoy communication range cone, as can be seen in Fig. 2.9b. Sample rate is another sensor characteristic to consider for the cell decomposition of a region. It is necessary to decompose the area in a manner that a sensor can obtain at least one sample of each of the resulting cells in a unit of time (see Figure 2.10). Then, the projected footprint area P(t) must guarantee that its size is proportional to the sample rate T and the UAV speed V . This can be achieved by either reducing speed or increasing altitude in order to grow the projection of the FoV. For most of the sensors, the sample rate is not an issue; however, this aspect has been pointed out for the sake of generality. 2.4.3 Lloyd optimisation grid enhancement The CDT decomposition strategy discussed in the previous sections, does not guarantee that all the segregation cells are going to be equilateral triangles, thus providing smaller triangles and larger sample overlapping. Its constraints are limited to either maximum
28 Chapter 2. Decomposition of complex coastal areas (a) (b) (c) (d) Figure 2.10 The appropriate cell decomposition is proportional to the velocity V , the sample rate T and the FoV projection footprint P . In ( a ), in every time step t0,t1,...,tn , F is not large enough for the sensor to take a complete sample, whereas in ( b ), T is not fast enough to obtain a sample from each area. In both cases, the problem could be solved be either reducing the speed, while increasing the sample rate if possible or increasing the altitude to increase the projection of the FoV. ( c ) shows the ideal solution in the limit, whereas ( d ) depicts the most usual case of the same portion of the sea being present in many samples. triangle side or minimum angle and the area itself. Hence the choice of optimising the resulting grid is merely a task and equipment defined choice: if a varying overlapping schema is not critical, then the initial grid can be maintained. On the other hand, if a normalised schema is preferred, a series of Lloyd optimization [ 28 ] iterations can be performed in the resulting grid. This algorithm improves the homogeneity of the mesh, by altering the angles of each cell, making them as close as possible to 60 degrees. The method is dependent on the iterations selected and since it moves all the points H of the CDT(H) triangulation in each iteration, it is computationally expensive [ 31 ]. Therefore, it is recommended mainly for offline computation; more information on the actual complexity of the Lloyd optimization can be found in [ 32 ]. Nevertheless, having the size of each cell of the triangulation as well as the FoV of the sensor, the magnitude of sample overlapping can be calculated in the initial coverage planning. It is clear that the more regular the triangles are, the less overlapping will occur. Since by increasing the Lloyd iterations the distribution of triangles closer to equilateral is also increased, this operation can be added to the path generation cost functions. Finally, there is another positive effect of the Lloyd
2.5 Results 29 Figure 2.11 A CDT of a non-convex polygon with an internal hole. The simulated shape is a 5km x 5km area, scaled down by a 1/10 factor for visualization and notation purposes . optimisation on the resulting grid, explained in detail in the following Chapters: since the overall angle is tending to reach 60o , the paths created by connecting every centroid of the triangular mesh have smoother angles, thus resulting in smoother turns; easier to be followed by non-holonomic, fixed wing UAVs. 2.5 Results In every case, the CGAL implementation of the triangulation has managed to produce a mesh consisting of homogeneous triangles. The triangulation and meshing functions were tested by giving several input constrain criteria for the lesser angle bound constraint of the produced triangles, as described in Section 2.4.1. The proposed method has been tested on various scenarios in order to obtain qualitative and quantitative results, iterating through angle bound constraints for the produced cells. Resulting meshes have been examined, obtaining information about triangle homogeneity. Initially, a relatively simple non-convex polygon with an internal no fly zone was used, as shown in Fig. 2.11. Produced vertices, total cells, minimum resulting edge and minimum resulting angle are presented in Table 2.1. In all trials, the input criteria for the edge size have remained the same; a predefined size for the FoV, as defined in Fig. 2.6, of 25 metres (250 metres in the actual area). It is noted that reaching the upper angle constraint limits of the Delaunay
30 Chapter 2. Decomposition of complex coastal areas Table 2.1 CDT performance based on lower angle bound constraint. Angle bound (deg) Vertices Cells Min edge length (m) Min angle (deg) 0.011-8.88 304 476 10 25 9.16-10.02 305 478 10 25 12.31 313 488 10 29 15.29 321 500 9.43 31 16.15 325 508 8.24 32 17.01 357 568 6.32 32 18.16 377 597 6.32 34 Table 2.2 Square grid and Constrained Delaunay Triangulation decomposition methods comparison. Footprint Size (m) Total cells (square mesh) Total cells (CDT) Partial coverage by square mesh 15 331 1395 125 20 211 784 99 25 138 476 73 30 104 352 64 triangulation, an increased number of sub-optimal cells were produced near the no fly zone. Moreover, the algorithm could not proceed further than the angle bound constraint of 18.9 degrees, a size which is consistent to the one being referred to in CGAL notes. Nevertheless, the FoV criteria are always met, as in all produced triangles at least one side has the edge size input constraint. Comparison with a square grid method The same shape was decomposed following a simple square grid decomposition strategy as in the studies presented in [ 4 ]. A footprint size of 15, 20, 25 and 30 metres has been tested in order to observe the amount of cells that partially cover the region of interest by that method. Note that using the strategy presented herein, there are no cells that partially cover an area, and the centroids of the triangular cells produced are always inside the domain. Results in Fig. 2.12 show that even though this method manages to produce a lesser amount of cells, its performance is highly dependent on the shape of the area, the restricted zones, as well as the resolution of the FoV. Comparison of the two methods in the 25 metres sized footprint as presented in Table 2.2, shows that half of the produced cells by the square grid method share a region outside the domain of interest, raising questions on the restrictions that have to be applied in order to either visit or not those specific cells. The same comparison has been performed in a more realistic scenario in an actual area,
2.5 Results 31 Figure 2.12 Grid decomposition performance of the area in Fig 2.11. In the 25 metres footprint case, out of the 138 total cells that cover the area, 73 cells are partially inside the domain of interest. Figure 2.13 A coastal area in Schinousa island on the right. On the left, the same area when a grid overlay is applied. Every orange shaded cell represents an ambiguous cell. as shown in Fig. 2.13. In a total of 218 square cells, 78 had partially a no-fly zone in them. As we will see in the next Chapters, these ambiguous areas raise critical questions in path generation. These preliminary results demonstrate the feasibility of the method which provides an initial tool for exact area decomposition for the next stages of area partition and the extraction of waypoint plans.
38 Chapter 3. Area partitioning for multiple vehicles Algorithm 1: Antagonizing wavefront propagation algorithm that computes the baseline area partition. Q is a queue list managed as a FIFO by functions insert and getFirst 1vIk: initial vertices/triangular cells for each UAV Uk, 2Sv: area size of triangular cells v, 3N(v): the set of neighbours of vertex v, 4A(v): the UAV identifier allocated to triangular cell v, 5Zk: area coverage capability of UAV Ukin square meters, 6SvMin: area size of the smallest triangular cell in CDT, 7C(vi): is the unitary transition cost of vertex vi, 8Cp: is a unitary transition cost place-holder 9Cp←1 10 foreach Uk∈CDT do 11 Q.insert(vIk) and mark vIk as visited; 12 Zk←Zk−SvIk ; 13 C(vIk)←Cp; 14 end 15 while Q not empty do 16 v←Q.getFirst(); 17 if C(v)6=Cpthen 18 Cp←Cp+1; 19 end 20 foreach vi∈N(v)do 21 k←A(v); 22 if vinot visited AND Zk>SvMin then 23 C(vi)←Cp+1; 24 A(vi)←A(v); 25 mark vias visited; 26 Q.insert(vi); 27 Zk←Zk−Svi; 28 end 29 end 30 end in a breadth-first manner [ 3 ], the strategy is called Antagonizing Wavefront Propagation (AWP). This strategy is presented in Algorithm 1 and works as follows. Let us consider each of the initial positions of the UAVs as the root node of a tree; each root node is given an initial step cost of one. In every recursion step, each vertex that has an edge connected to the parent vertex is given that cost plus one. In addition, vertex vi gets the same A(vj) attribute of its parent vertex vj , propagating the identifier of the UAV in that way. In case
3.3 Deadlock handling 39 the number of vertices for the Uk UAV meets its autonomy limit, denoted by the total area Zk it should cover, the algorithm for that UAV stops. Please note that these steps are not performed if a triangular cell already has any of these attributes. Since queue Q is accessed once for every i -th cell, the while iteration has a complexity of O(m) , where m is the number of vertices. The complexity of getting the first element is O(1) ; then, it inserts new elements according to the restrictions. The insertion in a stack has a complexity of O(1) . Hence, the complexity of Algorithm 1 is O(m) . The area partition computed is not sufficient for complex cases where deadlock situation occur after applying Algorithm 1. A further adjustment step is needed in order to assign regions where a deadlock has happened, as will be described in Section 3.3. 3.2.2 Reverse Watershed Schema By performing the previous step, each configuration space is either adjacent to another configuration space or to the borders of the whole region. In that manner, a second algorithm (see Algorithm 2) assigns to each vertex that already has an UAV identifier a unitary border-to-center cost attribute D(vi) of proximity from the borders to the center of the configuration space. The triangular cells that are adjacent to a border with another configuration space or to the whole area are given a high D(vi) cost and are considered as the root nodes of a tree. In each step of the algorithm, this cost is decreased and propagated to the adjacent triangular cells of these nodes. This function manages to create a border-to-center pattern resembling a watershed algorithm, and as such it is called the Reverse Watershed Schema (RWS) algorithm. Please note that this cost attribution will further facilitate the coverage plans to be produced, as will be described in the following Chapter. Finally, we have to note that the complexity of this algorithm is similar to the previous one. The first loop has a complexity of O(m) , O(m) for the initial foreach loop and O(1) for each insertion to the stack, since the second inner foreach has a maximum of three iterations. For the same reasons, the second while loop has also a O(m)complexity. 3.3 Deadlock handling The previous baseline partitioning algorithm is able to perform well in most cases where the area is simple or where the initial positions of the UAVs are evenly distributed in the area. Nevertheless, it may lead to several deadlock scenarios as the growing sub-areas meet each other, as can be seen in the example shown in Figure 3.2. Hence, a Deadlock Handling (DLH) algorithm that adjusts the initial partitioning in the non-allocated areas
40 Chapter 3. Area partitioning for multiple vehicles Algorithm 2: RWS algorithm for the generation of the border-to-center cost D(vi) attribute. Q is a queue list managed as an FIFO by functions insert and getFirst. 1N(v): the set of neighbors of vertex/triangular cell v, 2A(v): the UAV identifier for triangular cell v 3foreach v∈CDT do 4foreach vi∈N(v)do 5if A(v)6=A(vi)then 6mark vas visited; 7D(v)←∞; 8Q.insert(v); 9end 10 end 11 end 12 while Q not empty do 13 v←Q.getFirst(); 14 foreach vi∈N(v)do 15 if vinot visited then 16 D(vi)←D(v)−1; 17 mark vias visited; 18 Q.insert(vi); 19 end 20 end 21 end is needed, by exchanging UAV identifiers or assigning UAVs to the empty areas. Two different approaches have been tested, by applying the two algorithms presented before. As was stated before, each UAV Uk should cover an area of Zk . In a test case, after the initial partitioning, let us consider that Yk6=Zk space has been allocated to Uk . In the deadlock scenarios (see Figure 3.2), there are areas that do not belong to any UAV. These areas are allocated to a virtual UAV U−1 with area Z−1=0 . Let us consider a list LU that contains the results of Zk−Yk for each UAV and the area size of the smallest triangle in the CDT SvMin . Thus, each UAV can have an area surplus if Zk−Yk>SvMin or a shortfall if Zk−Yk<SvMin . The latter case always happens in the deadlock scenarios for U−1 . In each recursion of Algorithm 3, a pair of UAVs, for instance Ui and Uj , one having an area surplus and another with a shortfall, is chosen from the list in order to gradually exchange triangular cells between them to reach the desired area size, one area at a time. In order to do so, a feasible transition sequence Pi j must be found, as shown in Figure 3.3. The complexity of this algorithm is calculated as O(nm2) due to the findSequence function, which is actually a tree sort; n is the number of UAVs and m the number of
3.3 Deadlock handling 41 Figure 3.2 A deadlock scenario. Four UAVs U1 , U2 , U3 and U4 after the baseline partition Algorithm 1. UAVs U1 , U2 and U4 have met their autonomy capability of Zk by covering Yk area. Nevertheless, U3 was not able to overtake any more area, being “blocked” by the other UAVs and the borders of the whole region. Different colors/shades of areas indicate the configuration space of each UAV, while the numbers inside the cells indicate the isotropic cost, as has been assigned by Algorithm 1. The free or non-allocated areas belong to virtual UAV U−1. cells. The complexity of the move functions is displayed below, in each of the following transposition algorithms. Two algorithms called moveAWP and moveRWS have been implemented for the move function, which is used in Algorithm 3. In the former, the farthest vertex in the Ri sub-area of Ui UAV in sequence Pi j is chosen, by using the unitary transition cost information from Algorithm 1. This vertex has also to be adjacent to the second area in the transition sequence Pi j . Starting from that vertex, Algorithm 1 is applied again, overtaking the requested area size in the means of exchanging UAV identifiers between those triangular cells. Recursively, this operation is performed for every item of the sequence (see Algorithm 4). The complexity is O(k) , where k is the number of areas that are in the transition sequence. Then, since Algorithm 1 is used for the transposition function, the whole complexity is O(l∗m2) , where l is the number of cells to be transposed and w is the findSequence algorithm complexity. In the second approach, we apply the RWS algorithm in order to get a depth schema of
42 Chapter 3. Area partitioning for multiple vehicles Algorithm 3: Multi-UAV partitioning Deadlock Handling (DLH) algorithm. Baseline partitioning is performed by Algorithm 1, whereas this method is for the subarea size adjustment (if needed). Function getSurplusUAV (L) gets a UAV identifier from list L that has an area surplus, whereas function getShort f allUAV(L) gets the identifier of a UAV that has an area shortfall after Algorithm 1. Function findSequence finds a feasible transition sequence Pi j between UAV Ui and Uj , whereas the move function performs the transfer between triangular cells 1SvMin: area size of the smallest triangular cells in CDT 2while ∃U∈LU<SvMin do 3i←getSurplusUAV(L); 4j←getShortfallUAV(L); 5Pi j = findSequence(Ui,Uj); 6if Zi−Yi>Zj−Yjthen 7move(Zj−Yj,Pi j); 8end 9else 10 move(Zi−Yi,Pi j); 11 end 12 end Algorithm 4: MoveAWP algorithm. C(v) is the transition cost from the AWP algorithm (see Algorithm 1). Then, function findBiggestC(v)(Ui,Pi j) finds the triangular cell with the largest transition cost value of UAV Ui that is adjacent to UAV Ui+1 in the sequence. Then, function Awp takes as variables an initial cell v , the area size that needs to be exchanged and the UAV identifier that needs to be exchanged from. The growing function is similar to Algorithm 1 1Pi j : the transition sequence between Ui and Uj for triangular cells exchange, treated as a list 2vinit : initial triangular cell for identifier exchange 3S: area size to be moved 4foreach Ui∈Pi j do 5vinit = findBiggestC(v)(Ui,Pi j); 6Awp(vinit , S, Ui); 7end the adjacent areas, as was described in Section 3.2.2. In each recursion of the algorithm (see Algorithm 5), the amount of triangular cells that are in the borders of the first pair of the transition sequence exchange their UAV identifiers. If the area of these border triangular cells sum up less than the requested area, then the area size and the total amount of border triangular cells exchange their UAV identifiers. If not, then only the triangular cells in the
3.3 Deadlock handling 43 Figure 3.3 Transition sequence selection. After the initial partition process, U1 , U2 and U4 have met their sub-area size constraint and blocked the growth of U3 . As a result, three areas are not allocated ( U−1 ). The feasible transition sequences A (green(|)), B (red(||)) and C (blue(|||)) are used in order for U3 to obtain the requested total area, by gradually exchanging cells in every pair of the sequence. Sequence D (black) does not lead to a partition that has an area shortfall, and thus, it is a not feasible sequence. front (in the borders) exchange their UAV identifiers. This amount of triangular cells is then exchanged to the next UAV in the sequence and so on, maintaining the aforementioned restriction, until all of the requested area and associated triangular cells are transposed from the initial UAV in the sequence to the last. This algorithm’s complexity is O(nw2) due to the use of the findSequence function, as has been described in Algorithm 3. There are two main differences in these approaches. In the first approach, we have a wavefront pattern from a single triangular cell, whereas in the second, the exchange progresses as a kind of width sweep Morse function [4]. The second difference is that in the first approach, all of the triangular cells to be exchanged are the transposed UAV first, and in the second approach, only the amount of triangular cells that are in the adjacent
44 Chapter 3. Area partitioning for multiple vehicles Algorithm 5: MoveRWS algorithm. Function findSequence finds a valid transition sequence, as can be seen in Figure 3.3. This function is also called before the initial recursion of the MoveRWS algorithm. Function exchangeIdenti fiers makes use of the information of the RWS algorithm (see Algorithm 2), and it exchanges robot identifiers on two adjacent configuration spaces, by exchanging the amount of triangular cells that have the lowest coverage cost, but are adjacent. It also propagates and extends this cost. Function restO f Sequence returns the remaining sequence for the specific P[i]→P[i+1] transition, in order to initially transfer only the amount of triangular cells that are adjacent between i and i+1 until the final Uj UAV. In case this happens, the requested area has not been exchanged yet, so the algorithm runs recursively, and the last line takes a step back in sequence traversal 1S: area size to be moved 2Sad j(kl): the area size of adjacent triangular cells between UAV k and l 3Pi j the transition sequence between Uiand Ujfor triangular cell exchange, treated as a list 4foreach Ui∈Pi j do 5Pi j ←findSequence(Ui,Uj); 6if Sad j(P[i],P[i+1]) >Sthen 7exchangeIdentifiers(P[i], P[i+1], S); 8end 9else 10 Prest = restOfSequence(Ui); 11 MoveRWS(Sad j(P[i],P[i+1]),Prest ); 12 S=S−Sad j(P[i],P[i+1]); 13 i=i−1; 14 end 15 end borders are transposed in each step. In that manner, the triangular cells of the area are propagated respecting the total amount of cells that each UAV has each time, resolving overlapping UAV issues, as will be discussed in Section 3.4. 3.4 Results Two coastal areas in Greece have been selected for the experiments (Figure 3.4). The first is a broad and populated shore near the harbor of Piraeus, Salamina, the same as presented in Figures 2.2 and 2.5. The second is a remote island in the Aegean archipelago, Astipalea. The first area was used for evaluating and comparing the partitioning algorithms, whereas the second for evaluating the proposed strategy in various setups.
3.4 Results 45 (a) (b) Figure 3.4 Selected areas for testing: ( a ) Salamina area having narrow passages and complex shapes in shores; ( b ) the Astipalea area is used for testing the suitability of the proposed algorithm. In the first area, two square shaped FoV constraints with side sizes of 100 meters and 900 meters have been used, as has been defined in Section 2.4.1, in order to compare the behavior of the algorithms. Three test cases of different relative UAV capabilities have been tested: one has an equal relative percentage for the UAVs of 33,3%, another has a distribution of 10%,60% and 30% whereas the last test case has a distribution of 80%, 10% and 10% for the respective UAS (see Fig. 3.5). The partitioning strategies called MoveAWP and MoveRWS described in Section 3.3 have been compared. The former uses the transition cost of the AWP algorithm, and the latter applies the RWS algorithm for adjusting the baseline partition computed by Algorithm 1. The complexity of the area has managed to highlight some issues that were not evident for the majority of simple areas. The main problem occurs during cell exchange when the initial position of the UAV is close to the borders, because a sub-area could overtake the initial position of the UAV (see Figure 3.5a). Additional simulations have been performed to measure the performance of the different algorithms with respect to the metrics F and G explained in Section 3.2. In particular, the simulation environment shown in the second area of Figure 3.4 has been used for the metric F . Some results are detailed in Figure 3.6 for three and five UAVs and a FoV size of 30 m. In general, simulations have been executed for three and five UAVs, with initial locations evenly or randomly distributed in the area and different FoV sizes. The results for the sum of distances between each center of the triangular cell inside a sub-area and
46 Chapter 3. Area partitioning for multiple vehicles (a) - MoveAWP 250 m FoV (b) - MoveRWS 250 m FoV (c) - MoveAWP 2 km FoV (d) - MoveRWS 2 km FoV (e) - MoveAWP 250 m FoV (f) - MoveRWS 250 m FoV Figure 3.5 Partitioned area for three UAVs (indicated by the white cells) and visualized by using the ROS rviz node. Each four cases represent the results for different relative capabilities: cases a-d is the 10% (red), 60% (blue), 30% (green) case; cases e-h depict the 33% (red), 33% (blue), 34% (green) case; whereas the cases show the 80% (red), 10% (blue), 10% (green) case. Each consecutive pair of images indicate the comparison of the two algorithms. Images ( a , b ) show how the MoveAWP and MoveRWS algorithms have performed with the small (250 m) FoV, whereas ( c , d ) show the results for the large (2 km) FoV case. Table 3.1 An even distribution of initial locations for three and five UAVs, with different relative capabilities. FoV (15 m) FoV (30 m) #UAVs moveRWS moveAWP moveRWS moveAWP Metric F(m) 3 333.861,84 333.909,67 82.768,44 84.979,76 Metric F(m) 5 437.988,74 439.642,85 129.879,24 131.516,96 the initial location of the UAV inside that sub-area (metric F ) are shown in Tables 3.1 and 3.2. In both cases, it can be seen that the moveRWS algorithm has a better performance than moveAWP, since the metric is lower. Regarding the other metric G considered in Section 3.2, simulations have been performed
3.4 Results 47 (g) - MoveAWP 2 km FoV (h) - MoveRWS 2 km FoV (i) - MoveAWP 250 m FoV (j) - MoveRWS 250 m FoV (k) - MoveAWP 2 km FoV (l) - MoveRWS 2 km FoV Figure 3.5 Cont. (a) (b) Figure 3.6 Area of Figure 3.4b selected for the comparison of the two partitioning algorithms. ( a ) is partitioned for three UAVs, whereas ( b ) for five UAVs. The depicted FoV size is 30 m in both cases. Both figures are computed with the deadlock moveRWS handling of Algorithm 5. Lighter areas show the borders of each sub-area.
54 Chapter 3. Area partitioning for multiple vehicles Table 3.5 Cont. UAV Capability % FoV (15 m) FoV (30 m) Lloyd Iterations Lloyd Iterations 20 30 60 20 30 60 20% 0,01/0,03 0,01/0,21 11,23/0,77 0/1,8 0,02/1,43 10,78/0,4 30% 0,53/0,02 12,85/0,24 0,02/0,68 0/0,27 13,45/0,42 0,02/1,4 20% 0,01/0,03 0,01/0,11 11,59/2,65 6,24/0,57 6,88/0,44 10,78/0,64 10% 0,01/0,07 0,01/0,19 3,382/1,19 0/0,46 0,04/0,21 3,97/0,04 20% 0,01/0,14 0,01/0,14 0,01/1,36 0/1,03 0,02/0,37 0,02/0,4 Average% 0,11/0,05 2,58/0,18 5,25/1,33 1,25/0,82 4,12/0,57 5,11/0,58 G/area(R)%0,57/0,27 12,89/0,89 26,24/6,64 6,24/4,12 20,59/2,87 25,57/2,88 10% 0,02/0,46 4,36/0,05 4,2/0,26 0,04/0,36 0,04/0,18 0,04/0,09 20% 0,02/1,02 4,23/0,26 0,01/0,06 0,05/0,3 0,05/0,36 0,05/1,49 10% 0,02/0,42 0,01/0,47 0,01/0,01 0,04/1,2 0/0,41 0,04/0,25 30% 19,04/0,5 0,01/1,22 0,01/0,01 20,32/0,74 20/1,48 19,85/0,31 20% 0,01/0,9 0,02/0,04 0,01/0,02 10,32/0,49 10,32/0,24 10,08/0,65 10% 0,02/0,66 0,01/0,58 0,01/0,2 0,04/0,028 0/1,01 0,04/0,22 Average% 3,18/0,66 1,44/0,43 0,71/0,09 5,13/0,52 5,07/0,61 5,02/0,5 G/area(R)%19,13/3,96 8,64/2,62 4,25/0,56 30,81/3,12 30,41/3,68 30,1/3,01
3.5 Conclusions 55 Figure 3.8 A graphical representation of Tables 3.4 and 3.5. Average difference of ideal allocation for every UAV versus final area allocation, after the baseline algorithm (orange) and after the deadlock moveRWS treatment algorithm (green). As expected, random initial positioning of UAVs creates more often deadlock scenarios for the baseline algorithm. The algorithm has been tested for 3–6 UAVs, evenly or randomly distributed in the area. FoV projections of 15 and 30 m have been tested, and in each case, a different Lloyd iteration setting (20, 30 and 60) has been set. The horizontal lines show the average difference. 3.5 Conclusions The proposed method shows that the produced mesh and associated graphs manage to be consistent with the area properties and the capabilities of the UAVs and their on-board sensors. The defined metrics and the results show that the approach manages to provide an adequate solution without sacrificing the precision which a continuous solution would achieve. The presented partitioning algorithms are capable of treating any UAV setup since they solve deadlock scenarios regardless of the initial positions or the complexity of the area. The latter, as will be described in the following Chapters, along with the attribute of a uniformly expanding partitioning schema, show that the proposed framework provides smooth complete coverage waypoint plans in a decentralized manner. The framework presented so far is actually a generic waypoint planner by treating every centroid of the triangular grid as a waypoint. However, it does not take into consideration the UAV platform dynamics regarding path tracking, even though the pitch and roll upper boundary is used calculate the maximum triangle cell. Nevertheless, a platform might, in the case
56 Chapter 3. Area partitioning for multiple vehicles of a holonomic vehicle like a multirotor, or might not, in the case of a non-holonomic fast moving fixed wing, be able to follow sharp turns that are produced. Hence, so far this solution does not account for waypoint to waypoint flight trajectories. These issues are usually addressed by the trajectory generator and flight controller, for example, by assuming that a waypoint has been visited if the UAV passes close by. This metric usually is task specific and, in real-world applications, user defined. However, an optimization of the coverage waypoint plans and an enhancement of a state of the art trajectory tracker will be presented in the following Chapters, thus leading to a complete top-to-bottom architecture.
4 Coverage plans I can go anywhere The Who Ha ving acquired the subregions for each UAV, we can once more utilise the gridas-graph approach used in the previous Chapters, in order to produce coverage waypoint flight plans. By treating each triangular cell centroid as a waypoint of a flight plan, complete coverage can be assumed in case the vehicles pass over all of these waypoints. In that manner, the vertices of the area graph should be prioritised according to some optimisation criteria, in order to define the order in which they are visited. There are several restrictions to be considered in this prioritization; initially the paths must have limited or no repeated coverage. Several studies show that repeated coverage cannot be avoided, since this issue resembles the Seven Bridges of Königsberg or the multiple Travelling Salesman (mTSP) problems. Moreover, since this study is focused mainly on non-holonomic vehicles, the resulting coverage plans should have paths as smooth as possible, in order to provide feasible trajectories for the vehicle to follow. As it is already stated, the latter is a critical consideration, since the main goal is to respect the Region Of Interest (ROI) of every UAV; as we note in the next Section, current research provides paths which pass over no-fly zones in order to proceed in the next boustrophedon parallel path. Finally, as also presented in the previous Chapters, area characteristics impose several restrictions on the design of such a flight plan. In this Chapter, a proper identification of these area properties as well as a novel algorithm is presented in order to provide smooth, complete coverage flight plans, with minimal repeated coverage. 57
58 Chapter 4. Coverage plans 4.1 Related work Two main coverage path planning approaches can be identified in literature: continuous and discrete. In both of them, usually a back and forth motion strategy is applied, commonly known as “lawnmower path”, “zigzag” or boustrophedon motions. Several derivatives exist, like the Zamboni flight pattern [ 45 ], but in general the simplicity of calculating these paths is preferred over more complex trajectories. There are also efforts for sensor-driven paths[ 10 ] in order for the robot to increase the information gain from its sensors, instead of trying to preplan a solution, proposing like so an online strategy for coverage. In the discrete, grid-like methods, every cell of a segmented space is considered as a node in a, usually, undirected graph[ 5 ]. The vertices between these nodes might or might not have weights, which describe a transition cost. In literature, cost or weights are applied depending on e.g. the distance from the target or from the information gain that might have [39]. Regarding repeated coverage of a known area by using a team of robots, the work presented in [ 46 ] tries to tackle the problem by evaluating the results by the metrics of total path lengths, total average visiting period, total worst visiting period and balance of workload distribution. In order to also tackle the art gallery problem of maximizing the visibility of the robots over the whole area, that algorithmic approach generates a number of points that serve as a kind of guards which maximize complete visibility. These points along with an application of a CDT in the whole area, produce a set of nodes which serve as a graph for creating coverage waypoint lists, by using cluster-based or cyclic coverage methods. Coverage with aerial vehicles is also a relevant topic of the previously mentioned survey of Galceran et al. [ 4 ], which includes several related papers. In many cases, grid decomposition strategies have been used [ 5 ], [ 47 ], whereas the enhanced exact cellular decomposition of Li et al. [ 8 ] manages to decompose an area in convex cells and by using boustrophedon motions taking to account the FoV of the vehicle, provide coverage trajectories with minimum turns. In the latter case of convex decomposition, Maza et al. in [48] also used a convex decomposition method to partition the area. However, complete coverage is not always achieved and usually, as we have shown in the previous chapters, simple non convex areas are chosen as a test case. In case of non convex areas, the usual strategy is decomposing that region into a sum of convex polygons. Unfortunately sometimes this leads to repeated coverage paths as seen in [ 11 ] when the vehicle traverses to the next polygon. Like so, spiral paths might give a better coverage in comparison with boustrophedon movements [ 8 ]. In this Thesis we chose a spiral path strategy which can be either outward, from the initial position of the aerial robot to the
4.2 Waypoint prioritization 59 borders of the area, or inward, starting from the closest border cell to the initial position. This is achieved by prioritizing each of the triangular cells centroids according to the current position of the vehicle, based on the weights which have been attributed in the decomposition phase. As we will show in the next Sections, several geometrical properties of the area impose specific restrictions for these paths, and as such a novel algorithm is proposed in order to treat those properly. 4.2 Waypoint prioritization Consider a list of waypoints W={w1,w2,...,wi} of a ROI as an undirected graph G= (V,E) , where the vertices V are the waypoints, and the edges E are the transitions between the waypoints. Also consider the resulting theoretical straight line path P , which connects each of the waypoints, by not visiting a waypoint twice. Let us define as θwi]wi−1wiwi+1 the angle between two consecutive edges defined by three consecutive waypoints of P . As discussed before and shown in several studies with aerial vehicles [ 49 ], the transitions and turns between the waypoints or between the lines of a lawnmower path, are outside the ROI. This is particularly true in fixed wing vehicles, where a constant forward motion cannot allow sharp maneuvers. Some control frameworks treat this issue by making a turn before the waypoint, heading towards the next one, or by passing outside the ROI [ 12 ] [ 10 ]. The proposed method of grid construction guarantees full coverage if the aerial robot passes exactly over the waypoints. Then, the coverage flight plan problem is defined by the following three constraints: • Construct a spiral path by using each of W just once, avoiding passing over a triangular cell twice. • Find a minimal P total path length and reduce the repeated coverage which might occur in complex areas. •Increase the average value of the θwiangles of P: θw= N−1 ∑ i=2 θwi/(N−2),(4.1) where Nis the number of waypoints. By combining these constraints, smoother trajectories are expected [50].
60 Chapter 4. Coverage plans 4.2.1 Spiral coverage plans The framework presented in the previous Chapter, and in particular the weight attribution schema of Algorithm 2, can be applied to generate waypoint lists for the UAV to achieve complete coverage of a complex coastal sub-area. By using the border-to-center cost described in Section 3.2.2, inward spiral-like waypoint lists W can be generated. Algorithm 6 performs a selection of vertices by initiating from the vertex that has the highest Reverse Watershed Schema (RWS) cost D(v) and is closer to the starting position of the UAV. In every recursion, the closest adjacent cell vj that has the same cost ( D(vj) = D(v) ) is inserted in the waypoint list. In case all of the vertices of the same Dv cost have already been inserted in the list, the algorithm reduces the cost Dv and chooses the cell that is closer to the previous step. It should be mentioned that the complexity of this algorithm is O(n2) ; as such, Figure 4.1 presents the relative time cost of decomposition, partition, Lloyd optimization and waypoint plan calculation of Table 3.3, of the previous Chapter results. Algorithm 6: Waypoint list computation for coverage. Dc is an auxiliary variable with the current border-to-center cost in each step, whereas vIk is the starting position of the UAV Uk . Function findClosest finds the closest vertex to the current one that has its same border-to-center cost. CDTk is the sub-CDT for UAV Uk.Wis the produced waypoint list 1Dc←∞; 2v←findClosest(vIk,Dc); 3W.insert(v); 4foreach v∈CDTkdo 5if ∃v,D(v) = Dcthen 6vj←findClosest(v,Dc); 7W.insert(vj); 8v←vj; 9end 10 else 11 Dc←Dc−1; 12 end 13 end In order to simulate how an actual coverage trajectory would be by the proposed algorithm, a scenario has been chosen in the same area of Figure 2.13. The area was partitioned into sub-areas respecting the different UAVs’ capabilities as described in the previous Chapters (see Fig. 4.2a) and their respective coverage paths connecting the waypoints can
4.2 Waypoint prioritization 61 Figure 4.1 A percentage bar chart depicting the average time needed for each of the operations of all the simulations of Table 3.3. This plot also includes the waypoint generation strategy. A comparison between a large FoV projection side size (800 metres) and small FoV side size (100 metres) is presented, while also comparing the difference between different Lloyd iteration settings. be seen in Fig. 4.2b. Finally, Fig. 4.2c shows the simulated trajectory for an UAV, in the simulation framework which will be described in Section 4.3. This method offers complete coverage in simple environments. Nevertheless, geometric properties in complex areas indicate that this method increases repeated coverage as it introduces back and forth motions in narrow passages, as described in the following section. 4.2.2 Spiral plans and area characteristics By using the aforementioned method, spiral paths create repeated coverage patterns in areas like ( C ) of Figure 4.3. We use the term valleys to refer to this kind of areas, while areas like ( A ) and ( B ) of the same figure are refered as inlets. As it is also described in [ 51 ], the "valley" areas are produced in non-simply connected environments and a boustrophedon or, in our case, an inward or outward spiral path will cover them at least twice. Like so, the algorithm of the previous section is extended and improved by introducing an inlet sensitivity coefficient I . The new algorithm (see Algorithm 7) identifies the existence of inlets in the decomposition schema by comparing the total path length for different I values. In that manner, the resulting path visits isolated areas first before continuing to the rest of the area. However, since our goal is also to produce smoother trajectories, θw is also considered in the evaluation of I.
62 Chapter 4. Coverage plans (a) (b) (c) (d) Figure 4.2 In ( a ), area partitioning for three UAVs on the same location as in Figure 2.13. White areas indicate the no-fly zones, whereas the black triangles show the initial positions of the UAVs. The FoV sized cell distribution is shown along with the centres of the triangles. In ( b ), a detailed view of the UAV 3 waypoint prioritization, while in ( c ) latitude and longitude information received during the simulated flight of that UAV. Finally, ( d ) shows the the total coverage of a squared FoV simulated sensor, working at a very slow rate of 1 Hz. Let Gi be the centroid of cell vi and Gij the centroid of its j neighbouring cell. Consider the maximum euclidean distance Lamong all neighbouring centroids of the CDT: kLk=max i∈Nd(Gi,Gij),1≤j≤3(4.2) where Nis the number of waypoints. This distance is chosen as the initial and minimal value of I . Please note that I is actually a multiplier for L . As such, the initial value I=1 equals I=1∗L . Then, different I values are actually multiplier values for L. Let us consider two sets θw(α,I) k i=I =θwi,θwi+1...θwk=Θ Θ Θ(4.3)
4.2 Waypoint prioritization 63 Figure 4.3 A waypoint list of a triangulated area computed with the coverage Algorithm 6. The regions in A and B create a valley in C in the borders-to-center cost attribution schema as their surrounding region has higher costs. Since the algorithm initially visits cells according to their distance from the borders, a lot of back and forth paths are created in region C. and P(α,I) k I=i =Pi,Pi+1...Pk=P,(4.4) where α is the area to be covered, k is an upper limit for the I value and Θ Θ Θ,P are the subsets of the original sets, for the evaluated I’s. Let us also consider the functions f1= Θ Θ Θ normΘ Θ Θ (4.5) and f2= P(−I) normP (4.6) where P(−I)is the mirror function of P, and their difference ε=f1−f2(4.7) Then, the path enhancement problem is to find the optimal I , which occurs in the intersection of f1 and f2 , where ε=0 . This problem is a multi-objective optimization problem, where the optimal I value is the one which minimizes the total path P and maximizes θw.
70 Chapter 4. Coverage plans (a) (b) (c) Figure 4.9 An area segregated in triangular cells by a FoV projection size of 50 metres. Figure (a) shows the produced mesh of 519 cells and the waypoint list with each vertex connected to the following with a line segment, according to Algorithm 6. The triangular white cell indicates the initial location of the UAV, while the white areas represent the obstacles or no-fly zones. Sub-figure (b) shows the simulated trajectory results for the Rascal110 model aeroplane from the latitude and longitude information received by the mavros node during flight. In (c) , an in-flight screenshot of the ground station during the simulation, where the behavior of the simulated autopilot can be observed; narrow turns produce paths outside of the expected trajectories. Table 4.1 Path metrics comparison when a classical grid decomposition strategy has been used for the areas shown in Fig. 4.10. Waypoint generation case (a) (b) (d) (e) (f) Path length (m) 8132,26 10580,86 6741,16 6988,05 10315,48 Min angle (deg) 0,94 26,98 18,43 0,01 0,01 Max angle (deg) 179,91 179,94 151,97 173,66 139,4 Avg angle (deg) 120,22 147,14 84,87 87,75 108,98 % of path over no fly zone 0,032% 0,026% 11,467% 9,521% 8,39% the results suggest that: (a) the algorithms provide smoother trajectories since the average angles of the paths are wider and (b) there are no paths produced over no-fly zones or outside the areas of interest. Inlet coefficient method In order to identify the geometrical properties and spot the qualitative differences in different inlet sensitivity coefficient settings of Algorithm 7, two areas have been chosen. The first is a sample area with several “valley” areas. The second is a scenario in an actual area presented in Fig. 4.11. On the latter, a simulated flight has been performed in order
4.3 Results 71 (a) (b) (c) (d) (e) (f) Figure 4.10 A comparison between the proposed method and a classic grid decomposition with boustrophedon coverage motions. The first two figures use the CDT schema as described in the previous Chapters, resulting in a mesh of 391 triangular cells. Sub-figures (a) and (b) present a different inlet sensitivity coefficient in order to highlight the resulting paths. The rest of the sub-figures show a classic grid overlay on top of the area, having a total of 210 square cells of which 78, in orange color, partially having a no-fly zone in them. The black square is the respective initial position. Figures (d) and (e) show two boustrophedon strategies which have been tested in order to compare them with the strategies of (a) and (b) whereas Fig. (f) uses the same border-to- center cost strategy of Section 3.2.2 in order to obtain a waypoint list plan by Algorithm 6, properly adjusted for a grid decomposition. to obtain results of the real trajectory followed. Please note that in order to show the difference of different sensitivity settings, the Ivalues were selected implicitly. Sample area This area has been segregated and partitioned for two UAVs; one of them visits three out of four “valley” regions, as it can be seen in Fig. 4.12. The two settings that have been chosen demonstrate the effects that this algorithm might have: on lower inlet coefficient values, some vertices might be identified as inlets while they are not. As a result, repetitive turns are produced and the paths are far from smooth.
72 Chapter 4. Coverage plans Figure 4.11 On the left, a region in Great Abaco island near Miami. The red polygons define the area constraints. The blue boxes indicate the initial positions of each UAS. On the right, the same region partitioned for the three UAS. Their initial positions are the black cells, whereas each shade of orange is the border- to-center cost. Their coverage percentage in relation to the whole area is 30%, 40% and 30%. The CDT is a result of 10 Lloyd iterations. Figure 4.12 Spiral paths for different inlet coefficient values. In the left case, the algorithm identifies a lot of next steps which are considered inlets, making a lot of turns. On the right, the paths are smoother. In the black circles, some points where the algorithm identified a far next movement and chose the closest cell instead. On higher inlet coefficient values, the paths are smoother and the inlets are identified correctly, but the total length of the path is larger, as shown in Table 4.2. In order to compare these results with the initial Algorithm 6, tests have been conducted for three different Lloyd smoothing settings: 10, 30 and 60 iterations. The algorithm has been tested with a inlet coefficient of 4 times the kLk value. The upper panel of Fig. 4.13 shows that
4.3 Results 73 Figure 4.13 A comparison between Algorithm 6 (depicted in red color) and Algorithm 7. In the upper panel the average distance between each cell might be larger but as the distribution curve shows, also clear in Table 4.3, the total path length is smaller. In the lower panel, the angle distribution is shown. While there is no big difference, the distribution curves show that with the new algorithm, there is a shift towards larger angles, thus smoother paths. Both panels are in columns which represent the different Lloyd smoothing settings (10, 30 and 60). Table 4.2 Metrics on inlet coefficient corresponding to Fig. 4.12. Inlet coefficient 2*kLk4*kLk Total Path(m) 7139,51 8515,79 Min distance(m) 9,87 10,63 Max distance(m) 226,51 191,72 Min angle(degrees) 20,60 5,09 Max angle(degrees) 179,889 179,99 the distribution dispersion of distances between each cell is wider. This results in longer total path lengths, as it is also shown in Table 4.3. In the lower panel, the distribution curve shows that in average the new algorithm also manages to produce smoother paths, by increasing the angles between each three waypoints.
74 Chapter 4. Coverage plans Table 4.3 Comparison between Algorithm 6 and Algorithm 7, corresponding to Fig. 4.13. Algorithm Lloyd Total path (m) Distance Min (m) Distance Max (m) 10 7736,95 5,77 255,15 Initial 30 7322,59 6,02 290,29 60 7502,25 12,99 233,97 10 6958,26 10,41 144,13 Current 30 7144,34 12,99 174,32 60 7098,88 12,99 94,11 Figure 4.14 The difference in coverage path planning can be identified mainly in the number of sharp turns, repeated paths and total path length. The image on the left represents an inlet coefficient value of 2 whereas the one on the right a value of 4. The total metrics of the example can be seen in Table 4.4. Actual area scenario The area of the upper right UAV in Figure 4.11 has been selected as an actual area scenario for a mission. In Fig. 4.14, the qualitative difference in the calculated waypoint list path can be identified and Fig. 4.15 shows the angle distribution in the two different cases; it shows that the dispersion of angles is greater in smaller inlet coefficient values but the produced total paths are smaller. On the other hand, bigger I values provides wider angles, thus smoother paths. Then, the total metrics of the produced paths can be seen in Table 4.4. As shown from Fig. 4.15, both increased Lloyd iterations and inlet sensitivity values manage to increase the average angle metric. The actual simulated trajectory can be seen in Fig. 4.16(a) and a screenshot of the simulated experiment in Fig. 4.16(b).
4.3 Results 75 Figure 4.15 The angle distribution in the produced path of the Abaco region. Red, blue and green colors represent the different inlet sensitivity coefficient values of 2,3 and 4 respectively. Each panel shows the different Lloyd iteration settings. The horizontal axis distributes the angles in 5 bins, whereas the vertical axis show the percentile rank of each bin. Vertical lines show the average angle of each case. (a) (b) Figure 4.16 In (a) , the actual trajectory when the simulated model followed the produced waypoint plan. In (b) , the APM Planner 2 ground control station showing the loaded waypoints. The mission had a flight setup of 100 metres altitude. Table 4.4 Metrics on inlet sensitivity values corresponding to Fig. 4.14. Inlet coefficient value 2*kLk3*kLk4*kLk Total Path(m) 3807,88 4185,54 4344,05 Min distance(m) 12,90 12,90 12,90 Max distance(m) 216,87 127,52 118,24 Min angle(degrees) 11,5629 11,69 0,8535 Max angle(degrees) 179,667 179,667 179,913
76 Chapter 4. Coverage plans Figure 4.17 Angle distribution in different Lloyd iteration trials. The red Gaussian curve which represents the lesser number of iterations tested, shows that there is a larger dispersion in angles, resulting in more “sharp” cells. While the best performance is achieved by the 100 iterations, where there are more angles reaching the 60o ideal case, the required computation time should be considered. Depending on the application, a solution closer to 50 iterations could be competent enough. Configuration comparison Regarding the use of Lloyd optimization for creating a smoother grid, results show that the angle distribution of the grid is more coherent with increased number of Lloyd iterations. Nevertheless, this operation is strictly dependent on the area input as well as the constraints of the CDT. In the distribution shown in Fig. 4.17, in the 100 iterations case, the algorithm did not manage to overcome the 65 iterations thus not been able to smooth the grid any more. The complexity effect of the algorithm is clear in Table 4.5, where the time each smoothing procedure needs significantly increases, thus exposing the trade off among sample overlapping and smoother or shorter paths. 4.3.3 Field experiments Two sets of experiments (M1 and M2) were conducted to test the proposed algorithm. The first set is a classical boustrophedon mission (M1) in the same area and it is used in order to obtain a comparative result of the classical method [ 4 ] with the proposed method. The second set of flights demonstrates the proposed method while it includes an optimised
4.3 Results 77 Table 4.5 Time in seconds for Lloyd optimization in the whole area and in each UAV configuration space. Lloyd Iterations CDT 10 60 100 All 0,2s 0,86s 1,008s UAS1 0,076s 0,36s 0,32s UAS2 0,064s 0,316s 0,392s UAS3 0,156s 0,552s 0,6s (a) (b) Figure 4.18 In (a) , the area which the experiments took place. The red external and internal borders represent the limit and the internal no-fly zones. In (b) , the Skysurfer fixed-wing platform used, operated in a maximum speed of 15m/s. tuning with a stringent acceptance radius (M2). Please note that the term acceptance radius, refers to the UAV-waypoint distance where the waypoint is considered visited. The experiments have been conducted in the agricultural area near Brenes, Spain (Figure 4.18a). A COTS fixed wing platform(Figure 4.18b) has been used which had the proposed framework on-board. Moreover, a square theoretical sensor with a FoV projection of 100 metres side is considered in order to provide a metric of coverage. The comparison of the two trajectories (see Figure 4.19) shows that the proposed method has managed to significantly reduce the flight over no-fly zones ( Qout ) and total flight length (Q), while at the same time performing a slightly better coverage (see Table 4.6). The missions were executed in different meteorological conditions during August and September in the Spanish province of Seville where the winds are relatively low and the temperatures are high, as they can go above 50o C. These extreme conditions affect
78 Chapter 4. Coverage plans (a) (b) Figure 4.19 In (a) , the trajectory of the boustrophedon path and the respective theoretical sensor coverage (in blue shade). The red external and internal borders represent the limit and the internal no-fly zones. As shown, several violations, in green ellipsis, of the no-fly zones have occurred in turns. In (b) , actual trajectory if the proposed method, with a calculated inlet sensitivity coefficient value of 2,35. As before, the blue shade represents the theoretical sensor coverage, while the ROI violations are shown in green circles. The flight over the no-fly zone is limited, whereas the sensor coverage is extended. Table 4.6 Comparative experimental results for the boustrophedon and proposed method missions. Mission ID Qout (m) Q(m) Coverage (%) Qout (%) M1avg 775,09 7204,8 94 10,75 M2avg 96,65 5925,84 97,2 1,63 the performance of the power sources, therefore, the appropriate tuning is necessary to increase flight endurance in order to finish the missions within appropriate safety levels. This method will be described in the following Chapters. 4.4 Conclusions The experimental results show that the proposed algorithm aids the process of reducing the overall flight path length and time, while reducing the flight over no-fly zones, in complex areas. The latter also supports the safe separation consideration, in context with the main contribution of this Thesis. Even though initial simulations showed that a boustrophedon P path was shorter, the actual experiments demonstrated that this is not actual the case. This occurs due to the aforementioned extra trajectories which occur during turns. Regarding
4.4 Conclusions 79 the average turn radius, our intuition was correct: as smoother paths are provided by creating wider turns, less energy consumption has been achieved. We have to note that except the non-holonomic nature of the considered vehicles, flight or vehicle dynamics have not yet been considered. The following chapter will describe in detail several considerations and cases on those issues. Nevertheless, it is clear that holonomic vehicles like multi-rotors would not demonstrate this behaviour, as they are able to perform sharp turns and follow the produced waypoint plans with more precision. However, these considerations are somehow application specific and heavily dependent on the task to be performed.
86 Chapter 5. Path tracking • The vehicle has to pass as close as possible over every waypoint of the planned path. • The overshoot after visiting a waypoint has to be as small as possible, in order to respect the area constraints as described in the previous Chapters. • The solution has to take into consideration vehicle dynamics and delays of every platform to be used at, especially in COTS systems. In that manner, area A which the UAV will fly over which is not over the desired path, has to be minimised (Fig. 5.4). Consider the total area integral error of the current position Uc(t) of vehicle U to the ideal position Ur(t) on the edge Ew,w+1 from waypoint w to the next, while Tis the duration between the two waypoint visits. minZT 1 [[(e(t−1) + e(t))h]/2]dt (5.6) where e(t) = Ur(t)−Uc(t) at time t and h is the height of the trapezoid formed by Uc(t−1),Uc(t),Ur(t−1),Ur(t). Initially the path tracking algorithm for every straight line segment can utilize a lookahead distance similar to Eq. 5.4. However, preliminary experiments showed that pure delays and latencies make this method not adequate for passing exactly over the waypoints and might produce cut corners and paths like the one in Fig. 5.1b, where the control loop reports goal points 2 and 3 consecutively, even though they belong to a different straight segment. In those areas, the equation of the circular path (Eq. 5.5), even for an instantaneous moment, seem more fit. The pure pursuit strategies mentioned in the previous Section constantly fit an arc between the vehicles current position and a goal point. These are extended as follows: the algorithm selects a goal point PG in the path to follow, at a lookahead distance L from the closest point in the desired path with respect to the current position of vehicle U . In order not to proceed to the next straight line segment prematurely, the vector of the straight line path is extended in the same direction until the waypoint is visited. Then, by combining the two strategies of straight path following and circular path in order to treat delays that might occur as discussed before, for every next turn a different lookahead distance is chosen. Since the size of that angle is correlated with the curvature, lower or higher Ls values can be expected in different turns. This means that L is expected to follow the pattern of R in relation with θ of Eq 5.2 and Fig. 5.2, decreasing on the interval (0,π/2) and increasing on the interval (π/2,π). The preliminary experiments described in the following Section, validated that indeed a fixed distance does not reduce the previously mentioned area A (Fig. 5.5); the reduction is
5.3 An offline-online hybrid approach 87 Figure 5.4 A pure pursuit algorithm. In a desired path of a W=w1,w2,...,wi set of waypoints, the closest point on the path from the current position of the vehicle U is chosen and a L lookahead distance is used in order to obtain the goal point PG . The proposed method tries to reduce the shaded area A , using each previous and current position, e(t−1) and e(t) accordingly, for several L values. In the conducted experiments, current vehicle position is obtained 3 times per second. Figure 5.5 The Lcoefficient which reduces area Ais different for different θangles. dependent on a dynamic distance coefficient, based on the angle θw between two segments of the waypoint plan. Hence the reduction of area A is dependent on the properties of the underlying UAV model, like its speed and angular velocity, the angle of the turn to be performed and the rate of which the next goal point is given. In that manner, a series of offline simulations are performed, by introducing a method which tests several lookahead distances L against different angles. This procedure produces a lookup table of a distance coefficient L in relation with the angles. Then, during flight, the pure pursuit algorithm uses this table in order to select the appropriate lookahead distance values for every next waypoint. Since this approach is dependent on the underlying vehicle dynamics model, it is characterized as model based - angle dependent.
88 Chapter 5. Path tracking Figure 5.6 On the left, the Odroid XU4 companion computer. On the right the PX4 autopilot. 5.4 Results As in the previous Chapter, the proposed algorithm has been tested on a simulation framework as also as in field experiments; the simulated setup also permitted the extraction of the aforementioned lookup tables. The same Software In The Loop (SITL) simulation framework of the previous Chapter has been used for simulations, whereas an Odroid-XU4 1 computer has been used as an onboard companion computer for the field experiments, along with a Pixhawk PX4 2 autopilot (see Fig. 5.6). The onboard computer has the same experimental framework as the simulations, except the visualisation tools. In both cases, the APM Planner ground station software has been used in order to observe the trajectories in real-time (see Fig. 5.7). 5.4.1 Simulations The UAV model used in the simulated flights is a fixed wing Rascal11 model airplane. Its on-board controller is the open source autopilot Ardupilot and has been combined with the JSBSim flight dynamics model simulator 3 . The Arduplane controller used is the Pixhawk Flight Management System [59]. Lookahead table extraction Two sets of simulations have been conducted: the first is intended for checking the validity of the strategy and the extraction of the aforementioned lookup table and has been tested against the onboard autopilot control model, which acts as a trajectory tracker based on 1http://www.hardkernel.com/main/products/prdt_info.php?g_code=G143452239825 2https://pixhawk.org/modules/pixhawk 3http://jsbsim.sourceforge.net/
5.4 Results 89 Figure 5.7 An APM Planner 2 screenshot during the simulated flights of the first experiment for the lookahead table extraction, showing different trajectories of the vehicle based on different L values. Note that the outer right trajectory belongs to the control one, without using the algorithm presented in this paper. Also note that the ground control software updates the map every second, hence the rough, orthogonal looking trajectories. the study presented in [ 67 ]. The second is a set of simulated test flights in order to observe the results by using the previously extracted look-up table. For the calculation of the lookup table, the UAV had a constant platform specific speed of 20m/s. In order to include the expected delays of the experimental framework for the calculation of L , we had to take into consideration the rate at which every new goal point is sent to the controller. By examining the results of Fig.9 in [ 68 ], we could expect that a rate of 3Hz (since this is the limit the experimental framework permits us to update the actual position of the vehicle) in giving a new command would result an Ls value of approximately 2,3. Considering that V=20m/s and Ls=2,5 , and since from Eq. 5.4 L=Ls∗VT a minimum of 45 meters has been tested as the L value for the lookahead distance. Then, a random incremental set of angles of 14,27,48,76,90,119,136 and 154 degrees has been tested against the proposed method (Table 5.1 and Fig. 5.8). The lookahead distance which reduced the overhead area A for these settings, indicated the lookup distance coefficient value to be used on the respective turns, which was used in the second set of simulated experiments and field tests. Accordingly, the comparison between using or not the proposed method for the individual turns can be seen in Fig. 5.9 Moreover, the relation between the different values for the lookahead distance coefficient and the tested angles, can be seen in Fig. 5.10.
90 Chapter 5. Path tracking Table 5.1 Lookup distance coefficient table based on different path angles. Angle (degrees) Best lookup distance coefficient L (m) Area (km2) 14 55 0,005448 27 56 0,005404 48 55 0,004911 76 53 0,003239 90 45 0,002321 119 51 0,001245 136 52 0,000717 154 65 0,000564 Figure 5.8 Angle vs Ldistance. Simulated test flight The second simulated experiment is a full scale scenario of an actual area of river Douro in Porto, Portugal, where a UAV flying at 20m/s, performs a coverage task (see Fig. 5.11). Please note that on the indicated path, the average turn angle is of 119,84o . Also note that there is a specific radius Rv around each waypoint to be considered as visited. A total of nine scenarios have been tested; the first is without the proposed method and by letting the onboard controller to pass by a waypoint on a loose trajectory, meaning that the UAV could pass at a Rv=40m distance, if this would produce a smoother trajectory. This method is often used in fixed wings in order not to perform intense manoeuvres and it is the default behaviour of the SITL framework. The second case was also without using the proposed method, but restricting the on-board controller that the aforementioned
5.4 Results 91 Figure 5.9 The blue line represents the total area outside of the straight line paths by using the autopilot’s default path tracking algorithm. The green line is after the proposed method. Figure 5.10 The effect of different L lookahead distances (horizontal axis) in relation with the integrated area outside the straight line paths (vertical axis) and the angle of the next segment (individual coloured lines). The vertical lines indicate the average of the lookahead distance, whereas the dotted lines on the left and right indicate the standard deviation.
92 Chapter 5. Path tracking Figure 5.11 On the left, the estuary of Douro river in Porto; the red lines indicate the borders of the area of interest. On the right, same area as has been decomposed for one UAV, having an initial position on the blue cell. A coverage task is performed by visiting all the waypoints, indicated by the centroids of each cell. The green line indicates this path. relaxation would have a radius of Rv=10m ; a restriction which our proposed framework considers. The next six test cases were performed by using the proposed method, but with a static lookahead distance for all turns, which was extracted from the minimum and maximum best L values of Table 5.1, in a 5 meter step: 40, 45, 50, 55, 60 and 65 meters. Finally, the last scenario is by performing each turn by looking up on the next turn and changing the lookahead distance accordingly, demonstrating the proposed algorithm. As the results indicate (Figs. 5.12 and 5.13), when imposing the same restrictions to the autopilot, all cases of the proposed method have performed better, while the definitive last case (see case "Adaptive L" in Fig. 5.12) has managed to reduce the overhead area against all other instances. Another simulated flight has been conducted, this time in the same area and waypoint plan as of Fig. 4.19b of the previous Chapter. The proposed method was tested in two different speed scenarios of 10m/s and 20m/s. Figure 5.14 shows the four trajectories where the default method is tested against the proposed method. The prioritisation of passing over the waypoints with the proposed method is evident, while the area of error was reduced by an average of 18.79%. 5.4.2 Field experiments The same area of the experiments of the previous Chapter has been used in order to conduct a set of comparative test flights. Two sets of experiments have been conducted, with two different UAVs (see Fig. 5.15). The first set uses the platform already presented in the previous Chapter and extends those experiments, by introducing a tuning control strategy
5.4 Results 93 Figure 5.12 The integrated area for each test case of the second simulated experiment. The "AUTO_40m" case refers to the default autopilot behavior, without always passing over the waypoints. On the other hand, the strict "AUTO_10m" case is restricted on passing over them at least 10 metres away. The rightmost "Adaptive L" scenario shows the result by making use of the proposed method, whereas the rest refer to the static Lcoefficient for all turns. of the platform. The second set is focused only on the path tracking method in order to compare it with the default one. Table 5.2 shows the characteristics of the field experiment framework A tuning control strategy The controller gains of the selected autopilot need to be tuned to meet flight duration and energy consumption requirements. Moreover, the system has to be able to execute its missions in different conditions, both platform-related and ambient-related, e.g. payload or weather. Hence, the performance should be tested in different conditions, with different payloads, and emphasizing the presence of wind disturbances, which may alter the mission performance. Two more sets of experiments were conducted (M3 and M4) to test the tuning strategy along with the rest of the methods, with M1 (boustrophedon) and M2 (proposed) being from the previous Chapters. The platform tuning was performed by a set of ten flights, by using a sharp turn pattern (M3) of Fig. 5.16. Then, mission M4 is using the same methods as already proposed, but with a relaxed waypoint acceptance radius, in contrast with strict acceptance radius of mission M2 as
94 Chapter 5. Path tracking Figure 5.13 As a qualitative result, this groundstation image shows the recorded flight trajectories on different test scenarios. Denoted by the blue circles, some cases where the overshoot of the autopilot is evident. presented in the previous Chapter. A tuning baseline is required in order to adjust the control parameters of the platform and is obtained with the aircraft in manual mode, i.e. with no assist, in which a series of manoeuvres are performed to find the desired characteristics with certain throttle Tin , airspeed Va , climb rate Vh , pitch θ and roll φ angular demands. The baseline calibration was obtained in steady conditions with low winds and a temperature of 25oC . The goal is to modify different parameters to achieve the optimal characteristics of the flight. These parameters are: the feed-forward roll gain Kφ f f , the proportional (roll-rate) gain Kφ P , the integral (roll trim offset) gain Kφ I , the feed-forward pitch gain Kθ f f , the proportional (pitch rate) gain Kθ P , the integral (pitch trim offset) gain Kθ I . To illustrate the attitude control laws that need to be tuned, the roll φ control loop is depicted in Fig. 5.17, which is equivalent to the pitch and yaw control laws. As discussed before, the path tracking strategy of the PX4 autopilot software is the same described in [ 67 ]. Its control module has a L1 controller which is implemented in two stages: one in which a reference waypoint L1R is calculated on the trajectory that has to be followed, and a second one computing the lateral acceleration needed to reach L1R , which
5.4 Results 95 Figure 5.14 Comparison between the onboard tracker trajectory of the autopilot and the proposed method, for two different speed profiles. The error is shaded in each turn. Please note that the overlapped regions show the smaller error first. The black line indicates the waypoint path as sent to the autopilot. (a) (b) Figure 5.15 In ( a ), the MarineUAS-1 platform used for the coverage and tuning flights. In (b), MarineUAS-2 used for the path tracking experiment. is equivalent to the roll angle command shown in Fig. 5.17. The main advantage of this method is the reduction of configurable gains, as the Roll PD and cross track errors are replaced with the reference distance L1R . In addition, this function presents improvements in disturbance rejections and circular path tracking in the presence of wind and this also gives the possibility of combining it with a higher level trajectory generation algorithm.
102 Chapter 6. A decentralized framework for coastal missions architecture for task assignment, each of the UAVs can change roles and serve as data relays or as data acquisition agents. Then, a multi-layered architecture is proposed by providing a test case of a buoy searching team of UAVs, as described in the previous Chapters. This Chapter presents the complete architecture for such a scenario, where the goal is to locate the positions of the buoys while maintaining connectivity. Moreover, the TOR is able to be reconfigured in terms of ROI and roles in case it is needed. The decomposition, partition, waypoint extraction and path tracking enhancement methods presented in the previous Chapters as a top-to-bottom architecture of path enhancement modules is presented in Fig. 6.1. In each level of this organisation, the following modules can be identified: •Complex area decomposition based on the FoV of the robots. • Partitioning for an arbitrary number of robots, based on their relative flight time autonomy capabilities, along with feasible Lloyd smoothing procedures on the resulting ROI. • Path enhancement through an enhanced waypoint list extraction, in terms of path length and overall turn angle. • A pure pursuit online path tracker which guarantees the passage over the waypoints while reducing the flight error. • A tuning control strategy which performs a COTS platform parameterization, as described in the field experiments of the previous Chapter. This framework is platform independent and applied on-board each UAV, providing a team-wide homogeneity in respect with the communication protocols and world representation. As described in the next Section, depending on the team set-up, each robot has a global knowledge of the map partition and local knowledge of its trajectory, position and gathered information. Then, depending on the evolution of the tasks, the team can be reconfigured. 6.1 A market based decentralized architecture In centralised approaches for multi-robot task assignment, a single robot or a coordinator is used in order to assign tasks to the rest of the team. In the UAS case, this role is usually held by a ground station inside the communication range of every vehicle or through a messagehop schema, where information and instructions are propagated through the network of vehicles. This approach can result in optimal solutions since the coordinator has global
6.1 A market based decentralized architecture 103 Figure 6.1 The top-to-bottom proposed framework architecture, with every path enhancement module denoted in yellow. Each operation is performed onboard each robot, thus providing an initial step for a decentralized solution. Please note that the airspace restrictions can be communicated to each individual vehicle, which will then reconfigure its ROI. knowledge of the vehicles’ states as also as any other relative information. Moreover, if the coordinator’s task is solely the planning of task assignment, the computational burden of each of the vehicles is solely on the mission execution. Nevertheless, centralised architectures are not robust, since they possess the single-point-of-failure disadvantage. Moreover, and especially in teams of aerial robots with limited communication range operating in large areas, the communication constraints imposed are restrictive, leading to critical latencies. On the other hand, in distributed approaches only local knowledge is available is available for the vehicles, which is exchanged with other vehicles in specific intervals or depending on their proximity. During those periods, decisions are made based on that information which may alter the task assignment, based on specific local and global utility functions. These approaches are more robust since there is no single-point-of-failure and can be characterised by an increased flexibility since every vehicle can be its own and its neighbours coordinator.
104 Chapter 6. A decentralized framework for coastal missions Nevertheless this also implies that the robots have an increased computational burden of planning the feasibility of these complex tasks and since they only possess local information, the overall solution can be sub-optimal. Considering the taxonomy proposed by Gerkey et al. [ 69 ] for task allocation in multirobot systems, a categorisation is proposed based on three categories: • Task type: S ingle R obot or M ulti- R obot, distinguishing the tasks that need one or more robots to be completed. • Robot type: S ingle T ask or M ulti- T ask, distinguishing types of robots which can perform only one or multiple tasks simultaneously. • Schedule type: I nstantaneous A ssignment or T ime-extended A ssignment, differentiating the ability to provide to the system information only for the current plan or to plan for future task allocations too. The authors point out that only the SR-ST-IA problem, which is an instance of the optimal assignment problem in combinatorial optimization, can be solved in polynomial time, while the rest of the combinations are strongly NP-hard. In our case and considering the approaches already presented in the previous Chapters, we can identify each UAV task of area coverage as an SR-ST-IA case. Each sub-area is a task which is assigned to a single robot and each robot can be assigned to this task only. Regarding scheduling, current task assignment can be altered only if new information is provided and no future scheduling is provided. Please note that this categorisation is assumed if, as noted, each task is the sub-area coverage. If the notion of task is defined as the visit of every waypoint, then the problem changes to SR-ST-TA, since a serial waypoint list has to be constructed, dictating in that way the future scheduling of single waypoint visit task for every robot. Following the aforementioned taxonomy of Gerkey et al. [ 69 ], Korsah et al. [ 70 ] proposed an extension in order to include interrelated utilities and task constraints. This taxonomy named iTax, has the following categorisation: •N o D ependencies: Problems which each sub-task for a robot is not affected by other tasks or robots. •I n-schedule D ependencies: The utility function of a robot for a task is affected by the other tasks this robot has to perform or is performing. •C ross-schedule D ependencies: The utility function of a robot for a task is affected by the other tasks this robot has as also as the schedules of other robots.
6.1 A market based decentralized architecture 105 •C omplex D ependencies: The utility function of a robot is affected by the schedules of other robots for that specific task decomposition, thus not decoupling the task decomposition and task allocation problems. As before, only the first category falls in the polynomial time/linear assignment problems. Our case of UAV-area assignment also falls into this category, since each sub-area and task of area coverage is not affected by the schedules of other agents. Once more, if we consider as a task each individual waypoint visit, then the categorisation is that of CD , since each partition solution results in different utility values. However, we have already treated this task decomposition and allocation in the previous steps. Nevertheless, this categorisation is ideal and in order to meet the aforementioned decentralized robustness and flexibility, some properties have to be identified, which will provide the system several fail-prone characteristics. Initially, it is assumed that the robots will not always be in a communication range to provide uninterrupted data exchange. Moreover, UAV loss mitigation is considered; this requires methods which will perform team reconfiguration and area reassignment operations. These methods also aid sub-area alterations, like the appearance of a new no-fly zone. This would change the utility function of a robot since it would have to cover less area, thus affecting the global utility function and optimality of the original partition. Following the definition of a market based architecture based on Diaz et al. [ 71 ], the aforementioned taxonomies, categorisation and pros and cons of the centralised and distributed solutions, a hybrid centralised/decentralised market-based architecture is proposed, which has the following characteristics: • The TOR has a global objective function of obtaining task related information through complete area coverage. This is achieved by performing a decomposition into sub-area individual tasks, depending on the number of robots and area restrictions. • The global utility function for the TOR, or for the auctioneer in a centralised marketbased architecture, is the maximisation of the information retrieval by utilising the resources of the vehicles. This is achieved by the capability aware partitioning algorithms as described in the previous Chapters. Like so, this function is dependent on the relative capabilities of the vehicles and their initial positions. It is also dependent on the number of the vehicles and the total area to be covered. • The local utility function for every vehicle is the minimisation of the path length and maximisation of the average path angle as described in the previous Chapters.
106 Chapter 6. A decentralized framework for coastal missions • The local utility function of every vehicle might change in case a new no-fly zone appears or the energy needs become higher, thus decreasing its relative capability with respect to the other vehicles. • The latter characteristic, dictates the need for a specific resources/needs schema, where the global and local utility functions affect each other. • This need requests of a mechanism which will be responsible of obtaining the individual utility function results and try to maximize its own utility, through a fail-prone communication schema. As such, initially a coordinator is chosen which performs the aforementioned initial partition of the whole area. This coordinator serves as an auctioneer, which is responsible of obtaining local information and compute an outcome which maximizes the global utility. Then, a team of communication relay vehicles is deployed in order to serve as a relay between the coverage UAVs and the coordinator. These relay vehicles configure their trajectory in order to communicate with the coverage vehicles, obtain local information and perform the communication with the auctioneer. The whole team of UAVs, relays and coverage, should be able to reconfigure depending on several user or mission criteria. These characteristics and methods are demonstrated in the next Section, where a test case of data acquisition from sea buoys in a marine biology scenario is required. 6.2 Test case: a multi-layer framework for buoy coordinated data acquisition Following the MarineUAS objectives, bringing attention to the challenging problem of using UAVs for marine sciences, a buoy data extraction scenario is chosen. Figure 6.2 shows the proposed 3-layer architecture for a Coordinated Data Acquisition (CDA) model: in the top layer, a coordinator runs the initial partitioning and waypoint extraction algorithms and monitors the whole operation. The middle layer is responsible for the communication between the UAVs and the coordinator using a feedback control schema, described in the following Sections. Finally, in the bottom layer, the individual sub-area UAVs are responsible for data extraction from the buoys. Please note that the UAVs used in either of the last two levels can be of the same type or not, but are denoted differently, based on their role: the buoy searching UAVs ( UAVSearch ), which search for the buoys and communicate with them and the UAV communication relays ( UAVComm ) which are used to enable communication between the searching UAVs and the coordinator.
6.2 Test case: a multi-layer framework for buoy coordinated data acquisition 107 Figure 6.2 A multi-layered framework for buoy data acquisition. Following the market-based paradigm as discussed before, each of the layers has to satisfy its own utility function: the coordinator has to make sure a fair partition is performed for every UAVSearch . The relay UAVs have to maintain the connectivity between the coordinator and the searching UAVs. Finally, the searching UAVs have to perform their coverage tasks, according with the path enhancement framework of Fig. 6.1. The coordinator and search UAV functions have been described before, while the way the UAVComm manages to satisfy its tasks, is described in the following sub-sections. 6.2.1 Data acquisition layer and buoy specifications Sea buoys are used in oceanic studies as a low consumption / large autonomy devices in order to gather data from either sea surface or from a specific depth. They are also used as data relays for undrwater sensors or vehicles. When a buoy radio-device is active and waiting to receive or transmit data, it wastes energy on idle listening. Since traffic loads are usually low in buoy networks, such idle listening can use a considerable amount of energy unless efficient communication mechanisms are employed. Therefore, in order to save energy, buoys do not transmit data through the radio device continously, but by utilising an interrogation technique. Common low-power RF devices periodically wake-up to check if another device is trying to make contact. However, a new type of low power device became available which only wakes-up when it is explicitly interrogated. In order to perform this interrogation, the UAV RF extraction transmitter has to send a RF signal to
108 Chapter 6. A decentralized framework for coastal missions Figure 6.3 Multiprocessor and wake-up RF circuitry structure. wake-up the buoy radio device in order to start the extraction. This wake-up process is performed by using wake detectors like the ones in [ 72 ] and [ 73 ]. These are simple circuits which include two nano-powered operational amplifier comparators in conjunction with a super regenerative RF receiver which wakes up the communication micro-controller from its low-consumption, deep sleep state. The buoy’s embedded central unit structure is shown in Fig. 6.3. In order to make contact and extract information from a given buoy ( Bn ) a UAVSearch has to be close enough to enable a wake-up, handshake and data transmission sequence. As described in Chapter 2 and shown in Fig. 2.9, the on-range area RA to enable transmission is determined by the transceiver specifications. In the given test case, the location of the buoys is unknown and the initial position of every UAVSearch is randomly generated. Then, for each of the searching UAVs, the area is decomposed, partitioned and coverage waypoint plans are extracted according with the methods described in the previous Chapters. An example is shown in Fig. 6.4. When a UAV locates a buoy, executes loitering manoeuvres around the buoy in order to gather the data. Once the transmission is finished UAVSearch goes back to the original plan. Then, the UAVComm team is used to recover the data. 6.2.2 Communication layer as a feedback communication control model This middle layer demonstrates a Feedback Communication Control (FCC) model, which actuates on a UAV when it has the role of a UAVComm . This control model actually replaces the waypoint path planner and path tracker of the proposed framework shown in Fig. 6.1. It is in charge of establishing a communication link between one or more UAVComm , multiple UAVSearch and the Coordinator . Depending on the restrictions which the aerial platforms will impose, the goal is to enable communications with a given desired Quality of Service (QoS) level, taking into consideration the energy consumption, number of vehicles, link
6.2 Test case: a multi-layer framework for buoy coordinated data acquisition 109 Figure 6.4 A buoy searching CDT example of a team of 5 UAVSearch searching vehicles, in Kobe, Japan. Each of the vehicles has different flight duration capabilities, while the communication cone of each buoy is the same, resulting in same size triangular cells for the decomposition. bandwidth and time. This involves control and real time scheduling as in [ 74 ], where the UAVComm , has the role of a relay node to enable communication. The feedback strategy employed in this layer has the dual task of commanding and optimizing the UAVComm motion, adapting the scheduler algorithm according to the desired QoS level. This strategy allows the controller to correct the discrepancies between a real-time simplified estimation and actual flying and network conditions. This model, presented in [ 75 ], converts the user requirements given by a QoS selection to motion instructions for the unmanned vehicles. First, the required network variables and the physical constraints serve as input variables for the controller. Then, periodically a set of estimated variables are calculated to predict the behaviour and create the UAV paths. After that, the actuation output variables are sent to the UAVs. The variables used are summarized in Fig. 6.5a while Fig. 6.5b shows the proposed FCC model. Once again, please note that in case a UAV takes the role of the UAVComm , then the aforementioned FCC model replaces the waypoint path planner and tracker modules of the framework shown in Fig. 6.1. In that manner, the flexibility of the system remains, as every UAV can take up any role, as described in the next Sections. 6.2.3 Single UAV dynamic role reconfiguration The dynamic nature and flexibility of the chosen architecture enables the robustness of the system against unpredicted conditions such additional energy consumption due to changing weather conditions. Moreover, the uncertainty of buoy positions may affect the
110 Chapter 6. A decentralized framework for coastal missions (a) (b) Figure 6.5 In (a), the FCC model variables and in (b), the FCC model. overall performance of the system and the deployment of the communication network. To tackle this situation in real-time, in addition to the feedback controller, the proposed market-based strategy evaluates the objective functions of each UAV and two strategies have been implemented to dynamically adapt the roles of the deployed UAVs. The first one is the role reassignment which allows a UAVSearch to act as an UAVComm and vice versa. In this way, the architecture is able to change the number of UAVSearch as the buoy locations are being discovered and the distances between vehicles increase. On the other hand, in order to balance the communication network load, meet the QoS requirements and avoid congestion issues, the proposed architecture is able to recruit idle UAVs and convert UAVSearch into UAVComm. Figure 6.6 shows the conversion diagram divided into the two steps: • Converting UAVSearch into a UAVComm . During the mission, the required network quality through the QoS scheme could not be achieved and more UAVComm may be needed, so the main trigger for this conversion is based on a network-event. This forces the reconfiguration for the whole team and a repartition of the remaining area into new sub-areas. This operation also occurs when one or more searching vehicles are missing. • Converting UAVComm into a UAVSearch . This is triggered based on events like a loss of a UAVSearch vehicle or heavy weather conditions which make the searching mission unfeasible. In both cases, the number of available searching vehicles changes. In the context of the market-based architecture, what changes is the available global resources for the searching task. Like so, the auctioneer, or the coordinator, has to provide the new areas or sub-tasks to each of the remaining vehicles. This process resembles the bidding process, where
6.2 Test case: a multi-layer framework for buoy coordinated data acquisition 111 Figure 6.6 Role reassignment diagram, based on events. every vehicle bids according to its current status: its position, the area that has already covered, its remaining relative coverage capability etc. Dynamic team reconfiguration The second strategy is the implementation of a decentralized team reconfiguration solution which increases the robustness of the system. It may occur that a group of UAVSearch and/or UAVComm lost connectivity with the Coordinator during the searching mission. It may also happen that during that time one or more of the isolated vehicles becomes inoperative. This scenario is shown in Fig. 6.7. After a period of isolation from the Coordinator , the UAVSearch and the UAVComm vehicles that are in the same area and are able to communicate between them, execute a regrouping algorithm. The algorithm workflow is as follows, 1. A local coordinator is chosen. Every vehicle acting as a network node and depending on its dynamic and communication properties, has a pre-assigned rank. When a disconnection from the Coordinator is detected by any vehicle, it starts to share data packages with its neighbours, propagating the rank of all vehicles. The behaviour is similar to the token-holder selection in token-ring local networks. After waiting a predefined time, the higher rank vehicle is selected as a local coordinator and it is the group leader (UAVGLeader). 2. Vehicle enumeration. The UAVGLeader starts to enumerate the available vehicles and the local objectives like the remaining area to search or the acquired buoys. Since the coordination capabilities of UAVGLeader may be weaker than the main Coordinator in a given scenario, the motion control of the vehicles is performed by using self-triggered and event-triggered controllers. The self-triggered controller allows to perform the mission while decreasing the amount of data shared between vehicles as in [76]. 3. Regrouping. It may occurs that one vehicle is again able to establish connection with the main group. In this scenario all communication relay vehicles ( UAVComm ) must
118 Chapter 6. A decentralized framework for coastal missions Table 6.2 QoS Simulations. Simulation number Energy Bandwidth No. Vehicles #1. Keep Connectivity No limit No minimum No limit #2. Long Endurance Efficient >=0 bounded Figure 6.11 Simulation scenario including 5 exploration UAVs and 2 relay UAVs. Enhancing link connectivity The results for this simulation are presented in Fig. 6.11. The QoS index selected for is enhancing the connectivity. Here, the relay UAVs try to keep network connectivity calculating the minimum number of vehicles needed to do so. There is no energy consumption limitation for the relay UAVs. They are only tied to the maximum values set for the dynamic model used. There is also no limit for the number of vehicles in use at the same time. In this case in order to keep connectivity with the searching vehicles the minimum number was two. Fig. 6.11 shows the path for the 7 vehicles. The dark blue and red lines are the paths of relay UAVs, switching from a go-to-location behaviour to a loitering one continuously as needed. The rest of the coloured lines are the UAVSearch paths. In Figs. 6.12 and 6.13 we present the simulated bandwidth between each relay node and their assigned exploration vehicles. The objective is keeping the attainable bandwidth always above 0 in order to maintain the communication link up.
6.4 Results 119 Figure 6.12 UAV Relay #1 link bandwidth to 3 exploration UAVs. 0 1000 2000 3000 4000 5000 6000 7000 8000 9000 Time (s) 0 10 20 30 40 50 60 Throughtput estimation Mb/s Estimated throughput to Relay2 against time Figure 6.13 UAV Relay #2 link bandwidth to 2 exploration UAVs. Saving total propulsion energy by reducing the number of UAVs and bounding energy consumption. In this simulated scenario, shown in Fig. 6.14a, only one vehicle is used as a relay (dark
120 Chapter 6. A decentralized framework for coastal missions blue line path). Figure 6.14c shows the relative power consumption of the relay UAV, while the measured attainable bandwidth of all exploration vehicles to the relay UAV is shown in Fig. 6.14e. Please note that in this case study, there is no connectivity with some of the exploration UAVs during certain periods of time. In addition to the reduced number of vehicles, we have also added a reduced power consumption QoS level. In order to demonstrate the difference when the UAV relay power consumption is limited, sub-figures (b) 1 ,(d) and (f) of Fig. 6.14 show the relay UAV path, the relative power consumption and the difference in attainable bandwidth connectivity and when limiting the relay UAV speed. 6.4.2 Role reassignment and team reconfiguration experiments The following simulation illustrates a role reassignment scenario and its connection with the market based architecture. A set of UAVSearch is deployed for individual area coverage and a set of UAVComm is deployed as a communication network. At a given moment, the number of available UAVComm to act as relays is not enough to fulfil the user network quality requirements. Thus, in order to avoid the network overload or a disconnection scenario, the Coordinator auctions the conversion of 2 UAVSearch into UAVComm . Like so, the searching layer is reconfigured to cover the non-searched area using the remaining UAVSearch . From this point of view, the behaviour of the system is the same as if those UAVSearch are missing. In those cases, when the global resources change, the Coordinator acts as the auctioneer and the UAVSearch vehicles as bidders. The bid of every vehicle is calculated from the area that has already covered, the area which remains on its initial plan and its remaining flight time. Then, the Coordinator performs the auction based on the local information, or bids, of every vehicle. This information also includes any area changes, like the appearance of a new no-fly zone. In Fig. 6.15, the green area is the already travelled paths, the red triangles with the numbers are the current positions of the remaining UAVSearch and the red X’s are the last seen positions of the other two converted UAVs. In Fig.6.16 the reconfigured scenario is shown. Team reconfiguration In order to test the dynamic reconfiguration using self-triggered control we run a simulation to test the behaviour of the lost vehicles when they are trying to recover contact. The objective of this simulation is to move all together toward a given position but keeping the connectivity during the movement and reducing the number of exchanged messages between them to avoid unexpected network congestion. In Fig. 6.17 1Fig.6.14b: https://vimeo.com/user6158089/fbfig614b
6.4 Results 121 (a) (b) 0 1000 2000 3000 4000 5000 6000 7000 8000 9000 Time (s) 0 200 400 600 800 1000 1200 1400 Power coefficient (W) Relay 1 Power consumption difference with cruise speed against Time 0 2000 4000 6000 8000 10000 Time (s) 0 200 400 600 800 1000 1200 Power coefficient (W) Relay 1 Power consumption difference with cruise speed against Time (c) (d) 0 1000 2000 3000 4000 5000 6000 7000 8000 9000 Time (s) 0 5 10 15 20 25 30 35 40 45 50 Throughtput estimation Mb/s Estimated throughput to Relay1 against time 0 2000 4000 6000 8000 10000 Time (s) 0 10 20 30 40 50 60 Throughtput estimation Mb/s Estimated throughput to Relay1 against time (e) (f) Figure 6.14 In ( a ) a long endurance scenario where the relay UAV path is denoted by the dark blue line. In ( c ), the unlimited speed UAV Relay relative power consumption, while in ( e ) the link bandwidth measurements to all 5 searching UAVs is presented. On the other hand, sub-figures ( b ), ( d ) and ( f ) show the respective limited speed case.
122 Chapter 6. A decentralized framework for coastal missions Figure 6.15 Searching algorithm waypoint plans. The green area is the already travelled paths, the grey triangles with the numbers are the current positions of the remaining UAVSearch and the black X’s are the last seen positions of the other two converted UAVs. The black cells indicate the initial positions of the team. Figure 6.16 Reconfigured waypoint plans after team reassignment. The black cells indicate the initial positions. Note that the already covered areas are not considered no fly zones but the repeated coverage is not preferred.
6.4 Results 123 0 20 40 60 80 100 X [m] -10 0 10 20 30 40 50 60 70 80 90 100 Y [m] Position of robots, position of base and moving vehicle Robot 1 Robot 2 Robot 3 Moving Vehicle Base Figure 6.17 Motion simulation for 3 Relay-UAVs platoon deployment. 0 200 400 600 800 1000 1200 Time [s] 5 10 15 20 25 30 35 40 45 Attainable data rate [Mb/s] Link data rate seen from Base to End Node Link 1: Base-Relay1 Link 2: Relay 1-2 Link 3: Relay 2-3 Link 4: Relay 3-End Node Max. Bandwidth Figure 6.18 Estimated link data-rate for 3 Relay-UAVs platoon deployment motion simulation. Non-dashed line is the effective attainable bandwidth. we present a simulated scenario with 3 UAVComm previously disconnected from the main Coordinator . They are moving towards a location where they can communicate with the Coordinator directly or through another connected UAVComm . In Fig. 6.18 the impact over the available bandwidth during motion is shown and how the effective bandwidth (green non-dashed line) always exists between extreme nodes as required.
124 Chapter 6. A decentralized framework for coastal missions 6.5 Conclusions In the multi-robot task allocation domain, the main question which describes the problem is "which robot handles which task". In this Thesis, we have tried to decompose this question, providing an architecture which analyses the roles of each part of a UAS, proposing a method which adds flexibility in a variety of applications. The architecture discussed in this Chapter addresses the challenges of role reassignment, reconfiguration and replanning of a TOR, by demonstrating how a control framework for communication purposes can fit in the framework presented in the previous Chapters. The presented test-case of a middle communication relay between a different parts of an UAS and the coverage UAVs shows that the single-UAV framework can be part of a flexible architecture and be used by several distributed systems. The results indicate that the application of a market-based approach is dependent on the application itself and the QoS properties which are user defined.
7 Conclusions and future development Through the sorrow all through our splendour Don’t take offence at my innuendo The Queen In this chapter a summary of the main contributions and results is presented while highlighting the future directions this line of research can follow. Moreover, several advantages and disadvantages are discussed, while proposing some methods to address the issues faced. 7.1 Conclusions The main contribution of this Thesis is the proposal a modular, multi-vehicle area partitioning and coverage framework for complex area coverage tasks. The design of the proposed framework takes into consideration area and platform characteristics which impose restrictions, in the context of modern team configurations and tasks over coastal areas. The proposed methods and algorithms address the issues of complex area decomposition and partition for an arbitrary number of robots, depending on their on-board sensors and relative capabilities. These properties increase the flexibility of the system since different robots can acquire any role in a given mission, and the partition can be recalculated accordingly; this is also true in case the area restrictions change over time. Then, coverage 125
126 Chapter 7. Conclusions and future development waypoint flight plans are derived, aiming at producing full coverage paths having in mind the non-holonomic nature of fixed-wing robots. The latter is a design choice derived from the nature of low-cost COTS platforms used in a vast majority of modern UAS for research or commercial purposes. This Thesis also proposes a path tracking method which can be applied in any autopilot hardware by the addition of a low cost companion computer. Moreover, with the platform parametrisation methods we have presented, shorter paths and energy consumption is observed. The modular and decentralised nature of the proposed framework implies that it can be applied in any coastal and off-the-coast task with little effort. Even though this Thesis is focused on teams of aerial non-holonomic robots operating on coastal areas, its application on other types of vehicles or areas is rather trivial. Moreover, the multi-layered treatment of an area as a graph, proposes that this framework can be extended to include any type of heterogeneous robot or sensor, as shown in the final Chapter. 7.1.1 Decomposition of complex coastal areas In Chapter 2, a triangulation method by the means of a CDT is presented. The constraint of this triangulation is the FoV projection of the respective sensor to the surface to be covered. In that manner, a coupling of an area approximation with the actual characteristics of the platform is achieved. This method manages to create an abstraction of the continuous space and extract a list of the future coverage waypoint plans, since every centroid of every triangular cell will serve as a waypoint. If every waypoint is visited then complete coverage can be considered. The produced mesh is an Apollonian network, while the resulting graph structure of connecting each waypoint is a third degree undirected graph; every waypoint can have at most three neighbours. This property implies a simplification of the algorithms as presented in the following chapters over a 6-neighbour method or a square grid decomposition which would have at least four (4) neighbours or at most eight (8). The results and test cases show that in comparison with a square grid decomposition, a CDT manages to perform better regarding complex shapes. This property is evident in non-square borders and no-fly zones, since a square grid decomposition cannot have variable orientation for every square and maintain graph connectivity, while a triangulation can, without necessarily produce smaller triangles in order to properly decompose an area. The results show that this method properly decomposes any complex, non-convex area with internal holes and provides a two-dimensional point cloud with an associated graph. Nevertheless we have to note that a triangulation produces more waypoints in comparison
7.1 Conclusions 127 with a grid mesh. In order to resolve this issue, some triangles could be omitted, especially in narrow areas with many cells can be observed. In some cases, a combination of a CDT and a square grid decomposition might be more appropriate, since if the area is relatively simple or has large internal uninterrupted areas, its borders could be triangulated and the rest of the space be decomposed in a grid like manner. This would require a proper method to identify those areas and proper algorithms to perform the transition between the different segments. These cases can be tested in a future research effort on the same basis. 7.1.2 Area partitioning for multiple vehicles In Chapter 3 a novel partitioning algorithm has been presented in order to tackle the heterogeneous, multi-robot area allocation problem. The deadlock handling capabilities of this algorithm make it a valuable addition to tackle complex coastal areas and the probability of new restrictions or team reconfiguration during task execution. Moreover, the proposed methods take into consideration the initial positions of the vehicles, proposing in that manner a flexible and robust strategy for a variety of tasks. This is particularly true in remote oceanic missions where the take-off and landing positions of the vehicles are dynamic. These methods are derived by the initial choice of decomposing and treating the whole area as a graph. The cell weights produced in this step which belong to a vehicle, are used later on in the next steps. The results indicate that the goal of producing uniform sub-areas, according to the relative capabilities of the vehicles is accomplished. However, the proposed solution is still an approximation, due to the nature of area discretisation. A continuous solution along with an effort to reduce the complexity of the deadlock handling algorithm is considered as a future direction on that matter. 7.1.3 Coverage Plans In Chapter 4, each of the triangular mesh centroids has been treated as a waypoint and by their prioritisation, full coverage waypoint lists have been produced. These lists are treated as a flight plan composed by the segments connecting consecutive waypoints. The spiral nature of these plans has been a concious decision, having in mind fixed-wing vehicle dynamics, which would permit smoother turns in order to further aid the cause of maintaining the flight path inside the ROI. As described before, this is also a mitigation in the safe separation issues between vehicles of the team or external aerial or ground restrictions. Moreover, an algorithmic technique has been developed in order to identify area characteristics which would affect the overall plan by introducing repeated coverage.