Landing sites detection using LiDAR data on manycore systems

Noname manuscript No. (will be inserted by the editor) Landing sites detection using LiDAR data on manycore systems Oscar G. Lorenzo ·Jorge Mart´ınez · David L. Vilari˜no ·Tom´as F. Pena · Jos´e C. Cabaleiro ·Francisco F. Rivera Received: date / Accepted: date Abstract Helicopters are widely used in emergency situations, where knowing if a geographical location is adequate for landing is a critical issue, and it is far from being a straightforward task. In this work, we present a method to detect and classify landing sites from LiDAR data in parallel on multi- and manycore systems using OpenMP. Load balancing was identified as the main cause of poor performance because the computational cost depends mainly on the input data. Results for a set of LiDAR point clouds that represent different real scenarios were used as case studies in this work. Balancing strategies for three different multi and manycore systems were analyzed. The proposed load balancing techniques increase performance up to three times from the unbalanced case. Keywords LiDAR ·landing zone detection ·load balancing ·Xeon Phi 1 Introduction Helicopters are widely used in emergency situations such as rescue operations or firefighting efforts, sometimes in difficult terrain circumstances. Under these situations, the recognition of suitable landing locations may become difficult even for expert pilots. The automatic detection of suitable landing sites may represent a valuable aid for pilots in environments where visual inspection is insufficient. In this work, we present a method to use LiDAR point clouds [17] for landing sites detection and classification, with time performance as a main issue. A preliminary introduction to this method was already presented in [8]. Oscar G. Lorenzo CiTIUS Centro Singular de Investigaci´on en Tecnolox´ıas da Informaci´on, Universidade de Santiago de Compostela Tel.: +34-881-816400 E-mail: [email protected] 2 Oscar G. Lorenzo et al. LiDAR has become an established method for collecting accurate information of the earth landscape. The record of elevation data and the ability to penetrate through the canopy, are some of the advantages over traditional acquisition methods like aerial imagery. The required equipment consists of an airborne laser scanner (ALS), an inertial measurement unit (IMU) and a GPS receiver. These devices are installed in an aircraft that flies over the ground surface obtaining a georeferenced point cloud. On land a GPS system is synchronized with the GPS receiver on the aircraft. For each point a set of features are recorded: X and Y coordinates, elevation (Z), intensity, pulse id, number of pulses, and scan angle. Due to the huge size of the collected information, parallelization is essential to achieve efficient solutions, specially on manycores as the Intel Xeon Phi [14]. Processing LiDAR information often leads to situations in which execution times depend on the specific input data, as is the case of the algorithm proposed in this paper. Also, load balancing is a recurrent issue when parallelizing point clouds. In this work, we analyze the importance of this unbalance workload, and propose ways to deal with this issue in different platforms, including multi- and manycore systems. The rest of the paper is organized as follows: in Section 2, related work is presented. In Section 3, the strategy for classification of landing zones based on their suitability is introduced. In Section 4, the solution to the landing problem is presented, describing the point cloud preprocessing in 4.1, the classification algorithm in 4.2 and the parallelisation in 4.3. In Section 5, some case studies are presented, including examples of the landing detection, an analysis of execution times in 5.2, and a study of the effect of load balancing strategies in 5.3. Finally, some conclusions are summarized in Section 6. 2 Related work LiDAR information has been proven useful for landing hazard and landing sites detection for different aircraft by using geometrical approaches. A number of conditions have to be satisfied to characterize each LiDAR point as part of a landing zone, for example, the incidence angle of the rays and the roughness of the terrain. In [5], an early LiDAR obstacle detection was used to safely guide the descent of ships on Mars. In [18], it was shown how an unmanned helicopter can detect and avoid obstacles with little human intervention. Hemispherical LiDAR was used to avoid obstacles in flight, and also classify possible landing sites by computing the plane of a set of points and determining their slope and roughness –for obstacle detection–, and also checking for the proximity of hazards and obstructions. In [15], LiDAR is used for autonomous landing of helicopters in real time. The algorithmic approach consists of a coarse evaluation based on slope and roughness criteria, a fine evaluation for skid contact, and body clearance of a location. Additionally the ground path to a goal location and approach paths to the landing site are considered. First, the LiDAR points are processed to create a regular grid for the coarse evaluation of po- LiDAR Landing 3 tential landing areas and, second, a separate triangulation is used for the fine evaluation. For the coarse evaluation, a hash table of grid cells is kept, storing the points, and, if a site is considered adequate, a 2D Delaunay triangulation of all available data is built. In contrast, in our proposal, an octree is used to order the LiDAR points in both phases, speeding up the neighbor search [11] and allowing our algorithm to evaluate each point individually with higher precision. In [9], a 3D Convolutional Neural Network is used to differentiate between vegetation and solid obstacles, in order to output a prediction, which is integrated into a more complex evaluation system. Given the region of interest and a point cloud stream, a volumetric density map is incrementally defined, similar in spirit to occupancy grids [10]. Obstacles are classified using trained 3D Convolutional Neural Networks on subdivisions of this volumetric map. In contrast, our geometrical approach uses the normals of each point to create a probability of it belonging to an obstacle or vegetation, without needing to subdivide the space or to train the algorithm. While most of these algorithms are focused on finding a landing site on flight time, using the data collected from the helicopter, our work also intends to fully classify large LiDAR clouds in terms of the terrain suitability for landing, to help in mission planning. To achieve real time performance processing LiDAR data, parallelization is needed. In [6], the parallelization of the decomposition of the full digital waveforms of the returned LiDAR energy is presented. Building recognition, road detection or, in general, object recognition have been shown to be tasks that can be efficiently parallelized. In [7], the parallelization of a road detection algorithm runs into a load balance issue, which is solved by increasing the grain size of the data partition each process computes. In [2], 3D boundary points for recognizing objects are efficiently executed on distributed memory with MPI; data communications are minimized as the point cloud is stored in a quadtree instead of being shared among remote nodes. In [19], reinforced learning is used to balance the load in hybrid computational systems, where a task scheduler decides the optimal execution path. In [16], an assessment of the particularities of OpenMP programs on the Xeon Phi was conducted, showing the importance of a good load balance. In our approach, a parallel algorithm was developed which allows each point to be treated individually, speeding up usual LiDAR operations, such as neighbor search, by using an octree to store the point cloud, and taking care of issues about load balance. 3 Landing zones classification An automatic tool for landing sites detection must be accurate to help pilots to make decisions in emergency situations or plan missions. This should be done as quick as possible, so it can be used in real time emergency operations. According with professional pilots landing sites should meet the following requirements: 4 Oscar G. Lorenzo et al. –A planar surface where the helicopter skids may land: It must be large enough for the skids, it must have a slope not steeper than a given threshold, and vegetation and other obstacles in the surface must not exceed a given height. –An area free of obstacles large enough for the helicopter rotors and fuselage: obstacles must not exceed a given height. –A larger area devoid of obstacles that could impede a clear approximation to the landing zone and/or an aborted landing: the helicopter must be able to approach the area on a given angle with the terrain. The degree of fulfillness of the aforementioned requirements allows a soft classification of landing sites based on their suitability; for example, to detect zones of a high landing difficulty, but acceptable in an emergency, or zones with some criteria near the threshold, or good landing zones, etc. All these requirements must be expressed in terms of the LiDAR point cloud. According with the pilots specifications, the final classification should be divided into five groups, depending on the landing obstacles, where each point would also have a suitability value due to the terrain around it. These five groups are: unsuitable,planar,risky,suitable or good approximation. Following these specifications, our approach considers the computations in a per-point basis. In this way, each point is classified considering three main features: the type of terrain around it, the available room for the helicopter skids to land, and the size of the area around it devoid of obstacles for the helicopter rotors and fuselage. This way, each point is classified according to the points in its neighborhood, defined, for a LiDAR point cloud, as a voxel centered in the point with a given length, or radius if it is circular. The identification of members of these neighborhoods is costly in terms of computing performance. In the classification process, the neighborhood of every point is successively increased and analyzed. However, if at any time, due to the terrain classification or the detection of landing obstacles, a point can be definitely classified, then the remaining neighborhoods do not need to be computed. This feature reduces the execution time. For example, for an unsuitable point in a forested area, the terrain classification would discard it as a landing zone, so the following larger neighborhoods would not be searched for obstacles. While this approach improves performance, it is a source of load imbalance in the parallel code. 4 Parallel landing zones detection In a LiDAR point cloud, each point is characterized by its 3D coordinates, and some attributes. Nevertheless, they are not ordered and follow no predefined structure, so, in order to analyze a feature or a geometry, each 3D point needs to be checked. This makes the process of locating neighbors a very costly one in terms of performance. In order to process the point cloud efficiently, both ordering and storing them in a suitable data structure are needed. This is LiDAR Landing 5 performed in a preprocessing step. Afterwards the classification is done by the Landing Site Detection Algorithm (LSDA). 4.1 Point cloud preprocessing Initially, the points are sorted according first to their X coordinates, second to their Y coordinates, using the quick sort algorithm, and saved to an array. While this array is an useful data structure for many operations, it is not adequate for finding the neighbors of a point. Therefore, we use a read-only octree structure to store the point cloud. Octrees are well-known tree structures in which each node has exactly eight children named octant. The space is partitioned recursively until the size of a octant is below a given threshold (1.5 meters in this work). The use of this octree structure speeds up considerably the search of neighbors of a point [11]. As a second stage in the preprocessing, normal vectors are computed for each point. The normal vector to a surface, simply called the “normal”, is a vector perpendicular to the surface at a given point. To compute the normal of each point, a small neighborhood around it, inside a voxel with a given side, is searched in the octree to define the plane using Principal Component Analysis (PCA) algorithm [13]. The normals to all points are computed in parallel, using the point array to guide it. Nevertheless, obtaining each point neighborhood takes a different execution time, because points may not be evenly distributed and the octree is not balanced. Points with more neighbors traverse more of the octree and do more computations on the PCA than those with fewer ones. Also, even when the number of neighbors is similar, the points spacial distribution affects the shape of the octree. When points are regularly spaced, for example in a plain area, with low Z variation, the corresponding octree is more regular and neighbors are placed in closer leaves, whereas when point spacing is irregular, for example in forested areas with higher Z variation, neighbors are placed in leaves farther apart. This fact implies that searching for the same number of neighbors may take different times, depending on the point cloud. These two factors create a load imbalance in the normals calculation. In Figure 1, an example of an imbalanced 2D quadtree is presented (for the sake of clarity, the figure shows a 2D quadtree, but similar discussion is valid for a 3D octree). At the top of the figure the point distribution in X and Z in a neighborhood is shown, with its resulting quadtree at the bottom; the points placed in each small square go to one single leaf. In this figure, thread 0 is gathering the neighborhood of a point in a plain zone, where there are 14 neighbors close in the quadtree, in two leaves. The point whose normal thread 1 is computing corresponds to a zone with low Z variation, where there are still 14 neighbors, but placed with less locality in the quadtree, in different leaves; due to this distribution it would take longer to access the neighbors. Finally, thread N−1 is calculating the normal of a point with fewer neighbors, due to the high variation on Z, which makes the points in the LiDAR cloud 6 Oscar G. Lorenzo et al. Fig. 1 Example showing imbalances in accessing to the quadtree for normals calculation. farther apart in 2D space, so gathering them and performing PCA would take less time. To deal with this imbalance, in this work, the parallelization is done dynamically, dividing the work in blocks of a same number of points and allocating to each thread a new block as soon as it finished with the previous one. Once the point cloud is stored in the octree and the normals have been computed and added as a new attribute of the points, the landing algorithm is ready to classify the points, considering each one, individually, by turn. 4.2 Landing Site Detection Algorithm Our Landing Site Detection Algorithm (LSDA) is composed of several stages, which are depicted in Figure 2. These stages classify each LiDAR point in one of the groups defined in Section 3. In the first stage of the LSDA, the suitability of the terrain around the considered point is evaluated. A planar surface is calculated centered in the point using again PCA, with a neighborhood just large enough to accommodate the landing skids of the helicopter. This would be the surface where it would land, so, according to the requirements for skid landing in terms of slope and vegetation, a quality value is assigned to the terrain. A minimum number of points in each surface can be set to avoid evaluating zones with insufficient information. If, at any moment, the terrain does not fit any of the conditions, the LSDA stops for this point, which is then classified as unsuitable. The slope of the plane is checked first, if it is steeper than a given value, the point is immediately classified as unsuitable, otherwise its steepness is considered for quantifying its quality. Then the dispersion for the normal vectors of all LiDAR Landing 7 Fig. 2 Landing Site Detection Algorithm flow diagram the points in the plane is calculated. If it is higher than a certain threshold, the terrain is considered too rough, and the point is classified as unsuitable, otherwise its value is used to determine if the terrain is rough vegetation –such as grass, where the dispersion is higher– or a plain surface –such as a road, where the dispersion is lower–. Next, the distance of each point in the neighborhood to the plane is used to look for obstacles that may impede the skids to land. If any neighbor point is considered an obstacle, the center point is classified as unsuitable. Due to the discrete nature of the sensor data, the height of some objects might be poorly estimated. Therefore, with the aim of improving the robustness of 8 Oscar G. Lorenzo et al. Fig. 3 Example of obstacle detection the algorithm, the normal vector of the neighborhood point is used to give a confidence value of the height of the obstacle. Consider the example in Figure 3. The blue point is the center point to be considered, and the rest are points inside its neighborhood. If only the height of the neighborhood points were measured, the red point would be considered an obstacle and both the green and yellow points, which have the same height, would not. Nevertheless, if the red point were not in the LiDAR scan data, the yellow point would be a clue to the presence of an obstacle. For this reason, a confidence value to the height is given according to the angle between a point normal and the plane normal. If this angle is small, as with the green point, whose normal is parallel to the ground plane normal, its height is considered to be representative of the highest point of an obstacle or the terrain. If the point normal is in angle with the terrain, there is a possibility it is placed in a local slope, which may be representative of an obstacle, as with the yellow point in Figure 3, so the algorithm detects a possible obstacle, with a certain confidence. This confidence depends also on the underlying terrain –vegetation or plain–, since the chance of normals being in angle with the plane normal is different in each case. With all these values, a final quality assessment of the terrain is defined. If the point was not classified as unsuitable in this step, the LDSA continues searching if the larger area around it can accommodate the rotors and fuselage of the helicopter. In the second stage of the LDSA, successively larger neighborhoods are searched for obstacles that could collide with the rotors or fuselage. Since the larger the neighborhood is, the longer it takes to compute, the algorithm starts considering small neighborhoods, corresponding to more restricted landing zones, with less strict requirements. If obstacles are found, the points are classified accordingly and the larger neighborhoods are not searched. In this way, a new neighborhood is searched, with radius equal to the rotor and fuselage, which is the minimum possible landing site. According to the points in this larger area, the point can either be classified as planar –when the skids may land, but not enough space for the rotors is accounted for–, or a larger area is searched to further decide the quality of the landing site –when there is enough space for the whole helicopter–. LiDAR Landing 9 Fig. 4 Approximation to target point Afterwards, the quality of the landing site is determined. To check if the landing site is suitable, a new neighborhood, this time considering a radius equal to the rotor and fuselage, plus a safety distance, must be searched. Points already considered in the previous neighborhood are not considered now. The distances of all these neighborhood points to the plane determine if the landing site is risky or suitable. Finally, possible obstacles for a safe approximation of the helicopter are detected. A new and larger neighborhood around the suitable landing is gathered. All points that fall inside an inverted cone, with apex in the landing point and an aperture dependent of the desired approximation angle, are considered obstacles (see Figure 4). If there are no obstacles, the point is classified as good approximation. The points in the neighborhoods being studied are classified temporally as belonging to, or accompanying, a planar,risky,suitable or good approximation site. This temporal classification may be overridden if they become the center of a landing site. The final classification can be further refined according to the terrain quality calculated on the former steps. In this way, each point is fully classified in a single iteration. Note that there is no need for computing the whole point cloud to find a good enough solution depending on the input data. 4.3 Parallelization and load balancing The process of gathering all these neighborhoods means that each point needs a different execution time to be classified, since if a point is detected to be unsuitable in one of the stages of the algorithm, it is not considered in next stages with larger neighborhoods. Also note that a planar point does not require calculating the last neighborhood. In addition, the number of points inside these neighborhoods may vary greatly. These facts are the cause of a load balancing issue when the LDSA is parallelized, even while all points are processed independently. Nevertheless, this algorithm, as has been previously stated, allows for any landing point to be classified in just one sweep, which is useful when working with streaming data. 16 Oscar G. Lorenzo et al. Table 2 Throughput per Watt for each system‘s processor. System Thermal Design Power Mean p/s Mean points/(s*W) (W) (points/joule) Core 95 33,299.1 350.5 Phi 300 79,407.9 264.7 Xeon 130 114,746.5 850,0 the same time the points are obtained, with the I/O being performed concurrently with the classification. In this case, for the program to work alongside the sensor, the points would have to be classified faster than the ALS may measure them. The p/s an ALS may obtain varies depending on different factors, such as flight altitude and speed, desired resolution or sensor quality [3], and may be configured by the user. The classification program, executed on the Xeon, attains between 62,000 p/s and 256,000 p/s. This rate is good enough for practical applications, considering the robustness of the algorithm, that detects landing zones on low point density maps if necessary. In terms of energy usage, an important issue in case the system has to be placed in the helicopter, the Phi is the least energy efficient system, as can be seen in table 2, where the throughput per Watt of each processor is stated. 5.3 Load balance Load balance was studied varying the block size parameter. The case when all threads execute just 1 block is performed using OpenMP without scheduling clauses, so no scheduling is used and the work is divided beforehand. To study load balance is useful to consider the number of blocks the work is divided into, for each thread, named blocks per thread. Therefore, blocks per thread =number of points block size ×number of threads (1) The results for the Core system with 4 threads are shown in Figure 11. These results show that, once the work is split in more than 20 blocks per thread, with block size 16346, the workload balances, although the best results are with a few thousand blocks. When the work is not balanced, execution times are often about 2 times higher than well balanced cases. Execution times for the computation of normals are shown in Figures 11(a) and 11(b). Note that the larger the point density the longer execution time, since more points are searched in each neighborhood, except in the case of the Guitiriz map. This is because the point distribution is influenced by the kind of Li- DAR sensor used, and the way the Guitiriz map was created may explain this behavior, because it considered some returns as differentiated points, causing to artificially increase its density. In Figures 11(c) and 11(d), which show the execution time of LDSA, the percentage of risky and, specially, of suitable points greatly influence the LiDAR Landing 17 (a) Normals 1,700,000 points. (b) Normals 5,000,000 points. (c) LDSA 1,700,000 points. (d) LDSA 5,000,000 points. Fig. 11 Execution times on the Intel Core with 4 threads load imbalance. This is because LDSA performs more computation for risky and suitable points, including looking for larger neighborhoods. For this reason, Vaihingen, Australia and Florida, the ones with a higher percentage of this class of points, are on top. Nevertheless, Florida takes longer than Australia, although it has fewer risky and suitable points, due to its higher point density, which affects the neighbor search. In Figure 12 the results for the Phi system are shown. In this system, the work balances from about 5 blocks per thread, around 980 block size. When the workload is unbalanced, results are from 2 to 3 times worse, due to the importance of a good parallelization in the manycore architecture. The load balance effects vary with the map, being more acute on the Vaihingen and Australia cases than on the Florida one. Results for the Xeon system can be seen in Figure 13. In this case the load balances from around 15 blocks per thread, a block size of 4048, and the gain against the unbalance case is about 2.5×. Taking into account the best load balanced cases, with block size of 64, the performance improvement when increasing the thread count per core, for the Phi system –with baseline in one thread per core and all cores used– is shown in Figure 14(a). The performance improves as more threads are added per core, although it does not reach a 4×improvement. In Figure 14(b) the speedup for the Xeon system is shown –with baseline in one thread and one core–. Again, the speedup does not reach the optimal of 12, but remains high, taking advantage in this case of the multithreading, 18 Oscar G. Lorenzo et al. (a) Normals 1,700,000 points. (b) Normals 5,000,000 points. (c) LDSA 1,700,000 points. (d) LDSA 5,000,000 points. Fig. 12 Execution times for large maps on the Xeon Phi with 244 threads (a) Normals 1,700,000 points. (b) Normals 5,000,000 points. (c) LDSA 1,700,000 points. (d) LDSA 5,000,000 points. Fig. 13 Execution times for large maps on the Xeon with 12 threads LiDAR Landing 19 (a) Xeon Phi. (b) Xeon. Fig. 14 Execution time varying thread number using the best block size when there are 2 threads per core. For the Guitiriz and Alcoy maps the speedup actually decreases on the Phi with 244 threads, but still increases on the Xeon with 12 threads. This is likely related to the high percentage of unsuitable points in these maps–88.28 % and 84.53 %, see Table 1–. To classify a point as unsuitable takes the least amount of computation, so these cases do not benefit from increased parallelism. 6 Conclusions In emergency situations, knowing whether a geographical position is adequate for a helicopter to land is a crucial task. In this paper, we propose a technique to classify the suitability of different locations for helicopter landing from LiDAR point clouds. In order to be used in practical situations, the approach has to be robust and fast. To this end, parallel algorithms for automatic landing sites detection are introduced. The reliability of the proposal has been tested on a large set of data with different terrain characteristics. Computation speed has been optimized by parallelizing the algorithms, making the proposal suitable for real-time applications. Load balance was detected as the main obstacle to achieve good performance and it was characterized. Three different computer systems were considered with 4, 6 and 61 processing cores, respectively. The load balancing techniques applied increase performance up to 3 times from the unbalanced case. Acknowledgements This work has been partially supported by the Ministry of Ministry of Economy and Competitiveness of Spain under project TIN2013-41129-P and Xunta de Galicia under projects GRC2014/008 and GRC GI-1638. It has been developed in the framework of the European network HiPEAC-2, the Spanish network CAPAP-H, the Galician network under the Consolidation Program of Competitive Research Units (TLIX Network ref. R2014/049). This work is also a result of a collaboration with INAER and the Laborate group of the University of Santiago de Compostela. 20 Oscar G. Lorenzo et al. References 1. Dagum L, Enon R (1998) OpenMP: an industry standard API for sharedmemory programming. Computational Science & Engineering, IEEE 5(1):46–55 2. Gong J, Lee HJ, Cho G (2008) Parallel algorithm for object extraction from ranging data. URL http://dcs.chonbuk.ac.kr/~jgong/papers/ Parallel%20Algorithm%20for%20Object%20Extraction%20from% 20Ranging%20Data.pdf, technical Report, online 3. Hodgson ME, Bresnahan P (2004) Accuracy of airborne lidar-derived elevation. Photogrammetric Engineering & Remote Sensing 70(3):331–339 4. ISPRS (2016) ISPRS test project on urban classification, 3D building reconstruction and semantic labeling. URL http://www2.isprs.org/ commissions/comm3/wg4/tests.html 5. Johnson AE, Klumpp AR, Collier JB, Wolf AA (2002) LiDAR-based hazard avoidance for safe landing on Mars. Journal of guidance, control, and dynamics 25(6):1091–1099 6. Jung JH, Crawford MM, Lee SH (2009) Complexity estimation based work load balancing for a parallel LiDAR waveform decomposition algorithm. Korean Journal of Remote Sensing 25(6):547–557 7. Li J, Lee HJ, Cho GS (2008) Parallel algorithm for road points extraction from massive LiDAR data. In: Parallel and Distributed Processing with Applications, 2008. ISPA’08. International Symposium on, IEEE, pp 308– 315 8. Lorenzo OG, Mart´ınez J, Vilari˜no DL, Pena TF, Cabaleiro JC, Rivera FF (2016) Parallel landing sites detection using LiDAR data on manycore systems. In: 16th International Conference on Computational and Mathematical Methods in Science and Engineering, CMMSE, pp 805–816 9. Maturana D, Scherer S (2015) 3D convolutional neural networks for landing zone detection from LiDAR. In: Robotics and Automation (ICRA), 2015 IEEE International Conference on, IEEE, pp 3471–3478 10. Moravec H, Elfes A (1985) High resolution maps from wide angle sonar. In: Robotics and Automation. Proceedings. 1985 IEEE International Conference on, IEEE, vol 2, pp 116–121 11. Mosa ASM, Sch¨on B, Bertolotto M, Laefer DF (2012) Evaluating the benefits of octree-based indexing for LiDAR data. Photogrammetric Engineering & Remote Sensing 78(9):927–934 12. OpenTopography (2016) LiDAR point cloud data distribution and processing. URL http://opentopo.sdsc.edu/lidar 13. Peason K (1901) On lines and planes of closest fit to systems of point in space. Philosophical Magazine 2:559–572 14. Ramachandran A, Vienne J, Van Der Wijngaart R, Koesterke L, Sharapov I (2013) Performance evaluation of NAS parallel benchmarks on intel Xeon Phi. In: Parallel Processing (ICPP), 2013 42nd International Conference on, IEEE, pp 736–743 LiDAR Landing 21 15. Scherer S, Chamberlain L, Singh S (2012) Autonomous landing at unprepared sites by a full-scale helicopter. Robotics and Autonomous Systems 60(12):1545–1562 16. Schmidl D, Cramer T, Wienke S, Terboven C, M¨uller MS (2013) Assessing the performance of openmp programs on the intel xeon phi. In: European Conference on Parallel Processing, Springer, pp 547–558 17. Shan J, Toth CK (2008) Topographic laser ranging and scanning: principles and processing. CRC press 18. Whalley M, Takahashi M, Tsenkov P, Schulein G, Goerzen C (2009) Fieldtesting of a helicopter UAV obstacle field navigation and landing system. In: 65th Annual Forum of the American Helicopter Society, Grapevine, TX 19. Zhang T, Li J (2015) Online task scheduling for LiDAR data preprocessing on hybrid GPU/CPU devices: A reinforcement learning approach. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing 8(1):386–397