Geometric Targets for UAS Lidar

: Lidar from small unoccupied aerial systems (UAS) is a viable method for collecting geospatial data associated with a wide variety of applications. Point clouds from UAS lidar require a means for accuracy assessment, calibration, and adjustment. In order to carry out these procedures, speciﬁc locations within the point cloud must be precisely found. To do this, artiﬁcial targets may be used for rural settings, or anywhere there is a lack of identiﬁable and measurable features in the scene. This paper presents the design of lidar targets for precise location based on geometric structure. The targets and associated mensuration algorithm were tested in two scenarios to investigate their performance under di ﬀ erent point densities, and di ﬀ erent levels of algorithmic rigor. The results show that the targets can be accurately located within point clouds from typical scanning parameters to < 2 cm σ , and that including observation weights in the algorithm based on propagated point position uncertainty leads to more accurate results.


Introduction
Airborne lidar, or laser scanning, from small unoccupied aerial systems (UAS) has gained a reputation as a viable mapping tool for both academic researchers and commercial users within the past decade (e.g., [1][2][3]). Accompanying the release of the lightweight Velodyne VLP-16 in 2014 was the emergence of highly-accurate, lightweight navigational hardware, eliminating the need for auxiliary sensor data, such as computer vision techniques used in the earliest UAS lidar studies [4]. These two core components form a payload that is functionally similar to the conventional airborne lidar payloads aboard piloted aircraft, but at a fraction of the weight and cost. Commercial outfits (e.g., LiDARUSA, YellowScan, Phoenix LiDAR Systems) have offered turnkey UAS lidar systems which utilize the VLP-16 since at least 2015 [2]. From an accessible platform, the UAS, users can now collect an established and familiar data product, the lidar point cloud, which is compatible with existing data processing workflows and analysis.
The estimated positioning accuracy of the elevated 2 m diameter targets was around 5-10 cm and 2-3 cm for the horizontal and vertical components, respectively, at a lidar point density of 5 pts/m 2 . The study in [12] used white L-shaped targets painted on asphalt to similarly locate 3D positions within lidar-derived intensity images using least-squares matching [13]. The study in [14] used retro-reflective targets in a hexagonal configuration to find a location within a relatively sparse point cloud (approximately 2 pts/m 2 ) to aid in boresight calibration, with reported location precision of 5 cm and 4 cm in the vertical and horizontal components, respectively. The dimensions of a proposed, but not yet implemented, foldable 3D target were illustrated by [15], designed to use solely the point positions, and not intensity.
The current applied methods for targets in UAS lidar primarily employ bright flat targets which make use of the intensity value of lidar returns available in most scanners. For example, [16] used ICP to align UAS lidar swaths, ground control points (GCPs) marked with flat sheet targets to apply an absolute correction, and targeted and untargeted ground-surveyed checkpoints for accuracy analysis. The method in [17] used the average location of points striking reflective ground targets for accuracy assessment. That study described critical issues with this approach, including assumed distribution of points striking the target and effects of the laser beam divergence leading to bias and/or lower accuracy of measured locations in the point cloud compared to expectations from stochastic modeling. Similarly, [17] used reflective ground targets manually measured in UAS point clouds, but details of the method were not described. There is not much exploration of mensuration uncertainty, the uncertainty of target location within the point cloud, in the literature, and thus a direct comparison of the various targeting methods is not possible.
Often, the use of intensity is not a robust solution for artificial targeting in laser scanner data. For example, the Velodyne VLP-16 is one of the most commonly-used UAS laser scanning systems. As reported by [18] and in the present authors' experience, the VLP-16 intensity reading quality is highly dependent on range and the incidence angle of the beam with the intercepting object. Additionally, intensity is recorded at a lower quantization level than larger and more-expensive laser scanners [19]. At the low altitudes associated with UAS platforms, incidence angles for small homogenous areas and objects in the scene can have very diverse incidence angles, and therefore objects with similar reflectivity can yield substantially diverse intensity readings when surveyed in the typical boustrophedonic overlapping-swath flight pattern. It is often difficult, for example, to locate points from overlapping swaths associated with a white paper plate on healthy grass based on intensity from the VLP-16 as illustrated in Materials and Methods section. Use of very bright objects may be one approach to addressing this issue. However, [18] indicates that high-intensity VLP-16 readings are associated with range bias on the order of a centimeter, which has also been observed by the present authors. In addition to issues associated with locating target points by intensity, the flatness of targets can also be problematic. As [4] mentioned, non-uniform point distribution on the target can lead to poor, biased resolution of the horizontal reference point of the target, often located in the center of extreme target points, e.g., x = (max x − min x )/2. This is the case for both automatic methods and manual methods, where an operator estimates the center location of the target and interpolates the elevation of the target center based on surrounding points. To address these issues, this study introduces an alternative to intensity-based artificial targets for UAS lidar. More specifically, this study presents: 1.
The description of a developed small-UAS laser scanning system; 2.
The design of geometric structure-based targets which can be used for accuracy assessment, calibration, and strip adjustment; 3.
Mensuration algorithms for precise target location within point clouds; 4. Two case studies illustrating the efficacy of the targets and associated algorithms for accuracy assessment.
While the accuracy results of the UAS system are presented and are of interest, they are not the focus of the paper. Instead, we focus on the methods to achieve the accuracy assessment, which can Remote Sens. 2019, 11, 3019 4 of 20 be used for other similar UAS laser scanning systems. It should be noted that UAS lidar accuracy may vary depending on site location, control configuration and accuracy, system characteristics, and mission planning.

Materials and Methods
The UAS laser scanning system used in this study ( Figure 1) comprised a DJI s1000 octocopter with a sensor payload consisting of a Velodyne VLP-16 Puck Lite laser scanner head, NovAtel OEM615 Global Navigation Satellite System (GNSS) receiver, STIM300 IMU, Garmin GPS18x-for time stamping outgoing lidar pulses, a NovAtel GPS-702-GG antenna (Experiment Site 1 only), and a Maxtena M1227HCT-A2-SMA antenna (Experiment Site 2 only). Although this system was designed and built by our team, it can be considered typical, with various commercially-available turnkey analogues. Similarly, processing algorithms to generate point clouds from raw laser scanner and navigation sensor data, calibrate the sensor, and perform strip adjustment were developed internally. They can also be considered standard, except that they were designed to use our developed targets, and had the capability to propagate and carry forward uncertainties from the navigation and scanner observations to the point cloud. This allowed the estimated standard deviations of the X, Y, and Z components of each point in the point cloud to be stored and accessed, a benefit that we explored in our target mensuration experiments. The UAS laser scanning system used in this study ( Figure 1) comprised a DJI s1000 octocopter with a sensor payload consisting of a Velodyne VLP-16 Puck Lite laser scanner head, NovAtel OEM615 Global Navigation Satellite System (GNSS) receiver, STIM300 IMU, Garmin GPS18x-for time stamping outgoing lidar pulses, a NovAtel GPS-702-GG antenna (Experiment Site 1 only), and a Maxtena M1227HCT-A2-SMA antenna (Experiment Site 2 only). Although this system was designed and built by our team, it can be considered typical, with various commercially-available turnkey analogues. Similarly, processing algorithms to generate point clouds from raw laser scanner and navigation sensor data, calibrate the sensor, and perform strip adjustment were developed internally. They can also be considered standard, except that they were designed to use our developed targets, and had the capability to propagate and carry forward uncertainties from the navigation and scanner observations to the point cloud. This allowed the estimated standard deviations of the X, Y, and Z components of each point in the point cloud to be stored and accessed, a benefit that we explored in our target mensuration experiments. A primary objective of this study was to create portable, inexpensive, and multimodal (measurable in both lidar point clouds and imagery) targets capable of establishing accurate 3D locations within a point cloud from associated point observations with practical densities and distributions associated with UAS-lidar systems. The resulting design is shown in Figure 2. The targets were foldable, corner-cube trilateral pyramids made from white corrugated plastic. The base of the pyramid was 1.1 m and the peak was 0.4 m above the ground when deployed. The reference location, to where ground-survey and point cloud coordinates were measured, was where the three planes intersected at the peak. The targets' foldability allowed for easy transport and storage. Figure  2b shows a stack of 22 targets, which could be easily carried by a single person. Black tape was placed on the pyramids to allow for measurement of the reference point, the intersection of the centerline of the three black lines, in photography taken from any perspective showing at least two sides of the target. A primary objective of this study was to create portable, inexpensive, and multimodal (measurable in both lidar point clouds and imagery) targets capable of establishing accurate 3D locations within a point cloud from associated point observations with practical densities and distributions associated with UAS-lidar systems. The resulting design is shown in Figure 2. The targets were foldable, corner-cube trilateral pyramids made from white corrugated plastic. The base of the pyramid was 1.1 m and the peak was 0.4 m above the ground when deployed. The reference location, to where ground-survey and point cloud coordinates were measured, was where the three planes intersected at the peak. The targets' foldability allowed for easy transport and storage. Figure 2b shows a stack of 22 targets, which could be easily carried by a single person. Black tape was placed on the pyramids to allow for measurement of the reference point, the intersection of the centerline of the three black lines, in photography taken from any perspective showing at least two sides of the target. The mensuration algorithm for locating targets in the lidar point cloud was an iterative leastsquares process to fit a known target template, three intersecting planes based on the design dimensions of the target, to points associated with it using rigid-body six-parameter transformation (3D rotation followed by translation) [20]. The template's untransformed coordinate system origin was located at the reference point at the tip of the pyramid, the base was parallel with the horizontal plane (close to what will likely be the fitted pose, hastening convergence), and one of the edges was aligned such that y = 0 ( Figure 3). This edge orientation was arbitrarily chosen, with any orientation leading to convergence of the mensuration algorithm, meaning that actual targets need not be aligned with any particular orientation when placed in the field. The first step was to identify and isolate points associated with the target. In this study, segmentation of target points was done manually, based on the relative height of the points above the ground. Automatic extraction of these points was feasible by using a random sample consensus (RANSAC) [21] variant of the mensuration algorithm, which is proposed for future work. Nevertheless, manual segmentation was straightforward due to the structure of the targets. Figures  4 and 5 illustrate the identifiability of target points based on height compared to that based on intensity. The mensuration algorithm presented here allowed the use of any subset of the points on the target, as long as some points on each facet were found. The distribution of the points did not bias the estimated target location, a quality that is explored in the appendix. This allowed the user to be conservative in the confident identification of target points, reducing the possibility of blundered selection of non-target points. The mensuration algorithm for locating targets in the lidar point cloud was an iterative least-squares process to fit a known target template, three intersecting planes based on the design dimensions of the target, to points associated with it using rigid-body six-parameter transformation (3D rotation followed by translation) [20]. The template's untransformed coordinate system origin was located at the reference point at the tip of the pyramid, the base was parallel with the horizontal plane (close to what will likely be the fitted pose, hastening convergence), and one of the edges was aligned such that y = 0 ( Figure 3). This edge orientation was arbitrarily chosen, with any orientation leading to convergence of the mensuration algorithm, meaning that actual targets need not be aligned with any particular orientation when placed in the field. The mensuration algorithm for locating targets in the lidar point cloud was an iterative leastsquares process to fit a known target template, three intersecting planes based on the design dimensions of the target, to points associated with it using rigid-body six-parameter transformation (3D rotation followed by translation) [20]. The template's untransformed coordinate system origin was located at the reference point at the tip of the pyramid, the base was parallel with the horizontal plane (close to what will likely be the fitted pose, hastening convergence), and one of the edges was aligned such that y = 0 ( Figure 3). This edge orientation was arbitrarily chosen, with any orientation leading to convergence of the mensuration algorithm, meaning that actual targets need not be aligned with any particular orientation when placed in the field. The first step was to identify and isolate points associated with the target. In this study, segmentation of target points was done manually, based on the relative height of the points above the ground. Automatic extraction of these points was feasible by using a random sample consensus (RANSAC) [21] variant of the mensuration algorithm, which is proposed for future work. Nevertheless, manual segmentation was straightforward due to the structure of the targets. Figures  4 and 5 illustrate the identifiability of target points based on height compared to that based on intensity. The mensuration algorithm presented here allowed the use of any subset of the points on the target, as long as some points on each facet were found. The distribution of the points did not bias the estimated target location, a quality that is explored in the appendix. This allowed the user to be conservative in the confident identification of target points, reducing the possibility of blundered selection of non-target points. The first step was to identify and isolate points associated with the target. In this study, segmentation of target points was done manually, based on the relative height of the points above the ground. Automatic extraction of these points was feasible by using a random sample consensus (RANSAC) [21] variant of the mensuration algorithm, which is proposed for future work. Nevertheless, manual segmentation was straightforward due to the structure of the targets. Figures 4 and 5 illustrate the identifiability of target points based on height compared to that based on intensity. The mensuration algorithm presented here allowed the use of any subset of the points on the target, as long as some points on each facet were found. The distribution of the points did not bias the estimated target location, a quality that is explored in the Appendix A. This allowed the user to be conservative in the confident identification of target points, reducing the possibility of blundered selection of non-target points. Remote Sens. 2019, 11, x FOR PEER REVIEW 6 of 20   The target transformation parameters were initialized by setting the orientation angles equal to zero: ω = 0, ϕ = 0, and κ = 0, where ω, ϕ,. κ are sequential rotations about the target X, Y, and Z axes, respectively. Translation T = (T X , T Y , T Z ) T was initialized to be equal to the coordinates of the highest of the target points.
Each iteration of the least-squares solution begins with identifying the closest facet, f min , of the template in its current pose for each point. This is accomplished by finding the minimum dot product of the unit normal vectors of each facet f , is the position of the point p. Jacobian matrix J. is formed (Equation (1), where m is the number of points), containing the partial derivatives of the distance, d i , of each point i to its associated plane with respect to the six transformation parameters: Finally, a vector containing the negatives of the point-to-plane distances is formed, ε = −d, enabling the least-squares solution to the linearized observation equations using Equation (2): In Equation (2), ∆ contains corrections to the current approximations of the transformation parameters. The process is repeated until convergence of the standard deviation of unit weight, with the final resolved translation components, T, serving as the estimated position of the reference point.
If propagated uncertainties for point coordinates are available, the algorithm has the same steps, except that the solution is: In Equation (3), W = AΣA T −1 , where A is the Jacobean matrix with rows containing the partial derivatives of the distances of each point, d p , to their associated planes with respect to x p : Σ is the covariance matrix associated with lidar point coordinate uncertainty. The effect of including this weighting in the mensuration algorithm was explored in the experiments for Site 2.
The covariance matrix of the target template transformation parameters, Σ ∆∆ , may be calculated using Equation (5). The elements of Σ ∆∆ . associated with the translation components, T X , T Y , T Z , can be used as estimates of the target mensuration uncertainty: Note that in this study, point position uncertainties were based on propagated variance from navigation parameters and scanner-head hardware. Inter-state correlations among these components were not resolved (i.e., the weight matrix associated with point-to-plane distances, W, was a diagonal matrix), and thus it was expected that propagated uncertainties of the target positions would be optimistic, since the assumed independence of point positions, and therefore point-to-plane distances, was violated. Obtaining estimates of the inter and intra-state correlation among navigation and scanner head parameters was not readily achieved, as they are not fully available in commercial navigation processing software and from scanner hardware vendors, respectively. Still, the outcomes from this Remote Sens. 2019, 11, 3019 8 of 20 study presented in the Results and Analysis section suggest that the propagated uncertainties used, even without accounting for correlation among parameters, were valuable.
Since lidar points near target edges may have convoluted echoes and lead to spurious observations, each point is tested for edge proximity. If a point is closer to an edge than a selected tolerance it is removed from the solution. In this study, points were culled if they were closer than one fourth of the beam diameter to an edge, a tolerance found empirically to improve the fidelity of the observations without overly decreasing redundancy. Similarly, points were removed if their orthogonal projection onto their associated facet-planes fell outside of the facet triangle. Points were tested in each iteration of the solution after the first two iterations were completed (the algorithm requires at least three iterations before convergence is established). The distance of a lidar point p to a target edge e i,j associated with The projection of a point p onto the plane containing target vertices v i , v j , and v k . associated with the closest target facet f min is: If x p is inside of the triangle formed by v i , v j , and v k , the barycentric coordinates of x p with respect to the triangle will all be positive. This condition is used to test if the orthogonal projections of points fall on their associated target facets.
Two separate study areas were used to illustrate the efficacy of the pyramidal target methods. The first (Site 1) was used to illustrate precision and how lidar point density affects target mensuration within the point cloud, and the other (Site 2) to illustrate the benefit of using error propagation to locate targets.
Site 1 was located in Jonesville, Florida, at a recreation/sports park. Site 2 was located within Ordway Swisher Biological Station (OSBS), a 9500-acre University of Florida research facility in Melrose, Florida, which is primarily used for ecosystem monitoring and management research. The planned mission parameters are shown in Table 1. Site 1 consisted of three north-to-south flight paths ( Figure 6). Data were collected along these flight paths and associated turns twice at 40 m and once at 30 m flying height above ground level (AGL). The area was scanned with the full outgoing pulse rate of the VLP-16, 300 kHz, yielding an average combined-swath pulse density of 3000 pts/m 2 , or an approximate point spacing of 1.8 cm in the center of the project area. The scan rate was 10 Hz, although it should be noted that later investigation found maximizing the scan rate of the VLP-16 yielded better point distribution, reducing clustering of points, albeit with lowered point accuracy due to lower resolution of the encoded scan angles. In order to investigate the effect of point density on target mensuration, the original processed point cloud was decimated three times to generate target point clouds with 5 cm, 10 cm, and 15 cm minimum 3D point distances, by iteratively and randomly selecting points and removing all neighbors within less than the minimum 3D point distance. Initially, 20 cm point spacing was also attempted, but it sometimes (2 of 20 targets) led to failure of the target mensuration algorithm and was excluded from analysis. the center of the project area. The scan rate was 10 Hz, although it should be noted that later investigation found maximizing the scan rate of the VLP-16 yielded better point distribution, reducing clustering of points, albeit with lowered point accuracy due to lower resolution of the encoded scan angles. In order to investigate the effect of point density on target mensuration, the original processed point cloud was decimated three times to generate target point clouds with 5 cm, 10 cm, and 15 cm minimum 3D point distances, by iteratively and randomly selecting points and removing all neighbors within less than the minimum 3D point distance. Initially, 20 cm point spacing was also attempted, but it sometimes (2 of 20 targets) led to failure of the target mensuration algorithm and was excluded from analysis.   (Figure 7). The area was scanned at 30 m flying height above ground, with a pulse rate of 300 kHz, and 10 Hz scan rate. The average combined point cloud pulse density was approximately 800 pts/m 2 , or a point spacing of approximately 3.5 cm. Note that since the Site 2 flight paths were longer, the increase in the average point density associated with points collected on turns between primary lines was smaller than that for Site 1. Similarly, the longer lines for this flight led to higher uncertainty for some target-associated points, since uncertainty is largely dependent on laser range, making this dataset well-suited for the analysis of uncertainty-based weighting on target mensuration. To study this, the target locations were estimated using both the weighted and unweighted versions of the mensuration algorithm.
Boresight calibration of the laser scanner sensor suite, estimation of the angular alignment of the scanner head with respect to the platform axes corresponding to the filtered trajectory (IMU axes), was carried out for both missions, each time using three targets without using ground survey coordinate observations. The details of the boresight procedure are beyond the scope of this paper, but details will be given in a future article. Strip adjustment, alignment of the overlapping swaths, was not necessary to bring collection units into alignment for this study, since there was no observable discrepancy between flight strips after the boresight calibration procedure.
At each site, the targets were used to mark ground-surveyed checkpoints, which were well-distributed across the project area. There were 20 checkpoints within Site 1 that were nominally spaced at 7 m along three north-south rows. The eastern-most row was approximately 14 m from the central row, and the western-most row was 20 m from the central row. These can be seen in Figure 6. Site 1 target coordinates were established using post-processed kinematic (PPK) GNSS. The rover GNSS antenna was held above the targets so that the processed coordinates could be reduced to the target reference points. It is important to mention that the target survey used the same base station that was used to calculate the UAS trajectory, albeit at different times, and thus the two were not entirely independent. Site 2 had 18 widely distributed checkpoints, with a mean nearest neighbor distance of 22 m. Figure 7 shows the distribution of Site 2 targets. Site 2 target coordinates were established using a total station, with distances corrected for map projection distance scale. A reflective surveying prism was held above the targets so coordinates could be reduced to the target reference points. Similar to Site 1, since total station back sight observations were made to a (non-lidar-targeted) point with coordinates established by the same GNSS observations that were used for calculating the UAS trajectory, target coordinates were not completely independent from the estimated UAS trajectory.

Site 1: Target Point Precision and Accuracy, and Effect of Point Density
The extremely high point density of Site 1 allowed an analysis of the impact of point density on target point accuracy. Since the study area had flat, level terrain, the tilt parameters of the target, and , were constrained to zero. If they are solved for in the mensuration routine, a systematic error is observed at low point cloud densities. This effect is described in more detail at the end of this subsection. Example distributions of target points at different densities are shown in Figure 8. In the figure, note that the edges of the point clouds do not follow exactly the edges of the template, since they were manually segmented prior to mensuration without knowledge of the exact target location. This is entirely allowable with the mensuration algorithm, and enables aggressive segmentation to avoid the inclusion of non-target points. The average estimated target position mensuration uncertainties for various point densities, based on Equation (5), are given in Table 2. Note that based on the average number of points per target (given in Table 2), the standard deviations relative to each

Site 1: Target Point Precision and Accuracy, and Effect of Point Density
The extremely high point density of Site 1 allowed an analysis of the impact of point density on target point accuracy. Since the study area had flat, level terrain, the tilt parameters of the target, ω and φ, were constrained to zero. If they are solved for in the mensuration routine, a systematic error is observed at low point cloud densities. This effect is described in more detail at the end of this sub-section. Example distributions of target points at different densities are shown in Figure 8. In the figure, note that the edges of the point clouds do not follow exactly the edges of the template, since they were manually segmented prior to mensuration without knowledge of the exact target location. This is entirely allowable with the mensuration algorithm, and enables aggressive segmentation to avoid the inclusion of non-target points. The average estimated target position mensuration uncertainties for various point densities, based on Equation (5), are given in Table 2. Note that based on the average number of points per target (given in Table 2), the standard deviations relative to each other could be approximated by σ i ≈ σ j n j /n i where n i is the number of points on the target for spacing i, which is expected.  Mensuration uncertainty increased with an increase in point spacing, since fewer observations were included in the position estimates. The estimated standard deviations were 1 cm or less for each component for all point spacings less than 15 cm. All components had nearly identical mensuration uncertainty; however, it is important to note this was not the case when tilt angles were included as unknowns where the Z uncertainty was larger than that for the horizontal components. With a point spacing of 15 cm, the standard deviations were all 1.5 cm. This is near the upper-limit for practical application, approaching the level of independent GCP measurement uncertainty. Based on this estimate it is recommended that point spacing not exceed 10 cm for precise target mensuration.
The results for target position errors based on field-surveyed check-point coordinates are given in Table 3 and illustrated in Figure 9. The errors were not distributed significantly differently from a normal distribution based on the Shapiro-Wilk test for normality (p > 0.05). The results show that  Mensuration uncertainty increased with an increase in point spacing, since fewer observations were included in the position estimates. The estimated standard deviations were 1 cm or less for each component for all point spacings less than 15 cm. All components had nearly identical mensuration uncertainty; however, it is important to note this was not the case when tilt angles were included as unknowns where the Z uncertainty was larger than that for the horizontal components. With a point spacing of 15 cm, the standard deviations were all 1.5 cm. This is near the upper-limit for practical application, approaching the level of independent GCP measurement uncertainty. Based on this estimate it is recommended that point spacing not exceed 10 cm for precise target mensuration. The results for target position errors based on field-surveyed check-point coordinates are given in Table 3 and illustrated in Figure 9. The errors were not distributed significantly differently from a normal distribution based on the Shapiro-Wilk test for normality (p > 0.05). The results show that RMSE values increased with point spacing, and generally followed the trend of estimated mensuration uncertainty. There was not enough information to make a statistically-rigorous statement about the power of estimated mensuration uncertainty to predict target accuracy, but the results are consistent in this study. The RMSE increased with estimated mensuration uncertainty, and thus it could be considered an effective way to predict the expected accuracy of target locations within point clouds at various point densities. Mensuration uncertainty values were lower than the RMSE values, as expected due to imperfections in point location uncertainty estimation methods, as discussed in the Methods section, and due to errors in the surveyed coordinates of the check points not accounted for by mensuration uncertainty. The difference between mean square error and estimated mensuration variance was in the order of the expected field-surveyed coordinate mean squared error, approximately 1 to 2 cm 2 . Table 3. Accuracy assessment for various point spacings at Site 1.   Significant differences could not be found for the mean errors of the horizontal components based on t-tests of difference of means. The mean error in the vertical component for the 2 cm point spacing was found to be significantly different than all other point spacings (p < 0.01). The mean vertical error for the 5 cm point spacing was not significantly different to that for the 10 cm and 15 cm point spacings (p = 0.08 and p = 0.07, respectively), and the mean vertical error for the 10 cm point spacing was not significantly different to that for the 15 cm point spacing (p = 0.97). Coordinate error standard deviations always increased with point spacing. A significant difference was found between the X component variance of the 15 cm spacing coordinates and that from the 2 cm and 5 cm point spacing (p < 0.003, p < 0.02, respectively). The Z component variance was found to be significantly different for the 2 cm spacing compared to the 15 cm spacing (p < 0.05).

RMSE (m) Mean Error (m) Standard Deviation (m)
Based on these results, the developed targets provide a means to identify precise locations within UAS-lidar point clouds, and therefore are appropriate for accuracy assessment and use as conjugate features in data adjustment. For targets with the dimensions used in this study, it is recommended that average point spacing be around 5 cm and ensured to be less than 10 cm to obtain accurate Significant differences could not be found for the mean errors of the horizontal components based on t-tests of difference of means. The mean error in the vertical component for the 2 cm point spacing was found to be significantly different than all other point spacings (p < 0.01). The mean vertical error for the 5 cm point spacing was not significantly different to that for the 10 cm and 15 cm point spacings (p = 0.08 and p = 0.07, respectively), and the mean vertical error for the 10 cm point spacing was not significantly different to that for the 15 cm point spacing (p = 0.97). Coordinate error standard deviations always increased with point spacing. A significant difference was found between the X component variance of the 15 cm spacing coordinates and that from the 2 cm and 5 cm point spacing Remote Sens. 2019, 11, 3019 13 of 20 (p < 0.003, p < 0.02, respectively). The Z component variance was found to be significantly different for the 2 cm spacing compared to the 15 cm spacing (p < 0.05).
Based on these results, the developed targets provide a means to identify precise locations within UAS-lidar point clouds, and therefore are appropriate for accuracy assessment and use as conjugate features in data adjustment. For targets with the dimensions used in this study, it is recommended that average point spacing be around 5 cm and ensured to be less than 10 cm to obtain accurate within-cloud mensuration relative to check-point accuracy. This is readily achievable in practice with common scanners, such as the VLP-16, and UAS lidar products are often within this tolerance. However, it is expected that larger targets would allow for wider point spacing. The predicted uncertainty in target mensuration agreed with check-point accuracy results, and can be used in determining appropriate use of the targets, and guide mission planning to ensure adequate practical mensuration precision. As mentioned in the Materials and Methods section, 20 cm spacing led to failure of the mensuration algorithm with gross errors in 2 of the 20 targets. Removal of these two targets yielded RMSE values of 0.024 cm, 0.022 cm, and 0.044 cm for X, Y, and Z components for 20 cm spacing. Although these values may be acceptable for some applications, we did not include them in the analysis, since the algorithm had a high likelihood of failure.
It should be noted that results from initial experiments when including all three angular parameters as unknowns in the mensuration algorithm showed that when the point density decreased, the estimated centroids of the target templates fitted to the points tended to remain similar, but the angular orientations of the estimated target poses had substantially increased error, particularly for 15 cm spacing. Error in tilt of the target about the centroid of the points displaces all components, X, Y, and Z, of the reference point of the target (top of the pyramid). The Z component of the reference point is always shifted down with increased random tilt error about the centroid of the target, leading to a negative vertical bias increasing in magnitude as tilt error increases with point spacing. To allow for a fair comparison of higher point spacings in the experiments, the tilt parameters were constrained to zero for all trials at Site 1. This was viable since the ground was flat and level at the site. In practice, leveling the targets is achievable on moderately sloped ground.

Site 2: Effect of Uncertainty Estimation and Proper Weighting on Target Point Accuracy
All of the target coordinates for Site 1 were found using least squares with a weighting scheme based on the propagated uncertainty of point coordinates, as shown in Equation (3). Due to the larger area of Site 2, more points were associated with longer scanner ranges. Thus, points had larger average uncertainty, and a wider range of target point location errors and associated propagated uncertainties. This enabled an opportunity to compare the effects of including uncertainties in weighting the mensuration algorithm. The results of using the weighted and non-weighted schemes for Site 2, as presented in the Methods section, are shown in Table 4 and Figures 10 and 11. The target coordinate errors were not distributed significantly differently from a normal distribution based on the Shapiro-Wilk test for normality (p > 0.05). Table 4. Position error statistics based on ground-surveyed target coordinates at Site 2 for weighted and unweighted mensuration methods. Using the weighted method yielded lower RMSE values for each of the X, Y, and Z components of the estimated target reference point. Mean errors magnitudes (with the exception of the X component) and standard deviations of target point errors were also lower when using the weighted method. These results suggest it is advantageous to use the weighted scheme when measuring the target Remote Sens. 2019, 11, 3019 14 of 20 points and any other derived multi-point feature within the point cloud, such as corners or planar surfaces. Substantial statistical significance could not be established for the difference in variance for any components X (p = 0.12), Y (p = 0.17), Z (p = 0.11), and the mean errors were also not found to be significantly different for X (p = 0.96), Y (p = 0.47), and Z (p = 0.10). In addition to the apparent increase in precision and accuracy from using the weighted algorithm, inspection of the normalized residuals from each solution indicated that observation errors were better modeled using the weighting scheme. Figure 12 shows histograms of all normalized residuals for all the weighted and unweighted solutions. In addition to the apparent increase in precision and accuracy from using the weighted algorithm, inspection of the normalized residuals from each solution indicated that observation errors were better modeled using the weighting scheme. Figure 12 shows histograms of all normalized residuals for all the weighted and unweighted solutions. In addition to the apparent increase in precision and accuracy from using the weighted algorithm, inspection of the normalized residuals from each solution indicated that observation errors were better modeled using the weighting scheme. Figure 12 shows histograms of all normalized residuals for all the weighted and unweighted solutions. Weighted observations led to more normally distributed residuals and ignoring observation uncertainty led to higher kurtosis and skew. Higher kurtosis indicates heavy-tailed distributions and observations presenting as potential outliers relative to the expected normal distribution of their errors. While removing points with high residuals will decrease the kurtosis and potentially improve the solution, a better approach is to properly weight them to retain observation redundancy, as presented here.

Discussion
In Site 1, the estimated mensuration uncertainty and standard deviations of ground-surveyed checkpoints were less than 2 cm for point spacings up to 15 cm in all components, X, Y, and Z. To guarantee mensuration within the expected error of ground-surveyed points, it is recommended to not exceed 10 cm. Although the experiments indicated a significant decrease in mensuration accuracy in the sparsest clouds, the approach enabled target location and yielded reliable measured coordinates up to a spacing of 15 cm. When comparing against field-surveyed check point coordinates, the target coordinate accuracy decreased with decreased point density. However, even for the sparsest point clouds tested in this study, the lowest accuracy measured using RMSE was in the order of what is expected for the ground survey and the theoretical limits of the UAS lidar system, ~2 cm, illustrating that mensuration accuracy for the targets and associated algorithm were sufficient for accuracy assessment. The estimated uncertainties associated with the mensuration algorithm generally agreed with the results of the check-point accuracy experiment when taking into account scanner and check-point error, and could be used for weighting in rigorous adjustment and fusion of UAS lidar data. For example, when combining lidar data from two separate collections, the targets can be used to register the point clouds with each other via a weighted coordinate transformation to bring them into coincidence. Further exploration of the value of target position uncertainty estimation is given in the appendix.
Although useable results were achieved with and without it, the Site 2 experiment illustrated Weighted observations led to more normally distributed residuals and ignoring observation uncertainty led to higher kurtosis and skew. Higher kurtosis indicates heavy-tailed distributions and observations presenting as potential outliers relative to the expected normal distribution of their errors. While removing points with high residuals will decrease the kurtosis and potentially improve the solution, a better approach is to properly weight them to retain observation redundancy, as presented here.

Discussion
In Site 1, the estimated mensuration uncertainty and standard deviations of ground-surveyed checkpoints were less than 2 cm for point spacings up to 15 cm in all components, X, Y, and Z. To guarantee mensuration within the expected error of ground-surveyed points, it is recommended to not exceed 10 cm. Although the experiments indicated a significant decrease in mensuration accuracy in the sparsest clouds, the approach enabled target location and yielded reliable measured coordinates up to a spacing of 15 cm. When comparing against field-surveyed check point coordinates, the target coordinate accuracy decreased with decreased point density. However, even for the sparsest point clouds tested in this study, the lowest accuracy measured using RMSE was in the order of what is expected for the ground survey and the theoretical limits of the UAS lidar system,~2 cm, illustrating that mensuration accuracy for the targets and associated algorithm were sufficient for accuracy assessment. The estimated uncertainties associated with the mensuration algorithm generally agreed with the results of the check-point accuracy experiment when taking into account scanner and check-point error, and could be used for weighting in rigorous adjustment and fusion of UAS lidar data. For example, when combining lidar data from two separate collections, the targets can be used to register the point clouds with each other via a weighted coordinate transformation to bring them into coincidence. Further exploration of the value of target position uncertainty estimation is given in the Appendix A.
Although useable results were achieved with and without it, the Site 2 experiment illustrated the potentially significant impact on mensuration that can be made by including the uncertainty of each point as input in the algorithm. The RMSE was lower for each component, X, Y, and Z, when incorporating weights into the mensuration algorithm. Similarly, mean error magnitudes and standard deviations were generally lower when weighting. While rigorous point coordinate uncertainty requires navigation uncertainty metadata, which may preclude use by some practitioners, surrogate generic methods based on manufacturer-specified scanner encoder uncertainty, range uncertainty, and estimates of platform exterior orientation uncertainty may be as effective as rigorous physical modeling. Users are encouraged to explore either option for their systems, since beyond the benefits of carrying forward propagated uncertainties for post-processing and analysis, they also increase the accuracy of general location measurement, and avoid the reduced point cloud density and target-placement restriction associated with more naïve methods, such as point removal based on, e.g., long range or across-track distance.
The targets had some, albeit minor, drawbacks. For sparse point clouds, the targets generated noticeable systematic error in the form of lower Z coordinates relative to denser point clouds (and check point coordinates). This may be remedied by constraining tilt to zero, the approach used in the presented Site 1 experiment, since the ground was flat and level at the study site. However, on sloped ground one may need to level the targets when placing them or record the tilt to enable constraint to non-zero values.
Black tape was placed on the targets for use in fusion with photogrammetry. While it was effective in its purpose, drop-outs (no point returned) were observed for some points falling on the tape (some of these voids are visible in Figure 4a). Although this did not seem to have a significant effect on the mensuration of the targets in the point clouds, it could be remedied by using a different colored tape or only placing tape on targets that must be seen from imagery. Another approach the present authors utilized is the placement of flat targets on the ground next to the target (e.g., adjacent to the north facet) to obtain more reliable height, the weakest-resolved component for the targets. The targets were very light weight, since they were made from corrugated plastic. In windy conditions, they must be staked to the ground, or, depending on the ground composition, secured with sandbags on the corners.
The design of our targets precluded viewing any monument, such as capped rebar, underneath the targets when placed. To place a target associated with a monumented point, an effective approach has been to place a tripod with a leveled and centered optical-plummet tribrach over the monumented point, then place the target above it, centering the reference point of the pyramid visually using the tribrach. Height was measured using tape and checked using the design specifications of the target (0.4 m). The target points were easily and rapidly manually extracted from point clouds, and non-target and dubious points were identified and culled using tests based on Equations (6) and (7). However, for a large number of targets, manual extraction can be somewhat time consuming. Future work will include the development of automatic identification and extraction of target points to include a RANSAC-approach to discriminate target points and identify outliers.

Conclusions
Laser scanning from unoccupied aerial systems (UAS) is gaining popularity due to improved safety and lower entry costs relative to conventional large piloted systems. Compared to Structure from Motion (SfM) photogrammetry from UAS, UAS lidar offers improved processing speed and robustness of generating 3D models, and penetration of vegetation gaps for capturing terrain and production of, e.g., canopy height models. Like all geospatial data, point clouds and 3D models often require artificial targets identifiable within the products to enable accuracy assessment, calibration, adjustment, and fusion. In SfM photogrammetry, easily-manufactured flat targets with photo-identifiable markings are adequate for precisely measuring coordinates within the imagery towards 3D measurement of locations within the derived production. In UAS lidar, analogous flat intensity-based targets are often not effective, since intensity values may not allow for accurate discrimination of target points within the scene, and the distribution of target points may bias horizontal component measurements. This study introduced a new structure-based target method and associated rigorous mensuration algorithm to address this issue. The experiments presented here, with collection carried out at different sites under different flight configurations, illustrate that the developed targets are a viable method for mensuration in UAS-lidar point clouds.
The targets in this study were developed based on necessity, since the vast majority of research carried out with the UAS laser scanner involves surveys over rural, often highly vegetated areas. There are typically no artificial structures to enable mensuration, and intensity-based targeting is not effective due to the characteristics of the lidar system intensity measurements, and the reflectivity of the scenes (e.g., Figure 5). The goal in development was to create portable, multimodal (i.e., useful for photography and lidar) markers to allow reliable full 3D mensuration within UAS lidar point clouds. Beyond the presented experiments, the targets have been used in various environments, including agricultural, coastal, wetland, and diverse forested areas around Florida, and with different scanners, including the Velodyne VLP-16 LW, HDL-32E, and VLP-32C. Under typical VLP-16 scanning mission parameters (40 m or less AGL, 5-10 m/s speed, 50% nominal overlapping swaths), the targets are readily identifiable and their 3D location measurable. The only environment precluding their use has been under dense canopy, where too few points are found on the target for reliable application of the mensuration algorithm. They have been used effectively for accuracy assessment (as reported here), in situ boresight calibration, and for aligning swaths via naïve strip adjustment.  The swaths in Figure A1 have variable distributions of points for the target. Note that swath a, for example, has relatively few points on the lower left third of the target, and target d has relatively few points on the upper third. Even with these clustered distributions of points, the standard deviation of the coordinates was 0.6 cm, 1.3 cm, and 1.3 cm for X, Y, and Z, respectively. This is a sufficient level of precision for calibration and strip adjustment. The swath with the largest target location difference from the mean was b. Also, note that the estimated mensuration uncertainty for this swath was the highest of all the points. In practice, the larger estimated standard deviation would lead a user to inspect the estimated target location for this swath. The target location could be omitted, or could be included and weighted based on the estimated uncertainty. Removal of swath target b leads to standard deviations of the differences of the mean of 0.5 cm, 0.9 cm, and 0.8 cm for the X, Y, and Z components, respectively.  The swaths in Figure A1 have variable distributions of points for the target. Note that swath a, for example, has relatively few points on the lower left third of the target, and target d has relatively few points on the upper third. Even with these clustered distributions of points, the standard deviation of the coordinates was 0.6 cm, 1.3 cm, and 1.3 cm for X, Y, and Z, respectively. This is a sufficient level of precision for calibration and strip adjustment. The swath with the largest target location difference from the mean was b. Also, note that the estimated mensuration uncertainty for this swath was the highest of all the points. In practice, the larger estimated standard deviation would lead a user to inspect the estimated target location for this swath. The target location could be omitted, or could be included and weighted based on the estimated uncertainty. Removal of swath target b leads to standard deviations of the differences of the mean of 0.5 cm, 0.9 cm, and 0.8 cm for the X, Y, and Z components, respectively.